Bean Curve

Author

Eric W. Weisstein
November 5, 2021

This notebook downloaded from http://mathworld.wolfram.com/notebooks/PlaneCurves/BeanCurve.nb.

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

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

Sources

Cundy and Rowlett

Properties

Equations

Cartesian

eqn=x^4+x^2 y^2+y^4-x(x^2+y^2);

Factor[Expand[eqn]]

AlgebraicDegree[eqn,{x,y}]

Plot

V6

ContourPlot[x^4+x^2 y^2+y^4-x(x^2+y^2)==0,{x,-.1,1},{y,-.7,.7},ContourStyle->Red,PlotPoints->20,AspectRatio->Automatic,ImageSize->{Automatic,200},Axes->True]

RegionPlot[x^4+x^2 y^2+y^4<x(x^2+y^2),{x,0,1},{y,-.7,.7},AspectRatio->Automatic,ImageSize->Small]

Solving for y(x)

ysoln=Union[y/.Simplify[Solve[x^4+x^2 y^2+y^4==x(x^2+y^2),y]]]

Plot[Evaluate[ysoln],{x,0,1},PlotStyle->List/@{Red,Yellow,Green,Blue},AspectRatio->Automatic]

Solving for x(y)

ysoln=x/.Simplify[Solve[x^4+x^2 y^2+y^4==x(x^2+y^2),x]]

ParametricPlot[Evaluate[{#,y}&/@xsoln],{y,-1,1},PlotStyle->List/@{Red,Yellow,Green,Blue},
AspectRatio->Automatic]

Solve[soln[[3]]==soln[[4]],y]

soln[[4]]/.y->2/3

Simplify[x^4+x^2 y^2+y^4==x(x^2+y^2)/.{x->2/3,y->2/3}]

ArcLength

Entity["PlaneCurve","BeanCurve"]["ArcLength"]

Perimeter

$Version

Perimeter[Region[ImplicitRegion[x^4+x^2 y^2+y^4<x(x^2+y^2),{x,y}]]]//Timing

Perimeter[Region[ImplicitRegion[x^4+x^2 y^2+y^4<x(x^2+y^2),{x,y}]],PrecisionGoal->120,WorkingPrecision->250]//Timing

Precision[%[[2]]]

Arc Length y(x)

Plot[Evaluate[(Sqrt[1+((D[((Sqrt[x-(x^2)+(Sqrt[(x^2)+2 (x^3)-3 (x^4)])])/(Sqrt[2])),x])^2)])],{x,0,1}]

$Version

Integrate[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]),{x,0,1}]//Timing

Integrate[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]),x]//Timing

2NIntegrate[(Sqrt[1+((D[((Sqrt[x-(x^2)+(Sqrt[(x^2)+2 (x^3)-3 (x^4)])])/(Sqrt[2])),x])^2)]),{x,0,1},WorkingPrecision->220,PrecisionGoal->110]

2NIntegrate[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]),{x,0,1}]

Extract the leading singular (but integrable) terms at the ends, integrate them exactly and NIntegrate the remaining singularity-free function (from mtrott 2020-07-23):

FullSimplify[Normal[Series[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]),{x,0,1}]],0<x<1]//Apart

Simplify/@(FullSimplify[Normal[Series[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]),{x,1,0}]],0<x<1]//Apart)

Integrate[(1/(4 (((1-x))^(3/4))))+(3/(16 (((1-x))^(1/4))))+(1/(2 (Sqrt[x]))),{x,0,1}]

Plot[(Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))])-(1/(4 (((1-x))^(3/4))))-(3/(16 (((1-x))^(1/4))))-(1/(2 (Sqrt[x]))),{x,0,1}]

(9/2)+NIntegrate[2 (Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))])-(1/(2 (((1-x))^(3/4))))-(3/(8 (((1-x))^(1/4))))-(1/(Sqrt[x])),{x,0,1},WorkingPrecision->220,PrecisionGoal->110,MaxRecursion->50]

(9/2)+Integrate[2 (Sqrt[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))])-(1/(2 (((1-x))^(3/4))))-(3/(8 (((1-x))^(1/4))))-(1/(Sqrt[x])),{x,0,1}]//Timing

