(* ::Package:: *)

(* ::Subsubsection:: *)
(*Initialization Code*)


straightFoliationLines[{lineDensityHorizontal_:1,lineDensityVertical_:1},
{tanHorizontal_:0.0,tanVertical_:0.0},transform_:(#&),offset_:{0,0},style_:Red]:=
{If[lineDensityHorizontal!=0,Style[Table[Line[transform/@{{-100+offset[[1]],
k-100tanHorizontal+offset[[2]]},{100+offset[[1]],k+100tanHorizontal+offset[[2]]}}],
{k,-100.5,100.5,1/lineDensityHorizontal}],style],{}],
If[lineDensityVertical!=0,Style[Table[
Line[transform/@{{k-100tanVertical+offset[[1]],-100+offset[[2]]},
{k+100tanVertical+offset[[1]],100+offset[[2]]}}],{k,-100.5,100.5,1/lineDensityVertical}],
style],{}]}


arrayPlotFoliationLine[layer_, totalSteps_] := 
Module[{horizontalCoordinates, verticalCoordinates}, 
horizontalCoordinates = 
Catenate[ConstantArray[#, 2] & /@ Join[Range[totalSteps + 1 - layer, totalSteps, 2], 
Reverse@Range[totalSteps + layer, totalSteps + 1, -2]]];
verticalCoordinates = 
Join[#, Reverse[#]] &@Rest@Most@Catenate[
ConstantArray[#, 2] & /@ Range[totalSteps + 1, totalSteps - Ceiling[layer/2] + 1, -1]];
Transpose[{horizontalCoordinates, ReplacePart[verticalCoordinates, 
{1 -> verticalCoordinates[[1]] + 1, -1 -> verticalCoordinates[[-1]] + 1}]}]]


foliationLine[coordinates_,style_:Red][extraPoints_,chosenVertices_]:=
Module[{chosenCoordinates=
Join[coordinates/@chosenVertices,extraPoints],nearest1,nearest2,coordinateBounds=
CoordinateBounds[coordinates]},{nearest1,nearest2}=
Nearest/@{chosenCoordinates,Complement[Values[coordinates],chosenCoordinates]};{
ContourPlot[EuclideanDistance[nearest1[{x,y}][[1]],{x,y}]-
EuclideanDistance[nearest2[{x,y}][[1]],{x,y}]==0,{x,coordinateBounds[[1,1]],
coordinateBounds[[1,2]]},{y,coordinateBounds[[2,1]],coordinateBounds[[2,2]]},
ContourStyle->style],Table[With[{x=x},Module[{func},func[y_?NumericQ]:=
EuclideanDistance[nearest1[{x,y}][[1]],{x,y}]-EuclideanDistance[nearest2[{x,y}][[1]],{x,y}];{x,y/.FindRoot[
func[y]==0,{y,coordinateBounds[[2,1]],coordinateBounds[[2,2]]}]}]],
{x,coordinateBounds[[1,1]],coordinateBounds[[1,2]],1}];{}}]


drawFoliation[graph_,vertexLists_,style___]:=
Quiet[Module[{vertexCoordinates=
Association[Thread[VertexList[graph]->(VertexCoordinates/.AbsoluteOptions[graph,
VertexCoordinates][[1]])]],lines},Show[graph,FoldPairList[
foliationLine[vertexCoordinates,style],{},vertexLists]]],NearestFunction::neard]


(* ::Subsubsection:: *)
(*Final Code*)


(* ::Input:: *)
(*(*https://www.wolframcloud.com/obj/wolframphysics/TechPaper-Programs/Section-05/FoliationLines.wl*)*)
(*CloudGet["https://wolfr.am/KXgcRNRJ"]*)