Package ModularGroup

Introduction

This is a demonstration for the Mathematica package ModularGroup.

About

The package ModularGroup has been developed in the course of the diploma thesis Computer Algebra and Analysis: Complex Variables Visualized, carried out at the Research Institute for Symbolic Computation (RISC) of the Johannes Kepler University Linz by Thomas Ponweiser under supervision of Univ.-Prof. Dr. Peter Paule. The thesis is available online at http://www.risc.jku.at/publications/.

The package has been developed and tested with Wolfram Mathematica 9 (Student Edition).

The current version (v0.1) contains prototypical implementations of the algorithms described in the above mentioned thesis. Check https://github.com/ponweist/ModularGroup for newer versions of the package with improved and extended features.

License

The MIT License (MIT)

Copyright (c) 2014 Thomas Ponweiser

Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the “Software”), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED “AS IS”, WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

Installation

Prior to viewing this Mathematica Notebook (evaluating its cells), the package ModularGroup needs to be loaded.

It is recommended to install the package ModularGroup such that Mathematica can load it via:

Needs["ModularGroup`"];

Package installation can be done in (at least) two ways:

Variant A (works in particular for Mathematica 9)

Select File → Install...

In the install dialog, for “Type of item to install” select “Package” and for “Source” select “From file...”

In the file selection dialog navigate to the file “ModularGroup.m” and open it.

Back in the install dialog, click OK.

Variant B

Evaluate the following cell:

Print[$UserBaseDirectory]

Copy the file “ModularGroup.m” into the subdirectory “Applications”, which is located in the directory listed above.

Alternatively, you can of course open “ModularGroup.m” with Mathematica and run the package manually (without installing) prior to viewing this Notebook.

Feature overview and some demonstrations

Moebius Transformations

?MoebiusTransformation

phi=New[MoebiusTransformation,{{I,1},{1,I}}];

Note that phi is the “modified Cayley transformation”, which is also predifined:

?mtModCayley

Evaluation for complex numbers:

{phi["z"], phi[I], phi[-I], phi[ComplexInfinity]}

Evaluation for homogeneous representation of complex numbers:

?Inhom

hom={phi@{"z",1},phi@{I,1},phi@{-I,1},phi@{1,0}}
Inhom/@hom

Composition of transformations:

psi=New[MoebiusTransformation,{{1,-1},{0,1}}];
composition=phi@psi;
{composition[1+I],phi[psi[1+I]]}

composition["z"]==phi[psi["z"]]//Simplify

?Mat

Mat[phi]//MatrixForm

?Inv

{phi["z"],Inv[phi]["z"]}

Inv[phi]@phi@"z" ==phi@Inv[phi]@"z"=="z"//Simplify

Modular Transformations

?ModularTransformation

ModularTransformations are MoebiusTransformations:

phi=New[ModularTransformation,{{13,5},{-8,-3}}];
InstanceQ[MoebiusTransformation][phi]

Some MoebiusTransformations are ModularTransformations:

