Heart Surface

Author

Eric W. Weisstein
June 3, 2008

This notebook downloaded from http://mathworld.wolfram.com/notebooks/Surfaces/HeartSurface.nb.

For more information, see Eric's MathWorld entry http://mathworld.wolfram.com/HeartSurface.html.

©2008 Wolfram Research, Inc. except for portions noted otherwise

Additional Authors

Nice 3-D rendered version by
Michael Trott
January 25, 2002

Nordstrand's Heart

ContourPlot3D[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3,{x,-2,2},{y,-2,2},{z,-2,2},MaxRecursion->6,Mesh->None,Contours->{0},Boxed->False,PlotRange->All,Axes->False,ViewPoint->{2,.1,.5},ContourStyle->Red]

Extent

Through[{MinValue,MaxValue}[{#,(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3==0},{x,y,z}]]&/@{x,y,z}

N[%]

Equation

Factor[Numerator[Together[eqn=(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3]]]

Properties

ByteCount/@(k=GaussianCurvature[eqn,{x,y,z}])//Timing

k[[-1]]

FullSimplify[k[[-1]]]

ByteCount[%]

Volume

Cartesian

NIntegrate[Boole[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3<0],{x,-2,2},{y,-2,2},{z,-2,2}]//Timing

NIntegrate[Boole[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3<0],{x,-∞,∞},{y,-∞,∞},{z,-∞,∞}]//Timing

Integrate[Boole[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3<0],{x,-∞,∞},{y,-∞,∞},{z,-∞,∞}]

Cylindrical

FullSimplify[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3/.{x->r Cos[θ],y->r Sin[θ]},r>0&&0<θ<2π]

Through[{MinValue,MaxValue}[{#,(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3==0},{x,y,z}]]&/@{(Sqrt[x^2+y^2]),z}

NIntegrate[Boole[-(8/45) (r^2) (z^3) ((Cos[θ])^2)+(((-1+((25 (r^2))/9)+(z^2)+(7/9) (r^2) Cos[2 θ]))^3)-(r^2) (z^3) ((Sin[θ])^2)<0]r,{θ,0,2π},{r,0,(Sqrt[-Root[64+384 #1+769 ((#1)^2)+512 ((#1)^3)&,1]])},{z,-1,-Root[27-54 ((#1)^2)+2 ((#1)^3)+27 ((#1)^4)&,1]}]//Timing

Integrate[Boole[-(8/45) (r^2) (z^3) ((Cos[θ])^2)+(((-1+((25 (r^2))/9)+(z^2)+(7/9) (r^2) Cos[2 θ]))^3)-(r^2) (z^3) ((Sin[θ])^2)<0]r,{θ,0,2π},{r,0,∞},{z,-∞,∞}]//Timing

Integrate[Boole[-(8/45) (r^2) (z^3) ((Cos[θ])^2)+(((-1+((25 (r^2))/9)+(z^2)+(7/9) (r^2) Cos[2 θ]))^3)-(r^2) (z^3) ((Sin[θ])^2)<0]r,{θ,0,2π},{z,-∞,∞},{r,0,∞}]//Timing

Spherical

FullSimplify[(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3/.{x->r Cos[θ]Sin[ϕ],y->r Sin[θ]Sin[ϕ],z->r Cos[ϕ]},r>0&&0<θ<2π&&0<ϕ<π]

Through[{MinValue,MaxValue}[{#,(2(4/3x)^2+2y^2+z^2-1)^3-(4/3x)^2z^3/10-y^2z^3==0},{x,y,z}]]&/@{(Sqrt[x^2+y^2+z^2])}

N[%]

NIntegrate[Boole[-(8/45) (r^5) ((Cos[θ])^2) ((Cos[ϕ])^3) ((Sin[ϕ])^2)-(r^5) ((Cos[ϕ])^3) ((Sin[θ])^2) ((Sin[ϕ])^2)+(1/729) (((-9+17 (r^2)+(r^2) (-8 Cos[2 ϕ]+7 Cos[2 θ] ((Sin[ϕ])^2))))^3)<0](r^2)Sin[ϕ],{θ,0,2π},{r,0,1.21},{ϕ,0,π}]//Timing

Integrate[Boole[-(8/45) (r^5) ((Cos[θ])^2) ((Cos[ϕ])^3) ((Sin[ϕ])^2)-(r^5) ((Cos[ϕ])^3) ((Sin[θ])^2) ((Sin[ϕ])^2)+(1/729) (((-9+17 (r^2)+(r^2) (-8 Cos[2 ϕ]+7 Cos[2 θ] ((Sin[ϕ])^2))))^3)<0](r^2)Sin[ϕ],{θ,0,2π},{r,0,∞},{ϕ,0,π}]//Timing

Kuska's Heart

comp.soft-sys.math.mathematica posting
Jens-Peer Kuska (kuska@informatik.uni-leipzig.de)
2004-02-13 19:09:18 PST

ContourPlot3D[(2*x^2 + y^2 + z^2 - 1)^3 - (1/10)*x^2*z^3 - y^2*z^3, {x, -0.9, 0.9}, {y, -1.2, 1.2}, {z, -1.2, 1.4},Boxed->False,Contours->{0},ContourStyle->Red,ViewPoint->{2,.1,.5},Axes->False,Mesh->None,BoxRatios->{.6,.6,.8}]

Taubin's Heart

With[{r=1},ContourPlot3D[(x^2+(3/2)^2 y^2+z^2-1)^3-x^2z^3-(3/2)^2/20 y^2z^3==0,{x,-1.5,1.5},{y,-1,1},{z,-1,1.5},Boxed->False,ContourStyle->Red,Axes->None,PlotRange->All,ViewPoint->{1,1,.2},Mesh->None,BoxRatios->{1,.5,1}]
]

Export to POV-Ray

<<Utilities`POVray`

POVCoefficients[N[(x^2+(3/2)^2 y^2+z^2-1)^3-(x^2)z^3-(3/2)^2/20 y^2z^3],{x,y,z}]

The (He)art of Michael Trott

Michael Trott
January 25, 2002

Let's load some relevant packages.

Needs["Graphics`ContourPlot3D`"]
Needs["Graphics`Graphics3D`"]

Off[General::spell]; Off[General::spell1];

This is the function. (I prefer algabraic functions in canonical form).

f[{x_,y_,z_}]=320 Factor[(x^2+(3/2)^2 y^2+z^2-1)^3-(x^2)z^3-(3/2)^2/20 y^2z^3]

The y=0 cross–section shows the heart shape nicely.

ContourPlot[Evaluate[f[{x, 0, z}]], {x, -1.2, 1.2}, {z, -1.1, 1.3}, 
            PlotPoints -> 120, Contours -> {0}];

ContourPlot can be used for a first quick visualization.

ContourPlot3D[Evaluate[f[{x, y, z}]], {x, -3/2, 3/2}, {y, -1, 1}, {z, -3/2, 3/2},
             PlotPoints -> {20, 5}, MaxRecursion -> 0, Axes -> True];

For a more detailed plot, having say y=y(x,z) is easier to draw. In the y=0 plane the x–extension is maximal. So let's get there the bounds for the projection of the heart in the x,z–plane.

solZ  = Solve[f[{x, 0, z}] == 0, z];

Plot[Evaluate[{Re @ solZ[[1, 1, 2]], Re @ solZ[[2, 1 ,2]]}], {x, 10^-8, 2},
     AspectRatio -> Automatic, PlotStyle -> {Hue[0], Hue[0.78]}];

So we can define the following limits for z=z(x).

zMin[x_] = solZ[[1, 1, 2]];
zMax[x_] = solZ[[2, 1, 2]];

Now we need the largest possible x–value. We get this from implicit differentiation.

eqD = D[f[{x[z], 0, z}], z] /. x'[z] -> 0 /. x[z] -> x

(* eliminate z *)
res = Resultant[eqD, f[{x, 0, z}], z] // Factor

xMax = Select[x /. Solve[res[[-1]] == 0, x], Im[#] == 0 ∧ # > 0&][[1]] // Simplify

N[xMax]

It remains to express y as a single valued function over the x,z-plane.

yP[{x_, z_}] = Solve[f[{x, y, z}] == 0, y][[2, 1, 2]]

We we can calculate some points on the surface.

ppx = 120; ppz = 120; ε = 10.^-6;
data = Table[{x, Re @ yP[{x, z}], z},
             {z, zMin[x], zMax[x], (zMax[x] - zMin[x])/ppz},
             {x, ε, 1. xMax, (xMax - ε)/ppx}];

polys = Cases[ListSurfacePlot3D[data], _Polygon, ∞];

We now mirror in all equaivalent positions and have the heart.

heartGraphics3D = 
Show[Graphics3D[{SurfaceColor[Hue[0], Hue[0.12], 2.3], EdgeForm[],
     {polys, Map[{1, -1, 1}#&, polys, {-2}],
             Map[{-1, 1, 1}#&, polys, {-2}],
             Map[{-1, -1, 1}#&, polys, {-2}]}}], Boxed -> False];//Timing

Here the heart is cut.

Show[heartGraphics3D, PlotRange -> {All, {0, 3/4}, All}, ViewPoint -> {0, -4, 0}];//Timing

We could now color it according to Gauss curvature (would need more plotpoints).

GaussCurvature[F_, {x_, y_, z_}] := With[{d = Derivative},
(-(d[0, 0, 2][F][x, y, z] d[0, 1, 0][F][x, y, z] d[1, 0, 0][F][x, y, z] + 
   d[0, 0, 1][F][x, y, z] (-d[0, 1, 1][F][x, y, z] d[1, 0, 0][F][x, y, z] - 
   d[0, 1, 0][F][x, y, z] d[1, 0, 1][F][x, y, z] + 
   d[0, 0, 1][F][x, y, z] d[1, 1, 0][F][x, y, z]))^2 + 
  (d[0, 0, 2][F][x, y, z] d[0, 1, 0][F][x, y, z]^2 + 
   d[0, 0, 1][F][x, y, z] (-2 d[0, 1, 0][F][x, y, z] d[0, 1, 1][F][x, y, z] + 
   d[0, 0, 1][F][x, y, z] d[0, 2, 0][F][x, y, z])) 
  (d[0, 0, 2][F][x, y, z] d[1, 0, 0][F][x, y, z]^2 + 
   d[0, 0, 1][F][x, y, z] (-2 d[1, 0, 0][F][x, y, z] d[1, 0, 1][F][x, y, z] + 
   d[0, 0, 1][F][x, y, z] d[2, 0, 0][F][x, y, z])))/
  (d[0, 0, 1][F][x, y, z]^2 (d[0, 0, 1][F][x, y, z]^2 + 
   d[0, 1, 0][F][x, y, z]^2 + d[1, 0, 0][F][x, y, z]^2)^2)]

Timing[κ[{x_, y_, z_}] = 
 GaussCurvature[Function[{x, y, z}, Evaluate[f[{x, y, z}]]], {x, y, z}] // Simplify]

κC = Compile[{x, y, z}, Evaluate[N[κ[{x, y, z}]]]];

Show[heartGraphics3D /. Polygon[l_] :> {SurfaceColor[#, #, 2.3]&[
                      Hue[ArcTan[κC @@ (Plus @@ l/4)]/(Pi/2)]], Polygon[l]},
     ViewPoint -> {3, -6, -1}];//Timing

Or we could color it according to direction.

φCoord[{x_, y_, z_}] := ArcTan[x, y];
ϑCoord[{x_, y_, z_}] := ArcCos[z/Sqrt[x^2 + y^2 + z^2]];

Show[heartGraphics3D /. Polygon[l_] :> {SurfaceColor[#, #, 2.3]&[
                      Hue[Re[SphericalHarmonicY[7, 3, ϑCoord[#], φCoord[#]]]]&[
                      Plus @@ l/4]], Polygon[l]}];//Timing

Here is a semi–transparent heart (we get a neat Moire pattern).

Show[heartGraphics3D /. Polygon[l_] :> {SurfaceColor[Hue[0.78], Hue[0.22], 2.8], 
                      Polygon[Plus @@@ Partition[Append[l, First[l]], 2, 1]/2]}];

Here is another possibility to make the heart graphics: We calculate r=r(ϕ,θ) in spherical coordinates and use parametric plotting. Because the resulting r is the solution of a sixth degree equation this has to be done numerically careful.

solϕθ = Solve[f[{r Cos[ϕ] Sin[θ], r Sin[ϕ] Sin[θ], r Cos[θ]}] == 0, r];

solϕθ /. {ϕ -> 0.4, θ -> 1.5}

frSpherical[{ϕ_, θ_}] = solϕθ[[2, 1, 2]]

frSphericalN[{ϕ_, θ_}] := frSpherical[{SetPrecision[ϕ, 20], SetPrecision[θ, 20]}]

ε = 10^-6;
p3d = ParametricPlot3D[frSphericalN[{ϕ, θ}] {Cos[ϕ] Sin[θ], Sin[ϕ] Sin[θ], Cos[θ]},
                 {ϕ, 0, Pi/2}, {θ, ε, Pi - ε}, PlotPoints -> 90, Compiled -> False];

polysϕθ = Cases[p3d, _Polygon, ∞];

Show[Graphics3D[{SurfaceColor[Hue[0], Hue[0.12], 2.3], EdgeForm[],
     {polysϕθ, Map[{1, -1, 1}#&, polysϕθ, {-2}],
               Map[{-1, 1, 1}#&, polysϕθ, {-2}],
               Map[{-1, -1, 1}#&, polysϕθ, {-2}]}}], Boxed -> False];