1--Copyright The Numerical Algorithms Group Limited 1994.
2-- Drawing conformal maps.
3
4-- The functions in this file draw conformal maps both on the
5-- complex plane and on the Riemann sphere.
6
7-- Compile, don't interpret functions.
8)set fun comp on
9
10C := Complex DoubleFloat                -- Complex Numbers
11S := Segment DoubleFloat                -- Draw ranges
12R3 := POINT DoubleFloat                         -- points in 3-space
13
14-- conformalDraw(f, rRange, tRange, rSteps, tSteps, coord)
15-- draws the image of the coordinate grid under f in the complex plane.
16-- The grid may be given in either polar or cartesian coordinates.
17-- parameter descriptions:
18--   f:  the function to draw
19--   rRange: the range of the radius (in polar) or real (in cartesian)
20--   tRange: the range of theta (in polar) or imaginary (in cartesian)
21--   tSteps, rSteps: the number of intervals in each direction
22--   coord: the coordinate system to use.  Either "polar" or "cartesian"
23
24conformalDraw: (C -> C, S, S, PI, PI, String) -> VIEW3D
25conformalDraw(f, rRange, tRange, rSteps, tSteps, coord) ==
26  transformC :=
27    coord = "polar" => polar2Complex
28    cartesian2Complex
29  cm := makeConformalMap(f, transformC)
30  sp := createThreeSpace()
31  adaptGrid(sp, cm, rRange, tRange, rSteps, tSteps)
32  makeViewport3D(sp, "Conformal Map")
33
34-- riemannConformalDraw(f, rRange, tRange, rSteps, tSteps, coord)
35-- draws the image of the coordinate grid under f on the Riemann sphere.
36-- The grid may given in either polar or cartesian coordinates.
37-- parameter descriptions:
38--   f:  the function to draw
39--   rRange: the range of the radius(in polar) or real (in cartesian)
40--   tRange: the range of theta (in polar) or imaginary (in cartesian)
41--   tSteps, rSteps: the number of intervals in each direction
42--   coord: the coordinate system to use. either "polar" or "cartesian"
43
44riemannConformalDraw: (C -> C, S, S, PI, PI, String) -> VIEW3D
45riemannConformalDraw(f, rRange, tRange, rSteps, tSteps, coord) ==
46  transformC :=
47    coord = "polar" => polar2Complex
48    cartesian2Complex
49  sp := createThreeSpace()
50  cm := makeRiemannConformalMap(f, transformC)
51  adaptGrid(sp, cm, rRange, tRange, rSteps, tSteps)
52  -- add an invisible point at the north pole for scaling
53  curve(sp, [point [0,0,2.0@DoubleFloat,0], point [0,0, 2.0@DoubleFloat,0]])
54  makeViewport3D(sp, "Conformal Map on the Riemann Sphere")
55
56-- Plot the coordinate grid using adaptive plotting for the coordinate
57-- lines, and drawing tubes around the lines.
58adaptGrid(sp, f, uRange, vRange,  uSteps, vSteps) ==
59  delU := (high(uRange) - low(uRange))/uSteps
60  delV := (high(vRange) - low(vRange))/vSteps
61  uSteps := uSteps + 1; vSteps := vSteps + 1
62  u := low uRange
63  -- draw the coordinate lines in the v direction
64  for i in 1..uSteps repeat
65    -- create a curve 'c' which fixes the current value of 'u'
66    c := curryLeft(f,u)
67    cf := (t:DoubleFloat):DoubleFloat +-> 0
68    -- draw the 'v' coordinate line
69    makeObject(c, vRange::Segment Float, colorFunction == cf, space == sp, _
70               tubeRadius == 0.02,  tubePoints == 6)
71    u := u + delU
72  v := low vRange
73  -- draw the coordinate lines in the u direction
74  for i in 1..vSteps repeat
75    -- create a curve 'c' which fixes the current value of 'v'
76    c := curryRight(f,v)
77    cf := (t:DoubleFloat):DoubleFloat +-> 1
78    -- draw the 'u' coordinate line
79    makeObject(c, uRange::Segment Float, colorFunction == cf, space == sp, _
80               tubeRadius == 0.02,  tubePoints == 6)
81    v := v + delV
82  void()
83
84-- map a point in the complex plane to the Riemann sphere.
85riemannTransform(z) ==
86  r := sqrt norm z
87  cosTheta := (real z)/r
88  sinTheta := (imag z)/r
89  cp := 4*r/(4+r^2)
90  sp := sqrt(1-cp*cp)
91  if r>2 then sp := -sp
92  point [cosTheta*cp, sinTheta*cp, -sp + 1]
93
94-- convert cartesian coordinates to cartesian form complex
95cartesian2Complex(r:DoubleFloat, i:DoubleFloat):C == complex(r, i)
96
97-- convert polar coordinates to cartesian form complex
98polar2Complex(r:DoubleFloat, th:DoubleFloat):C == complex(r*cos(th), r*sin(th))
99
100-- convert a complex function into a mapping from (DoubleFloat,DoubleFloat) to R3 in the
101-- complex plane.
102makeConformalMap(f, transformC) ==
103  (u:DoubleFloat,v:DoubleFloat):R3 +->
104    z := f transformC(u, v)
105    point [real z, imag z, 0.0@DoubleFloat]
106
107-- convert a complex function into a mapping from (DoubleFloat,DoubleFloat) to R3 on the
108-- Riemann sphere.
109makeRiemannConformalMap(f, transformC) ==
110  (u:DoubleFloat, v:DoubleFloat):R3 +-> riemannTransform f transformC(u, v)
111
112-- draw a picture of the mapping of the complex plane to the Riemann sphere.
113riemannSphereDraw: (S, S, PI, PI, String) -> VIEW3D
114riemannSphereDraw(rRange, tRange, rSteps, tSteps, coord) ==
115  transformC :=
116    coord = "polar" => polar2Complex
117    cartesian2Complex
118  grid := (u:DoubleFloat , v:DoubleFloat): R3 +->
119    z1 := transformC(u, v)
120    point [real z1, imag z1, 0]
121  sp := createThreeSpace()
122  adaptGrid(sp, grid, rRange, tRange, rSteps, tSteps)
123  connectingLines(sp, grid, rRange, tRange, rSteps, tSteps)
124  makeObject(riemannSphere, 0..2*%pi, 0..%pi, space == sp)
125  f := (z:C):C +-> z
126  cm := makeRiemannConformalMap(f, transformC)
127  adaptGrid(sp, cm, rRange, tRange, rSteps, tSteps)
128  makeViewport3D(sp, "Riemann Sphere")
129
130-- draw the lines which connect the points in the complex plane to
131-- the north pole of the Riemann sphere.
132connectingLines(sp, f, uRange, vRange, uSteps, vSteps) ==
133  delU := (high(uRange) - low(uRange))/uSteps
134  delV := (high(vRange) - low(vRange))/vSteps
135  uSteps := uSteps + 1; vSteps := vSteps + 1
136  u := low uRange
137  -- for each grid point
138  for i in 1..uSteps repeat
139    v := low vRange
140    for j in 1..vSteps repeat
141      p1 := f(u,v)
142      p2 := riemannTransform complex(p1.1, p1.2)
143      fun := lineFromTo(p1,p2)
144      cf := (t:DoubleFloat):DoubleFloat +-> 3
145      makeObject(fun, 0..1, space == sp, tubePoints == 4, tubeRadius == 0.01,
146                 colorFunction == cf)
147      v := v + delV
148    u := u + delU
149  void()
150
151riemannSphere(u,v) ==
152  sv := sin(v)
153  0.99@DoubleFloat*(point [cos(u)*sv, sin(u)*sv, cos(v),0.0@DoubleFloat]) +
154    point [0.0@DoubleFloat, 0.0@DoubleFloat, 1.0@DoubleFloat, 4.0@DoubleFloat]
155
156-- create a line function which goes from p1 to p2 as its parameter
157-- goes from 0 to 1.
158lineFromTo(p1, p2) ==
159  d := p2 - p1
160  (t:DoubleFloat):Point DoubleFloat +-> p1 + t*d
161