psi=New[MoebiusTransformation,1/Sqrt[2]{{1,I},{-I,-1}}];
{Mat[psi]//MatrixForm,InstanceQ[ModularTransformation][psi]}

psiSqr=psi@psi;
{Mat[psiSqr]//MatrixForm, InstanceQ[ModularTransformation][psiSqr]}

There are some predefined ModularTransformations:

?mtId
?mtU
?mtT
?mtR

Finding a group word representation for a ModularTransformation:

?TRWord

TRWord[phi]

TRList gives the factors of the above product as list of ModularTransformations:

?TRList

Module[{factors},
factors = TRList[phi];Row@{TRWord[phi]," = ",Sequence@@Riffle[MatrixForm@Mat@#&/@factors,"·"], " = ", MatrixForm@Mat@phi}
]

If only the rightmost or leftmost factor in the product representation is needed, TRLeft and TRRight may be used:

?TRLeft
?TRRight

TRLeft[phi]===Inv@mtR
TRRight[phi]===Inv@mtR

In general, the word representation in terms of the generators T and U is shorter than that in terms of T and R, but it is not unique any more:

?TUWord

Row@{TRWord[phi], " = ", TUWord[phi]}

?QuotientFunction

?CentralQuotient
?CQuotient

Grid[Table[{f,TUWord[phi,QuotientFunction->f]},{f,{Quotient,CentralQuotient,CQuotient}}],Alignment->Left]

Sometimes the list of exponents of the U’s is more useful:

?TUExponents

Grid[Table[{f,TUExponents[phi,QuotientFunction->f]},{f,{Quotient,CentralQuotient,CQuotient}}],Alignment->Left]

TUEval can be used for evaluating group words:

?TUEval

Module[{factors,word,phi},
phi=New[ModularTransformation,{{13,5},{-8,-3}}];
Manipulate[
word=TUWord[phi,QuotientFunction->f];
factors=TUEval[phi,Mat[mtT],Mat[mtU],List,MatrixPower,QuotientFunction->f];
Row@{
word," = ", 
Sequence@@Riffle[MatrixForm[#]&/@factors,"·"], " = ",
MatrixForm[Dot@@factors]
},
{{f,Quotient,"Integer division function:"},{Quotient,CentralQuotient,CQuotient}}
]
]

For creating a List of ModularTransformations, ModularGroupList may be used:

?ModularGroupList

Creating a list of all ModularTransformations, where the absolute values of all U exponents add up to at most 3:

Module[{p,transformations},
Manipulate[
p=Function[t,Total@Abs@TUExponents[t,QuotientFunction->f]≤ 3];
transformations=ModularGroupList[p];
TUWord[#,QuotientFunction->f]&/@transformations,
{{f,Quotient},{Quotient,CentralQuotient,CQuotient}}
]
]

For searching certain ModularTransformations, it is more practical to iterate over an infinite sequence rather than working with finite lists:

?ModularGroupIterator

Find the “first” six ModularTransformations, whose coefficients add up to 10:

Module[{it,t,n,transformations},
it=ModularGroupIterator[PSL2ZGrading@Mat[#1]≤PSL2ZGrading@Mat[#2]&];
n=0;
transformations = Reap[
While[n<6,
t=GetNext[it];
If[Total[Mat@t,Infinity]==10,
Sow[t];
n++;
];
]
][[2,1]];
Column@Table[Row@{TUWord[t],": ", MatrixForm@Mat@t},{t,transformations}]
]

Generalized Circles

?GeneralizedDisk

There are some convenience functions for creating GeneralizedDisks:

?CompactDisk
?NoncompactDisk
?Halfplane

CompactDisks are defined by their radius and center, Noncompact disks by their “coradius” and “cocenter”:

?GDiskRadius
?GDiskCenter
?GDiskCoradius
?GDiskCocenter

Predefined GeneralizedDisks are:

?gdUnitDisk
?gdUpperHalfplane
?gdLowerHalfplane
?gdRightHalfplane
?gdLeftHalfplane

{GDiskCenter[gdUnitDisk],GDiskRadius[gdUnitDisk]}

Testing if a GeneralizedDisk contains a given point:

?InteriorQ

InteriorQ[gdUnitDisk,#]&/@{0,1,2,ComplexInfinity}

InteriorQ[gdUnitDisk,#]&/@{{0,1},{1,1},{2,1},{1,0}}

For drawing GeneralizedDisks, RegionPlot in conjunction with InteriorQ may be used:

Module[{disk},
disk["Unit disk"]=gdUnitDisk;
disk["Upper half"]=gdUpperHalfplane;
disk["Lower half"]=gdLowerHalfplane;
disk["Right half"]=gdRightHalfplane;
disk["Left half"]=gdLeftHalfplane;
Manipulate[
RegionPlot[InteriorQ[disk[Selection],x+I y],{x,-2,2},{y,-2,2},
Frame->False,Axes->True],
{{Selection, "Unit disk"},{"Unit disk","Upper half","Lower half","Right half","Left half"}}
]
]

A faster alternative for drawing GeneralizedDisks is provided by GDisk and GCircle:

?GDisk
?GDiskNPoints
?GCircle

The image of the unit disk under certain powers of the loxodromic Moebius transformation z ↦((16z + 16 + 32 i)/((1 + i)z + 16)):

Module[{transformations,m},
m={16{1,1+2I},{1+2I,16}};
transformations=Table[MatrixPower[m/Sqrt@Det@m,p],{p,-10,10,0.5}];
Graphics[{Lighter@Blue,Opacity[0.5],GDisk[gdUnitDisk,transformations]},
Axes->True,PlotRange->12]
 ]

Generalized circles can be transformed by Moebius transformations, e.g. the upper halfplane is transformed to the unit disk by z ↦ ((i z + 1)/(z + i)):

phi=New[MoebiusTransformation,{{I,1},{1,I}}];
disk=phi[gdUpperHalfplane];
InteriorQ[disk,"z"]//FullSimplify

Continuous animation of this transformation and its powers up to 4:

Module[{phi,disk},
Manipulate[
phi=New[MoebiusTransformation,MatrixPower[{{I,1},{1,I}},t]];
disk=phi[gdUpperHalfplane];
Column@{
Superscript["ϕ",t],
Graphics[{
{Lighter@Blue,GDisk[disk]},
{Darker@Blue,GCircle[disk]}
},Axes->True,PlotRange->2,ImageSize->Medium
]
},
{t,0,4,Animator,DefaultDuration->20,AnimationRunning->False}]
]

The above transformation can be thought of a rotation of the Riemann sphere:

Module[{phi,disk,cone,halfSphere,regionPlot},
halfSphere=ParametricPlot3D[{Cos[phi]Cos[theta],Sin[theta],-Sin[phi]Cos[theta]},{phi,-Pi,Pi},{theta,0,Pi/2},Mesh->False,PlotStyle->{Darker@Green}][[1]];
cone[disk_]:=Module[{c,r,color},
If[InstanceQ[Halfplane][disk],
{},
If[InstanceQ[CompactDisk][disk],
c=GDiskCenter[disk];r=GDiskRadius[disk];color=Yellow,
c=GDiskCocenter[disk];r=GDiskCoradius[disk];color=Black;
];
{color,
GeometricTransformation[Cone[{{0,0,-2},{0,0,0}},2 r],{{{1,0,-Re[c]},{0,1,-Im[c]},{0,0,1}},{0,0,1}}]
}
]
];
Manipulate[
phi=New[MoebiusTransformation,MatrixPower[{{I,1},{1,I}},Mod[t,4,-2]]];
disk=phi[gdUpperHalfplane];
regionPlot=RegionPlot[InteriorQ[disk,x+I y],{x,-2,2},{y,-2,2},PlotStyle->Lighter@Blue,Axes->True,AxesLabel->{"Re[z]","Im[z]"},Frame->False];
Column@{
Superscript["ϕ",t],
GraphicsRow[{
regionPlot,
Graphics3D[{
Rotate[halfSphere,Pi t/2,{-1,0,0}],
(regionPlot[[1]]/.{x_Real,y_Real}->{x,y,0}),
{PointSize->Large,Red,Point[{0,0,1}]},If[showProjection,{Opacity[0.25],EdgeForm[],cone[disk]},{}]
},PlotRange->2,Axes->True,AxesLabel->{"Re[z]","Im[z]",""},Lighting->"Neutral",ViewPoint->{5,-2,2}]
},ImageSize->Large]
},
{t,0,4,Animator,DefaultDuration->20,AnimationRunning->False},{{showProjection,True,"Show projection"},{True,False}}]
]

The fundamental domain and the tiling of the upper half-plane

The tiling of the upper complex halfplane associated to the modular group can easily plotted:

?ModularTiling

To produce a list of modular transformations used by the function ModularTiling, we define the “grading” of a modular Transformation to be the sum of its coefficients squared.

?PSL2ZGrading

Now, we list all transformations with a grading not greater than 100.

mtList=ModularGroupList[PSL2ZGrading@Mat@#≤100&];
Length@mtList

Graphics[ModularTiling[mtList,mtId,LabelMode->TUWord,Magnification->5],PlotRange->{{-#,#}&@GoldenRatio,{0,2}},Axes->True,AxesStyle->Gray]

Here is a more fancy demonstration of ModularTiling and its options:

Module[{mtListPreview},
mtListPreview=Select[mtList,PSL2ZGrading@Mat@#≤40&];
Manipulate[
Module[{phi0,phi,opts,mts,viewCenter,plotRange,labelFu},
mts=ControlActive[mtListPreview,mtList];
phi0=MatrixPower[{{I,1},{1,I}},t];
phi=phi0.{{1/Cos[s],Tan[s]},{0,1}};
viewCenter=Inhom@(phi0.{I,1});
plotRange={Re[viewCenter]+{-#,#}&@GoldenRatio,Im[viewCenter]+{-1,1}};
labelFu=Which[
label==="Function",Function[m,Inhom[If[MatrixQ[m],m,Mat@m].{"z",1}]],
label==="Matrix",Function[m,MatrixForm@If[MatrixQ[m],m,Mat@m]],
True,label
];
opts={};
If[!MemberQ[show,"Tessellation"],
AppendTo[opts,TilingMode->Off];
];
If[MemberQ[show,"Incircles"],
opts=Join[opts,{IncircleMode->GDisk,IncircleStyle->Directive[EdgeForm@Lighter@Blue,Lighter@Lighter@Blue]}];
];
If[MemberQ[show,"Ford circles"],
opts=Join[opts,{FordDiskMode->GDisk,FordDiskStyle->Directive[EdgeForm@Lighter@Brown,Lighter@Lighter@Brown]}];
];
If[!(label===Off),
opts=Join[opts,{LabelMode->labelFu,Magnification->labelSize}]
];
Graphics[{
If[MemberQ[show,"Fundamental region"],{LightGray,GDisk[gdUnitDisk,phi.#&/@{Mat@mtId,Mat@mtR,{{0,1},{-1,1}}}]}],
ModularTiling[mts,phi,Sequence@@opts]
},
PlotRange->plotRange,Axes->True,AxesStyle->Gray,ImageSize->imageSize
]
],
{{t,0,"Plane ↔ Disk"},0,1},
{{s,0,"Center: I ↔ ρ"},0,Pi/6},
{{show,{"Fundamental region", "Tessellation"},"Show"},{"Fundamental region", "Tessellation","Incircles","Ford circles"},TogglerBar},
{{label,Off,"Label"},{Off,TUWord,TRWord,"Matrix","Function"}},
{{labelSize,5,"Label size"},1,10},
{{imageSize,400,"Image Size"},200,1024},
ControlPlacement->Left
]
]

The action modular transformations on the modular tiling

We can animate the action of modular transformations on the modular tiling. 
Selection of a modular transformation A can be done by clicking on its corresponding tile (that is the image of the fundamental domain under A).
For a given point of the upper half-plane, the corresponding transformation A and its inverse can be found with the following routines:

?AssociatedMap
?InvAssociatedMap

We use the modified Cayley transform to map the tiling of the upper half-plane to the unit disk. 
In order to animate the action of the selected transformation A, we use Mathematica’s built-in MatrixPower.
For a parameter t ∈ [0,1], the “intermediate image” of the tiling is induced by the Möbius transformation MatrixPower[A,t].
Note the little trick that for drawing the tiling, we use a paramter t’ which cycles through the interval [-0.5,0.5] instead.
This ensures that we see the same tiling lines for t=0 and t=1.

Module[{mtListPreview,modCayley},
mtListPreview=Select[mtList,PSL2ZGrading@Mat@#≤40&];
modCayley={{I,1},{1,I}};
Manipulate[
Module[{mts,z,z0,fix,t,phipt,phitess,labelFu},
mts=ControlActive[mtListPreview,mtList];
z=Inhom[PSL2CInv[modCayley].{Complex@@P,1}];
t=If[z===ComplexInfinity,IdentityMatrix@2,AssociatedMap[z]];
z0=Inhom[PSL2ZInv[t].{z,1}];
fix=If[t==IdentityMatrix@2,Null,First@Select[Eigenvectors[t],Im[#[[1]]Conjugate@#[[2]]]≥0&]];
phipt=modCayley.MatrixPower[t,p];
phitess=modCayley.MatrixPower[t,Mod[p,1,-1/2]]; 
labelFu=Which[
label==="Function",Function[m,Inhom[If[MatrixQ[m],m,Mat@m].{"z",1}]],
label==="Matrix",Function[m,MatrixForm@If[MatrixQ[m],m,Mat@m]],
True,label
];

Column@{
Graphics[{
{
ModularTiling[mts,phitess],
ModularTiling[mts,phipt,
TilingMode->Off,BorderMode->Off,ExteriorMode->Off,
LabelMode->labelFu,LabelThreshold->1/labelThreshold,Magnification->labelSize]
},
{Red,PointSize@Medium,Point[{Re@#,Im@#}&@Inhom[phipt.{z0,1}]]},
If[t≠IdentityMatrix@2,{Blue,PointSize@Medium,Point[{Re@#,Im@#}&@Inhom[phipt.fix]]}]
},
PlotRange->1.05,(*Axes->True,AxesStyle->Gray,*)ImageSize->imageSize
],
Grid[{
{"Selected point:",Round[z,0.001]},
{"Associated transformation:", labelFu@t},
{"Fixed point in ℋ:",If[fix===Null,"(all points)",Simplify@Inhom@fix/.ComplexInfinity->∞]},
{"Equivalent point in ℱ:",Round[z0,0.001]}
},Alignment->{Left,Left}]
}
],
{{P,{0,0}},Locator},
{{p,0,"Animate Transformation"},0,1},
{{label,TUWord,"Label"},{TUWord,TRWord,"Matrix","Function"}},
{{labelThreshold,8,"Less labels ↔ More labels"},8,32},
{{labelSize,8,"Label size"},1,16},
{{imageSize,400,"Image Size"},200,1024},

ControlPlacement->Left
]
]

Hyperbolic geometry: Canonical fundamental regions

The following function finds a hyperbolic rigid motion A of the upper half-plane (that is a Möbius transformation with real coefficients)
which maps two points z’ and w’ on the unit disk with Im[z’]=Im[w’]
to the given points z and w on the upper half-plane.
* The hyperbolic line joing z and w is given by the image of the unit disk under A.
* The hyperbolic perpendicular bisector of z and w is given by the image of the imaginary axis under A.
Note: I’m sure that this can be done smarter and faster...

FindHypMotion[z_,w_]:=Module[{mat,γ,λ,μ,sol},
mat={{z,w},{1,1}}.{{λ,0},{0,μ}}.{{γ,1},{-1,γ}};
sol=FindInstance[{Flatten@mat∈Reals,Abs@γ==1},{γ,λ,μ}][[1]];
mat/.sol
];

ModSubgroupQ returns a membership predicate for the certain subgroups of the modular group:

ClearAll[ModSubgroupQ]

ModSubgroupQ["Γ",n_]:=True&;
ModSubgroupQ["Γ(n)",n_]:=(Abs@Mod[Mat@#,n,-n/2]==IdentityMatrix[2])&;
ModSubgroupQ["(Subscript[Γ, 0])(n)",n_]:=MatchQ[Mod[Mat@#,n],{{_,_},{0,_}}]&;
ModSubgroupQ["(Γ^0)(n)",n_]:=MatchQ[Mod[Mat@#,n],{{_,0},{_,_}}]&;
ModSubgroupQ["(Subscript[Γ, 0]^0)(n)",n_]:=MatchQ[Mod[Mat@#,n],{{_,0},{0,_}}]&;

A canonical fundamental domain for a subgroup G of the modular group can be found the following way:
* Fix a point z on the upper halfplane.
* For each g ∈ G\{1}, determine the hyperbolic perpendicular bisector (Subscript[L, g]) of the hyperbolic line through z and gz.
* (Subscript[L, g]) splits the hyperbolic plane (that is the upper half-plane) into two hyperbolic half-planes (Subscript[A, g]) and (Subscript[B, g]), where z ∈ (Subscript[A, g]) and gz ∈ (Subscript[B, g]).
* The points w ∈ (Subscript[A, g]) are characterized by the property d(w,z) < d(w,gz), where d is the hyperbolic distance metric.
* The canonical fundamental domain (with respect to the point z) is the intersection of all these half-planes (Subscript[A, g]), g ϵ G.

A canonical fundamental domain F is therefore characterized by the property d(w,z) < d(w,gz) for all w ∈ F and g ∈ G.

WARNING: The following demonstration is in a preliminary state and neither optimized for speed nor well-tested.

Needs["ModularGroup`"];
Module[{modgrp,GetHypMotion,oldz},
modgrp=ModularGroupList[PSL2ZGrading@Mat@#<300&];
modgrp=Sort[modgrp,(PSL2ZGrading@Mat@#1<PSL2ZGrading@Mat@#2)&];
Manipulate[
Module[{subgrp,z,ws,m,ts},
subgrp=Select[modgrp,ModSubgroupQ[subgroup,n]];

z=Complex@@Rationalize@P;
If[!TrueQ[oldz==z],
oldz=z;
ClearAll[GetHypMotion];
GetHypMotion[a_,b_]:=GetHypMotion[a,b]=FindHypMotion[a,b];
];

ws=DeleteCases[#[z]&/@subgrp,z];
If[Length@ws>max,ws=Take[ws,max]];
ts=GetHypMotion[z,#]&/@ws;
Graphics[{
{PointSize[Medium],Point[{Re@#,Im@#}&/@ws]},
{Opacity[0.1],Black,GDisk[gdRightHalfplane,ts]},
{Gray,Dashed,GCircle[gdUnitDisk,ts]},
{GCircle[gdRightHalfplane,ts]},
{White,GDisk[gdLowerHalfplane]}
},
Axes->True,PlotRange->{{-2,2},{0,3}},PlotLabel->z
]
],
{{P,{0,3/2}},Locator,ContinuousAction->False},
{{subgroup,"Γ","Subgroup"},{"Γ","Γ(n)","(Subscript[Γ, 0])(n)","(Γ^0)(n)","(Subscript[Γ, 0]^0)(n)"}},
{n,2,5,1,Appearance->"Open"},
{{max,9,"max. Transformations"},1,20,1,Appearance->"Open"}
]
]

The modular Group and continued fractions

MyTicks[{min_,max_},special_List,n_Integer]:=Module[{ticks,spInRange},
spInRange=Select[special,IntervalMemberQ[Interval[{min,max}],#]&];
Join[
N@Select[FindDivisions[{min,max},n],And@@(Thread[Abs[#-spInRange]≥(max-min)/n])&],
spInRange
]
];
MyContinuedFraction[x_,n_]:=
Module[{i,c,t=x,cs={}},
For[i=0;c=Quotient[t,1,-1/2],
i<n&&t-c≠0,
t=1/(t-c);c=Quotient[t,1,-1/2];i++,
cs={cs,c}];
If[i<n,cs={cs,c}];
Flatten@cs
];
ModularConvergents=Compile[{{quotients,_Integer,1}},
Module[{sgn,i,m},
sgn=-1;
i=1;
m={{0,-1},{1,0}};
Table[sgn*=-1;m=m.{{0,-1},{1,sgn quotients[[i]]}},{i,1,Length[quotients]}]
]
];

ModularAuxConvergents[quotients_List]:=Module[{exponents,list,factors},
exponents=MapIndexed[(((-1))^(First@#2+1))#1&,quotients];
list=Flatten@Riffle[
ConstantArray[Sign@#,Abs@#]&/@exponents,
0,2
];
If[Length[list>100],list=Take[list,100]];
factors=list/.{
-1->{{1,-1},{0,1}},
0->{{0,-1},{1,0}},
1->{{1,1},{0,1}}
};
FoldList[Dot,IdentityMatrix@2,factors]
];

The following is a demonstration for nearest integer continued fraction representations of arbitrary numbers (any Mathematica expression may be used as input). Be careful with too high zoom factors - errors may occur due to machine-precision calculations.

Module[{mtList,mtPreviewList,quotients,convergents,controlPoints,f,m,plotRange,zf,ts,viewCenter},
mtList=ModularGroupList[PSL2ZGrading@Mat@#<100&];
mtPreviewList=Select[mtList,PSL2ZGrading@Mat@#<25&];

Manipulate[
quotients=MyContinuedFraction[number,6];
convergents=ModularConvergents[quotients];
controlPoints={Re@#,Im@#}&/@GDiskMatCenter[GDiskMatMap[Mat@gdIncircle0,convergents]];
f=BSplineFunction[controlPoints,SplineDegree->2];

Manipulate[
ts=ControlActive[mtPreviewList,mtList];
m=convergents[[1+Floor[t (Length[convergents]-1)]]];
viewCenter=f[t];
zf= viewCenter[[2]];
plotRange={viewCenter[[1]]+GoldenRatio{-zf,zf},{0,2zf}};
Graphics[{
ModularTiling[ts,m,
FordDiskMode->GDisk,FordDiskStyle->LightGray,FordDiskThreshold->0,
IncircleMode->GDisk,IncircleStyle->Directive[White,EdgeForm[Gray]],IncircleThreshold->0,
TilingStyle->Gray, TilingThreshold->0,
ExteriorMode->Off,
LabelMode->((((-1))^(Length@#+1))Last@#&@TUExponents[m.Mat@#,QuotientFunction->CentralQuotient]&),
Magnification->8/zf,LabelThreshold->(2^(-4))zf,
FordLabelMode->(Inhom[m.(Mat@#).{1,0}]&),FordLabelThreshold->(2^(-3))zf
],
{Yellow,Opacity[0.5],GDisk[gdIncircle0,convergents]},
{White,Rectangle[Scaled@{-1,-1},Scaled@{2,0}]}
},
PlotRange->plotRange,Axes->True,ImageSize->Large,Ticks->{MyTicks[{#1,#2},{number},12]&,Automatic},PlotLabel->Row@{number,": ",quotients}
],
{{t,0,"Zoom"},0,1,ControlType->Animator,DefaultDuration->30,AnimationRunning->False}
],
{{number,GoldenRatio,"Number"},ControlType->InputField}]
]

Thesis figures

This section lists the Mathematica commands used to to create the figures of the diploma thesis Computer Algebra and Analysis: Complex Variables Visualized.

Circle inversion

With[{t=Exp[-1],phi=0.3},
Module[{c={0,0},p={t,0},q={1/t,0},n1,n2,m},
n1={t,Sqrt[1-((t )^2)]};
n2={{1,0},{0,-1}}.n1;
m=RotationMatrix[phi];
{c,p,q,n1,n2}={c,p,q,n1,n2}.Transpose@m;
Graphics[
{
{Darker@Green,Circle[c,1]},
{Darker@Green,Arrow@{c,c+m.{-0.6,0.8}},Text["r",c+m.{-0.6,0.8}/2,{0,1}]},
{Red,Line@{c,Scaled[m.{1,0},q]}},
{Dashed,Line@{Scaled[m.{0,-1},p],Scaled[m.{0,1},p]}},
{Dashed,Circle[(c+q)/2,Norm[q-c]/2]},
{PointSize@Medium,
Function[{pt,label},{Point[pt],Text[label,pt,{0,-1.5}]}]@@@{{c,"C"},{p,"P"},{q,"Q"}},
Point/@{n1,n2}
}
},PlotRange->{{-1,3},{-2,2}+0.25},ImageSize->200
]
]
]

Modified Cayley Transform, Stereographic Projection

Needs["ModularGroup`"];
Module[{phi,disk,cone,rays,halfSphere,regionPlot,plots,figure,plots2,figure2},
halfSphere=ParametricPlot3D[{Cos[phi]Cos[theta],Sin[theta],-Sin[phi]Cos[theta]},{phi,-Pi,Pi},{theta,0,Pi/2},PlotPoints->6,Mesh->False,PlotStyle->{Darker@Green}][[1]];
rays[{p1_,p2_},r_]:=Table[Line@{p1+{r Cos[alpha],r Sin[alpha],0},p2},{alpha,0,11Pi/6,Pi/6}];
cone[disk_]:=Module[{c,r,color},
If[InstanceQ[Halfplane][disk],
{},
If[InstanceQ[CompactDisk][disk],
c=GDiskCenter[disk];r=GDiskRadius[disk];color=Gray,
c=GDiskCocenter[disk];r=GDiskCoradius[disk];color=Black;
];
{color,rays[{{2Re[c],2Im[c],-1},{0,0,1}},2r]}
]
];
plots=Table[
phi=New[MoebiusTransformation,MatrixPower[{{I,1},{1,I}},Mod[t,4,-2]]];
disk=phi[gdUpperHalfplane];
regionPlot=RegionPlot[InteriorQ[disk,x+I y],{x,-2,2},{y,-2,2},PlotPoints->6,PlotStyle->{Opacity[0.8],Lighter@Blue},Axes->True,AxesLabel->{"Re","Im"},Frame->False];
{
Graphics3D[{
{Opacity[0.8],Rotate[halfSphere,Pi t/2,{-1,0,0}]},
(regionPlot[[1]]/.{x_Real,y_Real}->{x,y,0}),
(*{PointSize->Large,Red,Point[{0,0,1}]},*)
{Dashed,cone[disk]}
},PlotRange->2,Axes->True,AxesLabel->{"Re","Im",""},Ticks->{Range[-2,2],Range[-2,2],Range[-2,2]},Lighting->"Neutral",ViewPoint->{5,-1.8,2}],
regionPlot
},
{t,0,3,1/2}
];
plots2=Table[
phi=New[MoebiusTransformation,MatrixPower[{{I,1},{1,I}},Mod[t,4,-2]]];
disk=phi[gdUpperHalfplane];
regionPlot=RegionPlot[InteriorQ[disk,x+I y],{x,-2,2},{y,-2,2},PlotPoints->6,PlotStyle->{Opacity[0.8],Lighter@Blue},Axes->True,AxesLabel->{"Re","Im"},Frame->False];
{
Graphics3D[{
{Opacity[0.8],Rotate[halfSphere,Pi t/2,{-1,0,0}]},
(regionPlot[[1]]/.{x_Real,y_Real}->{x,y,0}),
(*{PointSize->Large,Red,Point[{0,0,1}]},*)
{Dashed,cone[disk]}
},PlotRange->2,Axes->True,AxesLabel->{"Re","Im",""},Ticks->{Range[-2,2],Range[-2,2],Range[-2,2]},Lighting->"Neutral",ViewPoint->{5,-1.8,2}],
regionPlot
},
{t,0,1,1/3.}
];
figure=GraphicsGrid[Table[{plots[[i,1]],plots[[i,2]],plots[[i+2,1]],plots[[i+2,2]]},{i,1,5}],
Spacings->{Scaled[.25],Scaled[.5]},Dividers->{{3->Black},{}},
ImageSize->500
];
figure2=GraphicsGrid[plots2,Spacings->{Scaled[.5],Scaled[.25]},ImageSize->400];
Column@{figure,figure2}
]

Minimal region

Graphics[{
{Lighter@Lighter@Blue,EdgeForm@Black,Rectangle[Scaled[{0,-1},{-1/2,0}],Scaled[{0,1},{1/2,0}]]},
{White,EdgeForm[Black],Disk[{#,0},1]&/@{-1,1}},
{Transparent,EdgeForm[Directive[Dashed,Black]],Rectangle[Scaled[{0,-1},{-1/2,0}],Scaled[{0,1},{1/2,0}]]},
{PointSize@Medium,Function[{phi,label},Module[{coords={Cos@phi,Sin@phi}},{Point[coords],Text[label,coords,-2coords]}]]@@@{{Pi/3,"Tρ"},{2Pi/3,"ρ"},{4Pi/3,"-Tρ"},{5Pi/3,"-ρ"}}
},
{Dashed,Circle[]}
},Axes->True,PlotRange->{{-2.2,2.2},{-#,#}&@GoldenRatio},ImageSize->300]

Fundamental region

Module[{gdA, gdB},
gdA=New[Halfplane,-1/2,-1];
gdB =New[Halfplane, 1/2,1];
Graphics[{
{Opacity[0.2],{Red,GDisk[gdA]},{Green, GDisk[gdB]},{Blue,GDisk[gdUnitDisk]}},
{GCircle[gdA],GCircle[gdB],GCircle[gdUnitDisk]},
{Dashed,Black, GCircle[gdIncircle0]},
{White,GDisk[gdLowerHalfplane]},
{
PointSize[Medium],Function[{phi,psi,label},Module[{coords={Cos@phi,Sin@phi}},{Point[coords],Text[label,coords,-2{Cos@psi,Sin@psi}]}]]@@@{{Pi/3,Pi/4,Tρ},{Pi/2,0,""},{2Pi/3,3Pi/4,ρ}}
},
Text[Magnify[a],{-0.5,1.75},{2,0}],
Text[Magnify[b],{0.5,1.75},{-2,0}],
Text[Magnify[c],{-0.25,0.85}],
Text[Magnify[d],{0.25,0.85}],
Text[i,{0,1},-{4,1}],
Text[Magnify[𝔸,2],{-1.25,1.5}],
Text[Magnify[𝔹,2],{1.25,1.5}],
Text[Magnify[𝔻,2],{0,0.25}],
Text[Magnify[ℱ,2],{0,2.25}],
Text[Magnify[ℐ,2],{0,1.5}]
},PlotRange->{{-#,#}&@1.25GoldenRatio,{0,2.5}},Axes->True,
ImageSize->300
]
]

Modified Cayley Transform and Modular Tiling

Needs["ModularGroup`"];
Module[{figure1,criteria,gds,minSize=(2^(-9)),maxDistance=3,mtList,phi,viewCenter,plotRange},
criteria=Function[{mt},
gds=GDiskMatMap[Mat@gdIncircle0,{Mat@mt,{{I,1},{1,I}}.Mat@mt}];
(GDiskMatRadius@gds[[1]]≥minSize&&Abs@GDiskMatCenter@gds[[1]]<maxDistance)||GDiskMatRadius@gds[[2]]≥minSize
];
mtList=ModularGroupList[criteria];
SetOptions[TUExponents,QuotientFunction->CQuotient];
plots1=Table[
phi=MatrixPower[{{I,1},{1,I}},t];
viewCenter=Inhom[phi.{I,1}];
plotRange={{{-GoldenRatio,GoldenRatio},{0,2}},1}[[t+1]];
Graphics[ModularTiling[mtList,phi,
IncircleMode->GDisk, IncircleStyle->Directive[EdgeForm[Gray],LightGray],IncircleThreshold->(2^(-6)),
TilingThreshold->(2^(-6)),
LabelMode->TUWord,LabelStyle->Black,LabelThreshold->({1/12,1/16}[[t+1]]),Magnification->{5,8}[[t+1]]
],
PlotRange->plotRange,Axes->True,AxesStyle->Gray,ImageSize->{400,Automatic}
],
{t,{0,1}}
];
plots2=Table[
phi=MatrixPower[{{I,1},{1,I}},t];
viewCenter=Inhom[phi.{I,1}];
plotRange={Re@viewCenter+{-GoldenRatio,+GoldenRatio},Im@viewCenter+{-1,+1}};
Graphics[ModularTiling[mtList,phi,
IncircleMode->GDisk, IncircleStyle->Directive[EdgeForm[Gray],LightGray],IncircleThreshold->(2^(-4)),
TilingThreshold->(2^(-4)),
LabelMode->TUWord,LabelStyle->Black,LabelThreshold->1/7,Magnification->3
],
PlotRange->plotRange,Axes->True,AxesStyle->Gray,ImageSize->{200,Automatic}
],
{t,0,1,1/3.}
];
figure1=Module[{h, w,w2,h2, m=0.2,u,l},
h=m+GoldenRatio;
h2=GoldenRatio/4;
w2=GoldenRatio h2;
w=1+m+w2;
u=3h2+m;
l=h2/2;
Graphics[{
Inset[plots1[[1]],{1/2,1+m},{0,0},1],
Inset[plots1[[2]],{1/2,1/2},{0,0},1],
Inset[plots2[[#]],{1+m+w2/2,u(4-#)/3+l(#-1)/3},{0,0},1/GoldenRatio]&/@Range[4],
Line@{{1+m/2,0},{1+m/2,h}}
},PlotRange->{{0,w},{0,h}},
ImageSize->500
]
]
]

The action of Modular Transformations

Needs["ModularGroup`"];
Module[{mtList,mtSublist,zs,modCayley,gx,labelpt,mytiling,myttiling,mytfigure},
mtList=ModularGroupList[PSL2ZGrading@Mat@#<200&];
modCayley=New[MoebiusTransformation,{{I,1},{1,I}}];
mtSublist=Select[mtList,(GDiskMatRadius[GDiskMatMap[Mat@gdFord0,Mat@modCayley . Mat@#]]>1/10&&ModularGroup`Private`fordCriteria@Mat@#)&];
zs=#[ComplexInfinity]&/@mtSublist;
gx[phi_,psi_,t_]:=ModularTiling[mtList,phi.MatrixPower[psi,N@t],
IncircleMode->GDisk, IncircleStyle->Directive[EdgeForm[Gray],LightGray],IncircleThreshold->(2^(-6)),
TilingThreshold->(2^(-6)),
LabelMode->TUWord,LabelStyle->Black,LabelThreshold->1/12,Magnification->3];
(*mytiling[psi_,t_]:=Graphics[gx[IdentityMatrix@2,Mat@psi,t],
Axes->True,AxesStyle->Gray,PlotRange->{{-2,2},{0,2GoldenRatio}},
PlotLabel->(Superscript[psi,t]/.{mtT->T,mtU->U,mtR->R})];*)
labelpt[m_,z_]:=Module[{v,pt},
v=If[z===ComplexInfinity,{1,0},{z,1}];
pt={Re@#,Im@#}&@Inhom[m.v];
{Point[pt],Text[z/.{ComplexInfinity->∞,1/2->0.5,-1/2->-0.5},pt,-2*pt]}
];
mytiling[psi_,t_,fix_]:=Module[{m},
m=N@Mat[modCayley].MatrixPower[Mat@psi,t];
Graphics[{
gx[Mat@modCayley,Mat@psi,t],
{PointSize[Medium],labelpt[m,#]&/@zs},
{Red,PointSize[Medium],Point[{Re@#,Im@#}&@Inhom[m.fix]]}
},
(*Axes->True,AxesStyle->Gray,*)PlotRange->1.25,
PlotLabel->(Superscript[psi,t]/.{mtT->T,mtU->U,mtR->R})
]
];
mytfigure[psi_,fix_]:=GraphicsGrid[{
(*{myttiling[psi,0],myttiling[psi,1/6]},
{myttiling[psi,1/3],myttiling[psi,2/3]},
{myttiling[psi,5/6],myttiling[psi,1]}*)
{mytiling[psi,0,fix],mytiling[psi,1,fix]},
{mytiling[psi,1/5,fix],mytiling[psi,4/5,fix]},
{mytiling[psi,2/5,fix],mytiling[psi,3/5,fix]}
},ImageSize->400];
Column@{
mytfigure[mtU,{1,0}],
mytfigure[mtT,{I,1}],
mytfigure[mtR,{Exp[2 Pi I/3],1}]
}
]

Modular Tiling and Exponential Transformation

Needs["ModularGroup`"];
circ=Function[{m,r,t},m+r Exp[2 I Pi t],Listable];
halfCirc=Function[{m,r,t},m+r Exp[I Pi t],Listable];
arc=Function[{m,r,α,β,t},m+r Exp[I(α+(β-α)t)],Listable];

expHomotopy=Function[{p,z},(((I/(z+I)))^(1-p))Exp[2 Pi I p z],Listable];
expCoords=Compile[{{p,_Real},{z,_Complex}},
Evaluate@Simplify[{Re@#,Im@#}&@expHomotopy[p,z]],
RuntimeAttributes->Listable
];

expTilingPlot[mtList_,p_,opts:OptionsPattern[]]:=Module[{ms,disks,arcs},
ms=Mat/@mtList;
disks=GDiskMatMap[Mat@gdUnitDisk,ms];
arcs=Evaluate[Join[
{
1/2-I+I/#,-1/2-I+I/#,
arc[0,1,Pi/3,2Pi/3,#],
arc[1,1,2Pi/3,Pi,#],arc[-1,1,0,Pi/3,#],
arc[1/3,1/3,Pi/3,Pi,#], arc[-1/3,1/3,0,2Pi/3,#],
arc[2/3,1/3,2Pi/3,Pi,#],arc[-2/3,1/3,0,Pi/3,#]
},
halfCirc[GDiskMatCenter@disks,GDiskMatRadius@disks,#]
]]&;
ParametricPlot[expCoords[p,#]&/@arcs[t],{t,0,1},
Evaluate@FilterRules[{opts}, Options[ParametricPlot]]
]
];
expIncirclePlot[mtList_,p_,opts:OptionsPattern[]]:=Module[{ms,disks,circs},
ms=Mat/@mtList;
disks=GDiskMatMap[Mat@gdIncircle0,ms];
circs=Evaluate[circ[GDiskMatCenter@disks,GDiskMatRadius@disks,#]]&;
ParametricPlot[expCoords[p,#]&/@circs[t],{t,0,1},
Evaluate@FilterRules[{opts},Options[ParametricPlot]]
]
];
expIncircleLabels[mtList_,p_,labelFu_,displayFu_]:=Module[{ms,disks,centers,radii,north,south,labelSizes,labels},
ms=Mat/@mtList;
disks=GDiskMatMap[Mat@gdIncircle0,ms];
centers=GDiskMatCenter@disks;
radii=GDiskMatRadius@disks;
north=centers+I radii;
south=centers-I radii;
labelSizes=N@Abs[expHomotopy[p,north]-expHomotopy[p,south]];
labels=labelFu/@mtList;
Thread["f"[labels,N@expHomotopy[p,centers],labelSizes]]/."f"->displayFu
];

Module[{mtList,tilingList,incircleList,labelList,labelPlot,figure},
mtList=ModularGroupList[PSL2ZGrading@Mat@#<400&];
tilingList=Rest@Select[mtList,TRRight@#===mtT&&Abs@Re@#[I]≤1/2&&Im@#[I]≤1/6&];
incircleList=Select[mtList,Abs@Re@#[2I]≤1/2&];
labelList=Select[incircleList,GDiskMatRadius[GDiskMatMap[Mat@gdIncircle0,{{1,-I},{1,I}}.Mat@#]]≥1/50&];
SetOptions[TUExponents,QuotientFunction->CQuotient];
labelPlot[size_,p_]:=Graphics[
expIncircleLabels[labelList,p,TUWord,Text[Magnify[#1,Clip[size #3,{0.5,3.}]],{Re@#,Im@#}&@#2]&]
];
Show[{
expTilingPlot[tilingList,1,PlotStyle->Black],
expIncirclePlot[incircleList,1,PlotStyle->Gray],
labelPlot[4,1]
},ImageSize->500,PlotRange->1,Axes->True,AxesStyle->Gray]
]

Module[{nPlots=6,pStep,criteria,mtList,
tilingList,tilingPlots,
incircleList,incirclePlots,
labelList,labelPlot,zoomRectangles,
fanPlot,zoomPlot,zoomList,
figure1,figure2
},
pStep=1./(nPlots-1);
criteria=GDiskMatRadius[GDiskMatMap[Mat@gdIncircle0,{{I,1},{1,I}}.Mat@#]]≥(2^(-7))&;
mtList=ModularGroupList[criteria];

tilingList=Rest@Select[mtList,TRRight@#===mtT&&Abs@Re@#[I]≤1/2&&Im@#[I]≤1/6&];
tilingPlots=Table[expTilingPlot[tilingList,p,PlotStyle->Black],{p,0,1,pStep}];

incircleList=Select[mtList,Abs@Re@#[2I]≤1/2&];
incirclePlots=Table[expIncirclePlot[incircleList,p,PlotStyle->Gray],{p,0,1,pStep}];

labelList=Select[incircleList,GDiskMatRadius[GDiskMatMap[Mat@gdIncircle0,{{1,-I},{1,I}}.Mat@#]]≥(2^(-4))&];
SetOptions[TUExponents,QuotientFunction->CQuotient];
labelPlot[size_,p_]:=Graphics[
expIncircleLabels[labelList,p,TUWord,Text[Magnify[#1,Clip[size #3,{0.5,2.}]],{Re@#,Im@#}&@#2]&]
];

zoomRectangles=Table[Graphics[{EdgeForm[Directive[Red,Dashed]],Transparent,Rectangle[{-#,-#},{#,#}]}]&@Exp[-2Pi p],{p,0,1,pStep}];
fanPlot[i_]:=Show[incirclePlots[[i]],tilingPlots[[i]],labelPlot[3,(i-1)/(nPlots-1)],PlotRange->1,AxesOrigin->{0,0},AxesStyle->Gray,PlotLabel->(t==Module[{n=(i-1)/(nPlots-1)},If[n∈Integers,n,N@n]])];
zoomPlot[z_,nextZ_,opts:OptionsPattern[]]:=Show[incirclePlots[[6]],tilingPlots[[6]],labelPlot[3/z,1],
Graphics[{EdgeForm[Directive[Red,Dashed]],Transparent,Rectangle[{-#,-#},{#,#}]&@nextZ}],
PlotRange->z,AxesOrigin->{0,0},AxesStyle->Gray,FilterRules[{opts},Options@Graphics]];

(*figure1=Graphics[{
Inset[fanPlot[1],{0,3.0},{0,0},1],
Inset[fanPlot[2],{1.4,2.4},{0,0},1],
Inset[fanPlot[3],{0,1.8},{0,0},1],
Inset[fanPlot[4],{1.4,1.2},{0,0},1],
Inset[fanPlot[5],{0,0.6},{0,0},1],
Inset[fanPlot[6],{1.4,0.0},{0,0},1]
},PlotRange->{{-.6,2},{-.6,3.6}}];*)
figure1=GraphicsGrid[
Table[{fanPlot[i],fanPlot[i+1]},{i,1,5,2}],
Spacings->{Scaled[.5],Scaled[.25]},
ImageSize->500
];

zoomList=Append[Exp[-Pi/2Range[0,5]],1];
(*figure2=Graphics[{
Inset[zoomPlot[zoomList[[1]],zoomList[[2]]],{0,3.0},{0,0},1],
Inset[zoomPlot[zoomList[[2]],zoomList[[3]]],{1.1,2.4},{0,0},1],
Inset[zoomPlot[zoomList[[3]],zoomList[[4]]],{0,1.8},{0,0},1],
Inset[zoomPlot[zoomList[[4]],zoomList[[5]]],{1.1,1.2},{0,0},1],
Inset[zoomPlot[zoomList[[5]],zoomList[[6]]],{0,.6},{0,0},1],
Inset[zoomPlot[zoomList[[6]],1],{1.1,0},{0,0},1]
},PlotRange->{{-.6,1.7},{-.6,3.6}}];*)
figure2=GraphicsGrid[
{
{zoomPlot[zoomList[[1]],zoomList[[2]]],zoomPlot[zoomList[[2]],zoomList[[3]]]},
{zoomPlot[zoomList[[3]],zoomList[[4]]],zoomPlot[zoomList[[4]],zoomList[[5]]]},
{zoomPlot[zoomList[[5]],zoomList[[6]],Ticks->{{-#,#},{-#,#}}&@0.001],zoomPlot[zoomList[[6]],zoomList[[7]]]}
},Spacings->{Scaled[.5],Scaled[.25]},
ImageSize->500
];

Column@{
figure1,
figure2
}
]

Hyperbolic geometry, Canonical fundamental domains

?InvAssociatedMap

FindHypMotion[z_,w_]:=Module[{mat,γ,λ,μ,sol},
mat={{z,w},{1,1}}.{{λ,0},{0,μ}}.{{γ,1},{-1,γ}};
sol=FindInstance[{Flatten@mat∈Reals,Abs@γ==1},{γ,λ,μ}][[1]];
mat/.sol
];
ClearAll[GetHypMotion];
GetHypMotion[a_,b_]:=GetHypMotion[a,b]=FindHypMotion[a,b];
GetHypMotions[z_,grpList_]:=Table[GetHypMotion[z,t[z]],{t,grpList}];
CanonicalFunDom[z_,grpList_,phi_,opts:OptionsPattern[]]:=Module[{ts,a,ainv,phim},
ainv=InvAssociatedMap[z];(* maps z to the fundamental domain *)
a=PSL2ZInv[ainv];
phim=If[MatrixQ[phi],phi,Mat@phi];
ts=GetHypMotions[Inhom[ainv.{z,1}],grpList];
ts=phim.a.#&/@ts;
Graphics[{
{LightGray,GDisk[gdRightHalfplane,ts]},
{Black,GCircle[gdRightHalfplane,ts]},
{Gray,Dashed,GCircle[gdUnitDisk,ts]},
{Red,PointSize[Medium],Point[{Re@#,Im@#}&@phi[z]]},
{Black,PointSize[Medium],Point[{Re@#,Im@#}&/@Table[Inhom[phim.a.Mat[t].ainv.{z,1}],{t,grpList}]]},
{White,GDisk[gdLowerHalfplane,{phi}]}
},opts]
];

Needs["ModularGroup`"];
Module[{modGrpList},
modGrpList=ModularGroupList[PSL2ZGrading@Mat@#≤3&];
CanonicalFunDom[GoldenRatio I,Rest@modGrpList,mtId,PlotRange->{{-2,2},{0,4/GoldenRatio}},Axes->True,ImageSize->400]
]

Needs["ModularGroup`"];
Module[{modGrpList},
modGrpList=ModularGroupList[PSL2ZGrading@Mat@#≤3&];
CanonicalFunDom[(2+2I)/3,Rest@modGrpList,mtId,PlotRange->{{-1,2},{0,3/GoldenRatio}},Axes->True,ImageSize->400]
]

Needs["ModularGroup`"];
Module[{modGrpList,gammat2},
modGrpList=ModularGroupList[PSL2ZGrading@Mat@#≤18&];
gammat2=Select[modGrpList,MatchQ[Mod[Mat[#],2],{{_,0},{0,_}}]&];
Column[MatrixForm@Mat@#&/@gammat2];
CanonicalFunDom[I,Rest@gammat2,mtId,PlotRange->{{-2,2},{0,4/GoldenRatio}},Axes->True,ImageSize->400]
]

Needs["ModularGroup`"];
Module[{modGrpList,gammat2},
modGrpList=ModularGroupList[PSL2ZGrading@Mat@#≤18&];
gammat2=Select[modGrpList,MatchQ[Mod[Mat[#],2],{{_,0},{0,_}}]&];
Column[MatrixForm@Mat@#&/@gammat2];
CanonicalFunDom[1/2+I,Rest@gammat2,mtId,PlotRange->{{-2,3},{0,5/GoldenRatio}},Axes->True,ImageSize->400]
]

The Modular Group and Continued Fractions

Module[{modCayley,mtList,mtSublist,pt,zs},
modCayley=New[MoebiusTransformation,{{I,1},{1,I}}];
mtList=ModularGroupList[PSL2ZGrading@Mat@#<150&];
mtSublist=Select[mtList,(GDiskMatRadius[GDiskMatMap[Mat@gdFord0,Mat@modCayley . Mat@#]]>1/30&&ModularGroup`Private`fordCriteria@Mat@#)&];
zs=#[ComplexInfinity]&/@mtSublist;
figure=Graphics[{
ModularTiling[mtList,modCayley,
InteriorMode->GDisk,
TilingMode->GCircle,TilingStyle->White,
FordDiskMode->GDisk,FordDiskStyle->LightGray,FordLabelMode->(Inhom[(Mat@#).{1,0}]/.ComplexInfinity->∞&),FordLabelThreshold->1/20,
IncircleMode->GDisk,IncircleStyle->Directive[EdgeForm[Gray],White],
LabelMode->TUWord,LabelThreshold->1/12,Magnification->8],
{PointSize[Medium],
Function[{z},
pt={Re@#,Im@#}&@modCayley[z];
{Point[pt],Text[z/.ComplexInfinity->∞,pt,-2*pt]}
]/@zs
}
},PlotRange->1.15,ImageSize->400];
]

MyTicks[{min_,max_},special_List,n_Integer]:=Module[{ticks,spInRange},
spInRange=Select[special,IntervalMemberQ[Interval[{min,max}],#]&];
Join[
N@Select[FindDivisions[{min,max},n],And@@(Thread[Abs[#-spInRange]≥(max-min)/n])&],
spInRange
]
];
MyContinuedFraction[x_,n_]:=
Module[{i,c,t=x,cs={}},
For[i=0;c=Quotient[t,1,-1/2],
i<n&&t-c≠0,
t=1/(t-c);c=Quotient[t,1,-1/2];i++,
cs={cs,c}];
If[i<n,cs={cs,c}];
Flatten@cs
];
ModularConvergents=Compile[{{quotients,_Integer,1}},
Module[{sgn,i,m},
sgn=-1;
i=1;
m={{0,-1},{1,0}};
Table[sgn*=-1;m=m.{{0,-1},{1,sgn quotients[[i]]}},{i,1,Length[quotients]}]
]
];
Module[{mtList,number,quotients,convergents,Frame,plotRange,figure},
number=Pi;
mtList=ModularGroupList[PSL2ZGrading@Mat@#<100&];
quotients =MyContinuedFraction[number,4];
convergents=ModularConvergents[quotients];
Frame[m_,viewCenter_,zf_]:=Graphics[
plotRange={viewCenter+(GoldenRatio^2){-zf,zf},{0,2zf}};
{
ModularTiling[mtList,m,
FordDiskMode->GDisk,FordDiskStyle->LightGray,FordDiskThreshold->zf/25,
IncircleMode->GDisk,IncircleStyle->Directive[White,EdgeForm[Gray]],IncircleThreshold->zf/50,
TilingStyle->Gray, TilingThreshold->zf/50,
ExteriorMode->Off,
LabelMode->((((-1))^(Length@#+1))Last@#&@TUExponents[m.Mat@#,QuotientFunction->CentralQuotient]&),
Magnification->2/zf,LabelThreshold->(2^(-3))zf,
FordLabelMode->(Inhom[m.(Mat@#).{1,0}]&),FordLabelThreshold->(2^(-3))zf
],
{Red,Thick,GCircle[gdIncircle0,convergents]},
{White,Rectangle[Scaled@{-1,-1},Scaled@{2,0}]}
},
PlotRange->plotRange,Axes->True,AxesOrigin->{viewCenter-zf GoldenRatio,0},
Ticks->{MyTicks[{#1,#2},{number},6]&,Automatic}
];
figure=GraphicsColumn[Frame@@@{
{convergents[[1]],2,1},
{convergents[[2]],Pi,1/50},
{convergents[[3]],Pi,1/10000},
{convergents[[3]],Pi,1/25000}
},
ImageSize->400
]
]

Klein J Invariant Function under modified Cayley Transform

Graphics[{
{Raster[Transpose@Array[List@@ColorConvert[Hue[Mod[Arg[#1+I #2]/(2Pi),1]],"RGB"]&,{207,207},N@{{-1,1},{-1,1}}],{{-1,-1},{1,1}}]},
{White,GDisk[gdUnitCodisk]}
},
Axes->True,PlotRange->{{-1,1},{-1,1}}
]

phi=mtModCayley;(*New[MoebiusTransformation,{{I,1},{1,I}}];*)
f[y_,x_]:=Module[{c,z},
c=Complex[x,y];
z=(Inv@phi)@c;
If[Abs@c≥1.,0.,N@KleinInvariantJ@z]
];
color=Function[c,Module[{lc,t,b},
List@@ColorConvert[
Which[
c===ComplexInfinity,White,
c==0.,Black,
True,lc=Log@N@c;t=Mod[Im@lc/(2 Pi),1];b=ArcTan@Re@lc/Pi+.5;Hue[t,1-(b^2),1-(((1-b))^2)]
],"RGB"
]
],Listable];

jdata=Quiet[With[{n=513},Array[f,{n,n},N@{{-1,1},{-1,1}}]],Power::infy];

mtList=ModularGroupList[PSL2ZGrading@Mat@#<160&];
ModularPlot[data_,mts_,opts:OptionsPattern[]]:=Graphics[{
Raster[color@data,{{-1,-1},{1,1}}],
ModularTiling[mts,{{I,1},{1,I}},TilingStyle->Gray,BorderStyle->Gray]
},PlotRange->1,opts];

ModularPlot[jdata,mtList,ImageSize->400]

mtList2=Select[mtList,PSL2ZGrading@Mat@#<100&];
jdata2=Table[jdata[[x,y]],{x,1,Length@jdata,2},{y,1,Length@jdata,2}];

Quiet[ModularPlot[1/jdata2,mtList2,PlotLabel->Superscript[J,-1],ImageSize->200],Power::infy]

Quiet[ModularPlot[jdata2-1,mtList2,PlotLabel->J-1,ImageSize->200],Power::infy]

Quiet[ModularPlot[(jdata2^2),mtList2,PlotLabel->(J^2),ImageSize->200],Power::infy]

Quiet[ModularPlot[Function[z,phi[z],Listable][jdata2],mtList2,PlotLabel->Φ[J],ImageSize->200],Power::infy]

Quiet[ModularPlot[Function[z,If[z===ComplexInfinity,0.,2z/(1-z-(z^2))],Listable]@jdata2,mtList2,ImageSize->200],Power::infy]

Quiet[ModularPlot[Function[z,If[z===ComplexInfinity,0.,2z/(1-z-(z^2))],Listable]@jdata,mtList,ImageSize->400],Power::infy]

ModularPlot2[data_,mts_,opts:OptionsPattern[]]:=Graphics[{
Raster[Transpose@color@data,{{-2,0},{2,3}}],
ModularTiling[mts,IdentityMatrix@2,TilingStyle->Gray,BorderStyle->Gray]
},PlotRange->{{-2,2},{0,3}},opts];

jdataH=Quiet[With[{n=513},Array[N@KleinInvariantJ[#1+#2I]&,{n,n},N@{{-2,2},{1/n,3}}]],Power::infy];

ModularPlot2[jdataH,mtList,Axes->True]

Nice Generalized Disks under the Modular Group

?TRIndicateRight

Module[{mtList,gdT,gdR,phi,ts},
gdT=New[CompactDisk, Sqrt[5]I/2,1/2];
gdR=New[CompactDisk,-1/2+I,1/2];
phi={{I,1},{1,I}};
mtList=ModularGroupList[GDiskMatRadius[GDiskMatMap[Mat@gdT,phi.Mat@#]]≥(2^(-8))&];
ts=phi.Mat@#&/@mtList;
Graphics[
{
{Blue,GCircle[gdT,Rest@Pick[ts,Thread[TRIndicateRight[Mat/@mtList]==0]]]},
{Red,GCircle[gdR,Rest@Pick[ts,Thread[TRIndicateRight[Mat/@mtList]==1]]]},
ModularTiling[mtList,phi,TilingMode->GCircle,IncircleMode->GCircle,FordDiskMode->GCircle,FordDiskStyle->Brown,IncircleStyle->Gray]
},PlotRange->1]
]

Self tests

On[Assert];

Module[{phi,result,inf},
phi=New[MoebiusTransformation,{{I,1},{1,I}}];
Assert[InstanceQ[MoebiusTransformation][phi]];
result=phi[#]&/@{I,-I,ComplexInfinity};
result=result/.ComplexInfinity->inf;
Assert[result=={0,inf,I}];
result=phi[#]&/@{{I,1},{-I,1},{1,0}};
result=Inhom[#]&/@result;
result=result/.ComplexInfinity ->inf;
Assert[result=={0,inf,I}];
];

Module[{phi,psi},
phi=New[MoebiusTransformation,{{I,1},{1,I}}];
psi=New[MoebiusTransformation,{{1,-1},{0,1}}];
Assert[(phi@psi)@"z"==phi@(psi@"z")//Simplify];
];

Module[{phi},
phi=New[MoebiusTransformation,{{I,1},{1,I}}];
Assert[InstanceQ[MoebiusTransformation][phi]];
Assert[!InstanceQ[ModularTransformation][phi]];
Assert[Mat[phi]=={{I,1},{1,I}}];
Assert[Inv@Inv@phi===phi];
Assert[Inv@Inv@Inv@phi===Inv@phi];
Assert[phi@Inv[phi]@"z"==Inv[phi]@phi@"z"=="z"//Simplify];
];

Module[{phi},
phi=New[MoebiusTransformation,{{3,2},{7,5}}];
Assert[InstanceQ[MoebiusTransformation][phi]];
Assert[InstanceQ[ModularTransformation][phi]];
];

Assert[Table[CentralQuotient[i,5],{i,-5,5}]=={-1,-1,-1,0,0,0,0,0,1,1,1}];

Assert[Table[CQuotient[i,5],{i,-5,5}]==
{-1,0,0,0,0,0,0,0,0,0,1}];

Assert[TRLeft[mtId]===mtT];
Assert[TRRight[mtId]===mtT];
Assert[TRLeft[mtT]===mtT];
Assert[TRRight[mtT]===mtT];
Assert[TRLeft[mtR]===mtR];
Assert[TRRight[mtR]===mtR];
Assert[TRLeft[Inv@mtR]===Inv@mtR];
Assert[TRRight[Inv@mtR]===Inv@mtR];
Assert[TRLeft[mtU]===mtT];
Assert[TRRight[mtU]===mtR];
Assert[TRLeft[Inv@mtU]===Inv@mtR];
Assert[TRRight[Inv@mtU]===mtT];

Assert[TRList[mtU@mtU]==={mtT,mtR,mtT,mtR}];
Assert[TRList[Inv@mtU@mtU]==={Inv@mtR,mtT,Inv@mtR,mtT}];

Assert[TUExponents[mtT@mtU]=={0,1}];
Assert[TUExponents[mtU@mtT@Inv@mtU]=={1,-1}];
Assert[TUExponents[Inv[mtU]@mtT]=={-1,0}];

Module[{disks},
disks={gdUnitDisk, New[NoncompactDisk,0,1],gdUpperHalfplane};
Assert[InteriorQ[#,0]&/@disks=={True,False,True}];
Assert[InteriorQ[#,2]&/@disks=={False,True,True}];Assert[InteriorQ[#,I]&/@disks=={True,True,True}];
Assert[InteriorQ[#,-I]&/@disks=={True,True,False}];
Assert[InteriorQ[#,ComplexInfinity]&/@disks=={False,True,True}];
];

Module[{phi,disk,expectedType},
expectedType={Halfplane,CompactDisk,Halfplane,NoncompactDisk,Halfplane};
Do[
phi=New[MoebiusTransformation,MatrixPower[{{I,1},{1,I}},i]];
disk=phi[gdUpperHalfplane];
Assert[InstanceQ[expectedType[[i+1]]][disk]],
{i,0,4}
];
];

Off[Assert];