Gibbs-Bagdasaryan-Baez reduction

Version 2:  10 February 2015:  Added checks that points W,X,Y are inside the triangle B'

This notebook analyses a “universal cover” devised by Gibbs, Bagdasaryan and Baez for Lebesgue’s Universal Covering Problem.

See http://arxiv.org/abs/1502.01251 for details.

Original, unrotated hexagon of width 1, circumscribed radius (1/(Sqrt[3])).  Vertices are listed in clockwise order starting from (0,(1/(Sqrt[3]))).

hrad=(1/(Sqrt[3]));

H1=Table[hrad{Sin[((π i)/3)],Cos[((π i)/3)]},{i,0,5}]

Rotated hexagon:  the original hexagon, rotated counterclockwise, by (π/6)+σ.  Here we use s=sin(σ) and c=cos(σ) rather than explicit trigonometric functions.

rotH[s_,c_]=Map[{{c,-s},{s,c}}.#&,TrigExpand[Table[hrad{Sin[((π i)/3)-(π/6)],Cos[((π i)/3)-(π/6)]},{i,0,5}]]]

csAssume=0≤ s<Sin[π/6]&&Cos[π/6]<c≤1

Routines to tidy up expressions containing cosines and sines

cst0[x_]:=(MapAll[ExpandAll,x]//.{(s^(n_))/;Mod[n,2]==0:>(((1-(c^2)))^(n/2)),(s^(n_))/;Mod[n,2]==1:>s(((1-(c^2)))^((n-1)/2))})

Map[cst0,Table[(s^i),{i,1,10}]]

cst[x_,n_]:=Nest[cst0,x,n]

cstNumerator[x_,n_]:=Module[{t},t=Together[x];Factor[cst[Numerator[t],n],Extension->Automatic]/Simplify[Denominator[t]]]

Point of intersection of two line segments, with endpoints a, b for first line and c, d for second.

ils[{ax_,ay_},{bx_,by_},{cx_,cy_},{dx_,dy_}]=
Module[{p1,p2},
p1=λ {ax,ay}+(1-λ){bx,by};
p2=μ{cx,cy}+(1-μ){dx,dy};
Simplify[p1/.Solve[p1==p2,{λ,μ}][[1]]]]

Points of intersection between the two hexagons.  Each side of the original hexagon intersects both the same side of the counterclockwise-rotated hexagon, and the “next” side in clockwise order.

pih[s_,c_]=Module[{ip,j,jp,rh},
rh=rotH[s,c];
Flatten[Table[
ip=Mod[i+1,6,1];
j=Mod[i+k,6,1];
jp=Mod[j+1,6,1];
Simplify[cst[Simplify[ils[H1[[i]],H1[[ip]],rh[[j]],rh[[jp]]],csAssume],2]],{i,1,6},{k,0,1}],1]];

The vertices of the unrotated hexagon are named (Subscript[A, 1]),(Subscript[B, 1]),..., (Subscript[F, 1]), while the first intersection point clockwise from each vertex is (Subscript[A, 2]),(Subscript[B, 2]),..., (Subscript[F, 2]) and the first intersection point counterclockwise from each vertex is (Subscript[A, 3]),(Subscript[B, 3]),...(Subscript[F, 3]).

Some of these are specially named in the paper:  point N = (Subscript[E, 3]) and point O=(Subscript[C, 2]) are the 8th and 5th of the intersection points.

Npt[s_,c_]=pih[s,c][[8]]

Opt[s_,c_]=pih[s,c][[5]]

letters={"A","B","C","D","E","F"};

H1Labels=Map[Subscript[#,1]&,letters]

pihLabels=Table[Subscript[letters[[Mod[Floor[(i+2)/2],6,1]]],Mod[i,2,1]+1],{i,1,12}]/.{(Subscript["E", 3])->(Subscript["N=E", 3]),(Subscript["C", 2])->(Subscript["O=C", 2])}

drawHexagons[σ_]:=Module[{s,c,pih0},
s=Sin[σ];
c=Cos[σ];
pih0=pih[s,c];
Graphics[
{Table[
{Line[{H1[[i]],H1[[Mod[i+1,6,1]]]}],
Line[{pih0[[2*i]],pih0[[Mod[2*i+1,12,1]]]}],
Text[H1Labels[[i]],1.1H1[[i]]]
}
,{i,1,6}],
Table[Text[pihLabels[[i]],1.1pih0[[i]]],{i,1,12}]
},BaseStyle->18,ImageSize->{500,500}
]]

drawHexagons[1.3 Degree]

Area of an arbitrary polygon

areaPoly[p_]:=Module[{n,p1,p2},
n=Length[p];
Sum[
p1=p[[i]];
p2=p[[Mod[i+1,n,1]]];
(p1[[2]]p2[[1]]-p1[[1]]p2[[2]])/2,
{i,1,n}]
]

Area between a chord of length c and a circle of radius 1

chordArea[c_]=ArcSin[(c/2)]-(c/2)(Sqrt[1-((((c/2)))^2)])

Find the points on the original hexagon where circles of radius 1 centred at O and N are tangent to the hexagon.

Otan[s_,c_]=Simplify[H1[[1]]+Opt[s,c]-H1[[3]]]

Ntan[s_,c_]=Simplify[H1[[1]]+Npt[s,c]-H1[[5]]]

Find the intersection of two circles of radius 1 with specified centres a and b.  We single out the point of intersection that lies above the line from a to b if a is on the left.

intC=icirc[{ax_,ay_},{bx_,by_}]=Module[{diff,midpoint,magdiffsq},
diff={bx,by}-{ax,ay};
midpoint=({ax,ay}+{bx,by})/2;
magdiffsq=diff.diff;
midpoint+(Sqrt[(1/magdiffsq)-(1/4)]){-diff[[2]],diff[[1]]}
]

Verify solution

Simplify[(intC-{ax,ay}).(intC-{ax,ay})]

Simplify[(intC-{bx,by}).(intC-{bx,by})]

Verify with a simple case that we have the correct solution, above the line from a to b

icirc[{-(1/2),0},{(1/2),0}]

Point M is the midpoint of the fourth side of the unrotated hexagon.

Mpt=(H1[[4]]+H1[[5]])/2

Normal to the mirror-line of symmetry that passes through point M

Mnorm=Normalize[H1[[5]]-H1[[4]]]

Reflection in the mirror-line of symmetry that passes through point M

refM=IdentityMatrix[2]-2 Outer[Times,Mnorm,Mnorm];

Point Q is the reflection of the 10th intersection point between the hexagons in the mirror-line.

Qpt[s_,c_]=Simplify[refM.pih[s,c][[10]],csAssume]

Point W is the intersection of circles of radius 1 centred at M and O

Wpt0[s_,c_]=cst[Simplify[icirc[Mpt,Opt[s,c]]],3]/.{Sqrt[x_]:>Sqrt[Together[x]]}

α0=√(((-37-2 (Sqrt[3]) c+26 (c^2)+10 s-34 (Sqrt[3]) c s)/(-11+2 (Sqrt[3]) c+6 (c^2)-10 s+2 (Sqrt[3]) c s)));

α1=Sqrt[Map[Collect[#,c]&,FullSimplify[Together[cst[((α0)^2),2]]]]]

Wpt1[s_,c_]=Collect[Wpt0[s,c]/.{α0->α},α,Simplify]

Wpt[s_,c_]=Wpt1[s,c]/.α->α1

Point Y is the intersection of circles of radius 1 centred at N and Q

Ypt0[s_,c_]=cst[Simplify[icirc[Npt[s,c],Qpt[s,c]]],3]/.{Sqrt[x_]:>Sqrt[Together[x]]}

β0=(Sqrt[((3+(Sqrt[3]) c-3 (c^2)+3 s-5 (Sqrt[3]) c s)/(3-(Sqrt[3]) c-(c^2)-3 s+(Sqrt[3]) c s))]);

β1=Sqrt[Map[Collect[#,c]&,FullSimplify[Together[cst[((β0)^2),2]]]]]

Ypt1[s_,c_]=Collect[Ypt0[s,c]/.{β0->β},β,Simplify]

Ypt[s_,c_]=Ypt1[s,c]/.β->β1

Point X is the intersection of circles of radius 1 centred at N and O

Xpt0[s_,c_]=cst[Simplify[icirc[Npt[s,c],Opt[s,c]]],3]/.{Sqrt[x_]:>Sqrt[Together[x]]}

γ0=(Sqrt[((12+5 (Sqrt[3]) c-44 (c^2)-4 (Sqrt[3]) (c^3)+32 (c^4))/(6-5 (Sqrt[3]) c-4 (c^2)+4 (Sqrt[3]) (c^3)))]);

γ1=Sqrt[Map[Collect[#,c]&,FullSimplify[Together[cst[((γ0)^2),2]]]]]

Xpt1[s_,c_]=Collect[Xpt0[s,c]/.{γ0->γ},γ,Simplify[cst[Simplify[#],3]]&]

Xpt[s_,c_]=Xpt1[s,c]/.γ->γ1

Dot product between the vectors YW and YQ to determine whether the angle ∠ WYQ is greater than 90 degrees, i.e. this needs to be negative for a valid cover.

dp1[s_,c_]=Collect[Collect[FullSimplify[cst[FullSimplify[Together[(Wpt1[s,c]-Ypt1[s,c]).(Qpt[s,c]-Ypt1[s,c])]],3]],{α,β},FullSimplify]/.{(β^2)->((β1)^2)},{α,β},FullSimplify]

dp[s_,c_]=dp1[s,c]/.{α->α1,β->β1};

Dot product between the vectors WM and WY to determine whether the angle the angle ∠ MWY is greater than 90 degrees, i.e. this needs to be negative for a valid cover.

dp2[s_,c_]=Together[(Mpt-Wpt[s,c]).(Ypt[s,c]-Wpt[s,c])];

Vertices of the reflected  triangle B'.  The original triangle B has vertices (Subscript[B, 1]),(Subscript[B, 2]),(Subscript[B, 3]) in clockwise order; these correspond to the second vertex of the unrotated hexagon, and intersection points 3 and 2.  So the vertices of the reflected triangle in counterclockwise order are:

BpVerts[s_,c_]=Simplify[Map[refM.#&,{H1[[2]],pih[s,c][[3]],pih[s,c][[2]]}]]

drawBp[σ_]:=Module[{s,c,bpv,marg},
s=Sin[σ];
c=Cos[σ];
bpv=BpVerts[s,c];
marg=0.02;
Graphics[
{
Table[
{Line[{bpv[[i]],bpv[[Mod[i+1,3,1]]]}],
Text[(Subscript["B'", i]),bpv[[i]],{{-0.5,-1.2},{1.5,0},{-1.5,0}}[[i]]]}
,{i,1,3}],
{PointSize[.0075],
{Point[Xpt[s,c]],Text["X",Xpt[s,c],{1,-1.2}]},
{Point[Wpt[s,c]],Text["W",Wpt[s,c],{0,-1.2}]},
{Point[Ypt[s,c]],Text["Y",Ypt[s,c],{-1,1.2}]},
Line[{Wpt[s,c],Ypt[s,c]}],
Circle[Opt[s,c],1],
Circle[Npt[s,c],1]
}
}
,BaseStyle->18,ImageSize->{600,300},PlotRange->{{bpv[[2,1]]-marg,bpv[[3,1]]+marg},{Min[bpv[[2,2]],bpv[[3,2]]-marg],bpv[[1,2]]+marg}}
]]

drawBp[1.3 Degree]

Outwards-pointing normals for the sides of triangle B'

BpNorms[s_,c_]=Simplify[cst[Simplify[Table[RotationMatrix[π/2].(BpVerts[s,c][[i]]-BpVerts[s,c][[Mod[i+1,3,1]]]),{i,1,3}]],3]]

Values on sides of dot products with these normals

BpDP[s_,c_]=Simplify[cst[Simplify[Table[BpNorms[s,c][[i]].BpVerts[s,c][[i]],{i,1,3}]],3]]

Check whether a point lies inside the triangle B'

insideBp[s_,c_,x_]:=Apply[And,Table[BpNorms[s,c][[i]].x<BpDP[s,c][[i]],{i,1,3}]]

Polygon for the third reduction:  this joins up all the vertices of the cover, but needs two chord areas added to it, between (Subscript[O, tan]) and W and between (Subscript[N, tan]) and Y.

red3[s_,c_]=Module[{rh,ph},
rh=rotH[s,c];
ph=pih[s,c];
{Ypt[s,c],Ntan[s,c],H1[[2]],ph[[4]],ph[[5]],H1[[4]],ph[[8]],ph[[9]],H1[[6]],Otan[s,c],Wpt[s,c]}
];

Find the two chord lengths

chord1[s_,c_]=Sqrt[(Wpt[s,c]-Otan[s,c]).(Wpt[s,c]-Otan[s,c])];

chord2[s_,c_]=Sqrt[(Ypt[s,c]-Ntan[s,c]).(Ypt[s,c]-Ntan[s,c])];

Total area of cover is area of polygon plus areas between chords and circles

coverArea[s_,c_]=areaPoly[red3[s,c]]+chordArea[chord1[s,c]]+chordArea[chord2[s,c]];

Nominal rotation angle from paper:  1.3 degrees

σNom=(13/10) Degree

N[coverArea[Sin[σNom],Cos[σNom]],100]

Check that dot products are negative

N[{dp[s,c],dp2[s,c]}/.{s->Sin[σNom],c->Cos[σNom]},100]

Check that points W,X and Y are inside the triangle B'

Map[insideBp[Sin[σNom],Cos[σNom],#]&,
{Wpt[Sin[σNom],Cos[σNom]],Xpt[Sin[σNom],Cos[σNom]],Ypt[Sin[σNom],Cos[σNom]]}]

Plots around this value

Plot[coverArea[Sin[σ*Degree],Cos[σ*Degree]],{σ,1.2,1.4}]

Plot[dp[Sin[σ*Degree],Cos[σ*Degree]],{σ,1.2,1.4}]

Plot[dp2[Sin[σ*Degree],Cos[σ*Degree]],{σ,1.2,1.4}]

Find a value closer to the minimum possible

N[Sin[σNom]]

s0=N[s/.FindRoot[Evaluate[dp[s,(Sqrt[1-(s^2)])]],{s,22/1000},WorkingPrecision->2000],20]

ArcSin[s0]/Degree

We round up to get a safe value

σ0=(1294389444703601012/(10^18)) Degree;

N[coverArea[Sin[σ0],Cos[σ0]],100]

Check that dot products are negative

N[{dp[s,c],dp2[s,c]}/.{s->Sin[σ0],c->Cos[σ0]},100]

Check that points W,X and Y are inside the triangle B'

Map[insideBp[Sin[σ0],Cos[σ0],#]&,
{Wpt[Sin[σ0],Cos[σ0]],Xpt[Sin[σ0],Cos[σ0]],Ypt[Sin[σ0],Cos[σ0]]}]

See if we can find an explicit formula for the critical angle where the first dot product is zero.

dpsc0=Map[Collect[#,{α,β},Simplify]&,Factor[cst[Collect[Numerator[Together[Collect[Numerator[Together[dp1[s,c]]],{α,β},Simplify]/.{(β^2)->((β1)^2)}]],{α,β}],3],Extension->Automatic]]

dpsc1=Collect[Collect[Collect[Numerator[Together[(α^2)-((α1)^2)/.Solve[dpsc0==0,α][[1]]]],β,cstNumerator[#,3]&],β]/.{(β^2)-> ((β1)^2)},β,cstNumerator[#,3]&]

dpsc2=Collect[cst0[Collect[cst[Numerator[Together[(β^2)-((β1)^2)/.Solve[dpsc1==0,β][[1]]]],3],c]],s,Simplify]

dpsc3=Factor[Numerator[Together[(s^2)+(c^2)-1/.Solve[dpsc2==0,s][[1]]]],Extension->Automatic]

The best we can do is express the cosine of the critical angle as a root of a degree-16 polynomial.  With a change of variable, we can give the polynomial integer coordinates.

dpsc4=Factor[Numerator[Together[dpsc3[[-1]]/.{c->χ/(Sqrt[3])}]]]

((ArcCos[χ/(Sqrt[3])])/Degree)/.NSolve[dpsc4==0,Reals,WorkingPrecision->100]