(HyperlinkTemplate[A193506, https://oeis.org/A193506None, https://oeis.org/A193506, HyperlinkActionRecycled, HyperlinkActive, , ])

ExpandAll[1+((((1-2 x+((1+3 x-6 (x^2))/(Sqrt[1+(2-3 x) x]))))^2)/(8 x (1-x+(Sqrt[1+(2-3 x) x]))))]//FullSimplify

2NIntegrate[(Sqrt[((1+(Sqrt[1+(2-3 x) x])+x (5+4 (Sqrt[1+(2-3 x) x])+2 x (1-8 (Sqrt[1+(2-3 x) x])+2 x (-5+3 x+3 (Sqrt[1+(2-3 x) x])))))/(8 (((-1+x))^2) x (1+3 x)))]),{x,0,1}]

2Integrate[(Sqrt[((1+(Sqrt[1+(2-3 x) x])+x (5+4 (Sqrt[1+(2-3 x) x])+2 x (1-8 (Sqrt[1+(2-3 x) x])+2 x (-5+3 x+3 (Sqrt[1+(2-3 x) x])))))/(8 (((-1+x))^2) x (1+3 x)))]),{x,0,1}]//Timing

dsx=Simplify[(Sqrt[1+((D[xsoln[[3]],x])^2)]),0<x<1]

FullSimplify[dsx]

FullSimplify[dsx,0<x<1]//Timing

dsx2=FullSimplify[Together//@dsx]

2NIntegrate[dsx,{x,0,1},WorkingPrecision->20,MaxRecursion->15,SingularityDepth->10]

Integrate[(1/2) √(((Sqrt[(x^2)+2 (x^3)-3 (x^4)])+x (1+5 (Sqrt[(x^2)+2 (x^3)-3 (x^4)])+x (6-4 (Sqrt[(x^2)+2 (x^3)-3 (x^4)])+x (-1+4 x (-7+6 x)))))/((-1+x) x (1+3 x) ((-1+x) x-(Sqrt[(x^2)+2 (x^3)-3 (x^4)])))),{x,0,1}]//Timing

V4.2

2Integrate[dsx,{x,0,1}]//Timing

V5.1

2Integrate[dsx,{x,0,1}]//Timing

V6

2Integrate[dsx,{x,0,1}]//Timing

2Integrate[dsx2,{x,0,1}]//Timing

ArcLength x(y)

dsy=(Sqrt[1+((D[#,y])^2)])&/@ysoln[[{3,4}]]

FullSimplify[dsy]//Timing

Numerical

NIntegrate[#,{y,-2/3,2/3}]&/@dsy

Plus@@%

2Plus@@(NIntegrate[#,{y,0,2/3}]&/@dsy)

V4.2

2Plus@@(Integrate[#,{y,0,2/3}]&/@dsy)//Timing

...hangs...

V5.1

2Plus@@(Integrate[#,{y,0,2/3}]&/@dsy)//Timing

V6

2Plus@@(Integrate[#,{y,0,2/3}]&/@dsy)//Timing

V13

$Version

2Total[Integrate[#,{y,0,2/3}]&/@{(Sqrt[1+(1/16) (((((4 y)/(Sqrt[1+4 (y^2)]))+((2 (Sqrt[2]) y (-3+(1/(Sqrt[1+4 (y^2)]))))/(Sqrt[1-6 (y^2)+(Sqrt[1+4 (y^2)])]))))^2)]),(Sqrt[1+(1/16) (((((4 y)/(Sqrt[1+4 (y^2)]))-((2 (Sqrt[2]) y (-3+(1/(Sqrt[1+4 (y^2)]))))/(Sqrt[1-6 (y^2)+(Sqrt[1+4 (y^2)])]))))^2)])}]//Timing

ResourceFunction["IntegrateAlgebraic"]

ResourceFunction["IntegrateAlgebraic"][#,y]&/@{(Sqrt[1+(1/16) (((((4 y)/(Sqrt[1+4 (y^2)]))+((2 (Sqrt[2]) y (-3+(1/(Sqrt[1+4 (y^2)]))))/(Sqrt[1-6 (y^2)+(Sqrt[1+4 (y^2)])]))))^2)]),(Sqrt[1+(1/16) (((((4 y)/(Sqrt[1+4 (y^2)]))-((2 (Sqrt[2]) y (-3+(1/(Sqrt[1+4 (y^2)]))))/(Sqrt[1-6 (y^2)+(Sqrt[1+4 (y^2)])]))))^2)])}//Timing

ResourceFunction["IntegrateAlgebraic"][(Sqrt[1+(1/16) (((((4 y)/(Sqrt[1+4 (y^2)]))+((2 (Sqrt[2]) y (-3+(1/(Sqrt[1+4 (y^2)]))))/(Sqrt[1-6 (y^2)+(Sqrt[1+4 (y^2)])]))))^2)]),y]

Area

Entity["PlaneCurve","BeanCurve"]["Area"][1]

N[%]

$Version

NIntegrate[Boole[x^4+x^2 y^2+y^4<x(x^2+y^2)],{x,-∞,∞},{y,-∞,∞},WorkingPrecision->20]

Integrate[Boole[x^4+x^2 y^2+y^4<x(x^2+y^2)],{x,-∞,∞},{y,-∞,∞}]//Timing

AreaMomentOfInertiaTensor

Entity["PlaneCurve","BeanCurve"]["AreaMomentOfInertiaTensor"][1]

j=NIntegrate[{{(y^2),-x y },{-x y,(x^2)}}Boole[x^4+x^2 y^2+y^4<x(x^2+y^2)],{x,-∞,∞},{y,-∞,∞},WorkingPrecision->20]

(j=Integrate[{{(y^2),-x y },{-x y,(x^2)}}Boole[x^4+x^2 y^2+y^4<x(x^2+y^2)],{x,-∞,∞},{y,-∞,∞}])//Timing

N[j]

Centroid

Entity["PlaneCurve","BeanCurve"]["Centroid"][1]

RegionPlot[-(x^3)+(x^4)-x (y^2)+(x^2) (y^2)+(y^4)<0,{x,0,1},{y,-.7,.7},AspectRatio->Automatic,ImageSize->150]

$Version

((Integrate[{x,y}Boole[x^4+x^2 y^2+y^4<x(x^2+y^2)],{x,-∞,∞},{y,-∞,∞}])/((7 π)/(12 (Sqrt[3]))))//Timing

(1/((7 π)/(12 (Sqrt[3]))))2Integrate[x((Sqrt[x (1-x+(Sqrt[1+(2-3 x) x]))])/(Sqrt[2])),{x,0,1}]//Timing

Factor[1+(2-3 x) x]

Expand[1+(2-3 x) x]

Expand[-3(((x-(1/3)))^2)+(4/3)]

(4/3)(1-(9/4)(((x-(1/3)))^2))//Expand

Solve[(9/4)(((x-(1/3)))^2)==Cos[t]^2,x]

FullSimplify[((Sqrt[x (1-x+(Sqrt[1+(2-3 x) x]))])/(Sqrt[2]))/.x->(1/3) (1+2 Cos[t])]

D[(1/3) (1+2 Cos[t]),t]

Solve[(1/3) (1+2 Cos[t])==0,t]

Solve[(1/3) (1+2 Cos[t])==1,t]

Plot[(1/3)(1/3) (1+2 Cos[t]) (Sqrt[(1+2 Cos[t]) (1-Cos[t]+(Sqrt[3]) (Sqrt[((Sin[t])^2)]))])((2 Sin[t])/3),{t,0,2π/3}]

(2/area)NIntegrate[(1/3) (1/3) (1+2 Cos[t])(Sqrt[(1+2 Cos[t]) (1-Cos[t]+(Sqrt[3]) (Sqrt[((Sin[t])^2)]))])((2 Sin[t])/3),{t,0,2π/3}]

FullSimplify[(2/area) (1/3) (1/3)(1+2 Cos[t])(Sqrt[(1+2 Cos[t]) (1-Cos[t]+(Sqrt[3]) (Sqrt[((Sin[t])^2)]))])((2 Sin[t])/3),0<t<2π/3]

((16(Sqrt[3]) )/(63 π))NIntegrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1- Cos[t]+ (Sqrt[3]) Sin[t]]),{t,0,2π/3}]

(xbar=((16(Sqrt[3]) )/(63 π))Integrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1- Cos[t]+ (Sqrt[3]) Sin[t]]),{t,0,2π/3}])//Timing

Integrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1- Cos[t]+ (Sqrt[3]) Sin[t]]),t]//Timing

TrigReduce[1- Cos[t]+ (Sqrt[3]) Sin[t]]

Subtract@@{-Cos[t]+(Sqrt[3]) Sin[t],-2Cos[t+ArcTan[(Sqrt[3])]]}//FullSimplify

((16(Sqrt[3]) )/(63 π))NIntegrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1-2Cos[t+ArcTan[(Sqrt[3])]]]),{t,0,2π/3}]

i=((16(Sqrt[3]) )/(63 π))Integrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1-2Cos[t+ArcTan[(Sqrt[3])]]]),t]

Plot[i,{t,0,2π/3}]

(xbar=((16(Sqrt[3]) )/(63 π))Integrate[(((1+2 Cos[t]))^(3/2)) Sin[t] (Sqrt[1-2Cos[t+ArcTan[(Sqrt[3])]]]),{t,0,2π/3}])//Timing

Curvature

SortBy[FullSimplify[Curvature[eqn,{x,y}],(x|y)∈Reals],ByteCount]

ByteCount/@%

InflectionPoints

InflectionPoints[eqn,{x,y}]

SingularPoints

pts=SingularPoints[eqn,{x,y}]

Reduce[eqn==D[eqn,x]==D[eqn,y]==0,{x,y},Reals]

Slope

FullSimplify[Slope[x^4+x^2 y^2+y^4==x (x^2+y^2),{x,y}]]

FullSimplify[Slope[x^4+x^2 y^2+y^4==x (x^2+y^2),{x,y}]]