COMPLEX ANALYSIS: Mathematica 4.1 Notebooks
(c) John H. Mathews, and
Russell W. Howell, 2002
Complimentary software to accompany our textbook:

COMPLEX ANALYSIS: for Mathematics and Engineering, 
4th Edition,  2001,  ISBN: 0-7637-1425-9
Jones & Bartlett Publishers, Inc.  
40 Tall Pine Drive, Sudbury, MA 01776  
Tele.  (800) 832-0034, FAX:  (508)  443-8000  
E-mail:  mkt@jbpub.com  
Internet:  http://www.jbpub.com/  
This free software is compliments of the authors.
John H. Mathews,   mathews@fullerton.edu 
Russell W. Howell,  howell@westmont.edu

CHAPTER 10  APPLICATIONS of HARMONIC FUNCTIONS

a = A = b = B = c = k = R = T = u = U = v = V = w = W = x = X = y = Y = z = Z = 
                                        0; 
Clear[a, a0, a1, a2, a3, A, b, B, 
        CR, Δx, Δy, eq1, eq2, eqns, f, f1, F, F1, 
        F2, F3, F4, g, G, lines, n, P, Ψ, ϕ, Φ, 
        r, r0, solset, θ, t, T, T1, T2, Q, u, u0, U, U1, v, v0, V, V1, 
        w, w0, W, W1, x, x1, x2, x3, y, z, z1, z2, z3, 
        z4]; 
Needs["Graphics`ComplexMap`"]; 
Needs["Graphics`ParametricPlot3D`"]
Needs["Graphics`Colors`"];

Section 10.1	Preliminaries

Section 10.2	Invariance of Laplace's Equation and the Dirichlet Problem

Section 10.3	Poisson's Integral Formula for the Upper Half Plane

Section 10.4	Two Dimensional Mathematical Models

Section 10.5	Steady State Temperatures

Section 10.6  Two Dimensional Electrostatics

A two-dimensional electrostatic field is produced by a system of charged wires, plates, and cylindrical conductors that are perpendicular to the z-plane.  The wires, plates, and cylinders are assumed to be so long that the effects at the ends can be neglected as mentioned in Section 10.4.  This sets up an electric field  E(x, y)  that can be interpreted as the force acting on a unit positive charge placed at the point (x,y).  In the study of electrostatics the vector field  E(x, y)  is shown to be conservative and is derivable from a function  ϕ (x, y),  called the electrostatic potential, as expressed by the equation  

	E(x, y)  = - grad ϕ (x, 
              y)   =  -  (Subscript[ϕ, x]) (x, 
                y)  -  i  (Subscript[ϕ, y]) (x, y).

	If we make the additional assumption that there are no charges within the domain  D,  then Gauss' Law for electrostatic fields implies that the line integral of the outward normal component of   E(x, y)  taken around any small rectangle lying inside  D  is identically zero.  A heuristic argument similar to the one for steady state temperatures with  T (x, y)  replaced by   ϕ(x, y)  will show that the value of the line integral is  -[(Subscript[ϕ, xx]) (x, y)  +   (Subscript[ϕ, yy]) (x, 
                  y)] Δx Δy.  Since this latter quantity is zero, we conclude that  ϕ(x, y)  is a harmonic function.  We let  ψ(x, y)  denote the harmonic conjugate, and 

	  F (z)  =  ϕ (x, y)  +  i ψ (x, 
              y)

is the complex potential (not to be confused with the electrostatic potential).

The curves  ϕ (x, y) = (Subscript[K, 1])  are called equipotential curves, The curves  ψ (x, y) = (Subscript[K, 2])  are called the lines of flux.   If a small test charge is allowed to move under the influence of the field  E(x, y),  then it will travel along a line of flux.

Example 10.18,  Page 425.  Consider two parallel conducting planes that pass perpendicular to the z plane through the lines  x = a  and  x = b,  which are kept at the potentials  (Subscript[U, 1]) and (Subscript[U, 2]), respectively.   Then according to the result of Example 10.1, the electrical potential is ϕ (x, y) = 
        (Subscript[U, 1]) + (((Subscript[U, 2])  -  (Subscript[U, 1]))/(b   -   a)) (x - a).

Solution 10.18.

This example is similar to Example 10.1, Page 310  and  Example 10.14,  Page 331.  For illustration we use the boundary values  
\(ϕ \((\(-1\), y)\) = 30, 
ϕ \((1, y)\) = \(-20. \)\)

Enter the formula for the electrical potential.

Remove[a, b, c, x, y, 
        U]; 
Clear[cont, ϕ]; 
ϕ[x_, y_]  =  
        (Subscript[U, 1])  +  (((Subscript[U, 2])  -  (Subscript[U, 1]))/(b  -  a)) (x  -  
                a); 
Print["ϕ[x,y] = ", ϕ[x, 
          y]];

Enter the boundary values  a = -1,   b = 
        1,   (Subscript[U, 1]) = 30  and  (Subscript[U, 2]) = -20,  and check out the formula.

a  = -1; 
b  =  1; 
(Subscript[U, 1])  =  30; 
(Subscript[U, 2])  =  -20; 
Print["ϕ[x,y] = ", ϕ[x, y]]; 
Print["ϕ[-1,y] = ", ϕ[-1, y]]; 
Print["ϕ[0,y] = ", ϕ[0, y]]; 
Print["ϕ[1,y] = ", ϕ[1, 
          y]];

Use Mathematica to make a contour plot of the electrical potential.

cont  =  
        Table[c, {c, -20, 30, 5}]; 
ContourPlot[ϕ[x, y], {x, -1, 1}, {y, 0, 
          2}, 
Contours  ->  cont, 
        PlotPoints -> 10, 
AspectRatio -> 1, 
        ColorFunction -> Hue]; 
Print["Contour plot for"]; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print[""]; 
Print[ϕ[x, y] == c]; 
Print["where  c = ", cont];

Hence the electrical potential  ϕ[x, y] = 30 - 25 (1 + x)  is harmonic in the horizontal strip  -1 < Re[z] < 1  and has the desired boundary values.

Example 10.19,  Page 425.  Find the electrical potential  ϕ (x, y)  in the region between the two infinite coaxial cylinders of radius  (Subscript[R, 1]) = a  and  (Subscript[R, 2]) = b, which are kept at potentials  (Subscript[U, 1])  and  (Subscript[U, 2]),  respectively.

Solution 10.19.

This example is similar to Example 10.3,  Page 311.  For illustration we use the boundary values  

\(ϕ \((x, y)\) = \ \ \(0\ \ \ when\ \ \  | z |  = 
          1\), 
ϕ \((x, y)\) = \(40\ \ \ when\ \ \  | z |  = 2. \)\)

Enter the formula for the electrical potential.

Remove[x, y, R, U]; 
Clear[cont, ϕ]; 
ϕ[x_, y_] = 
        (Subscript[U, 1]) + (((Subscript[U, 2])  -  (Subscript[U, 1]))/(Log[(Subscript[R, 2])] - 
                  Log[(Subscript[R, 1])])) (Log[(Sqrt[(x^2) + (y^2)])] - 
                Log[(Subscript[R, 1])]); 
Print["ϕ[x,y] = ", ϕ[x, 
          y]];

Enter the boundary values  a = 1,   b = 2,   (Subscript[U, 1]) = 0  and  (Subscript[U, 2]) = 40,  and check out the formula.

\(
\(R\_1\ \  = \ 1;\)\ \)
(Subscript[R, 2])   =  2; 
(Subscript[U, 1])  =  0; 
(Subscript[U, 2])  =  40; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print[""]; 
Print["ϕ[1,0] = ", ϕ[1, 0]]; 
Print["ϕ[0,1] = ", ϕ[0, 1]]; 
Print["ϕ[2,0] = ", ϕ[2, 0]]; 
Print["ϕ[0,2] = ", ϕ[0, 
          2]];

Use Mathematica to make a contour plot of the electrical potential.

cont  =  
        Table[0 + 10  k, {k, 0, 4}]; 
ContourPlot[ϕ[x, y], {x, -2, 2}, {y, -2, 
          2}, 
Contours  ->  cont, 
        PlotPoints -> 25, 
AspectRatio -> 1, 
        ColorFunction -> Hue]; 
Print["Contour plot for"]; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print[""]; 
Print[ϕ[x, y] == c]; 
Print["where  c = ", cont];

Hence the electrical potential  ϕ[x, y] = ((40 Log[(Sqrt[(x^2)  +  (y^2)])])/Log)[2]  is harmonic in the annulus  1< | z | <2  and has the desired boundary values.

Example 10.20,  Page 425.   Find the electrical potential  ϕ (x, y)  produced by two charged half-planes that are perpendicular to the z plane and pass through the rays  x < -1, y = 0  and  x > 1, y = 0  where the planes are kept at the fixed potentials  

\(ϕ \((x, 0)\) = \(-300\)\ \ \ for\ \ \ x < \(-1\), 
ϕ \((x, 
            0)\) = \ \ 300\ \ \ \ for\ \ \ x > 1. \)

Solution 10.20.

Enter the formula for the electrical potential.

Remove[u, x, y]; 
Clear[cont, ϕ]; 
ϕ[x_, y_]  =  
        (600/π) Re[ArcSin[x  +  I y]]; 
Print["ϕ[x,y] = ", ϕ[x, y]]; 
Print["ϕ[x,y] = ", 
        ComplexExpand[ϕ[x, y]]]; 
ϕ[x_, y_]  =  
        (600/π) ArcSin[(((Sqrt[(((x + 1))^2)  +  (y^2)])  - (Sqrt[(((x - 1))^2)  +  (y^2)]))/2)]; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]];

Check out the formula at some boundary values.

Print[" ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print["ϕ[-1,0] = ", ϕ[-2, 
          0]]; 
Print["ϕ[-1,0] = ", ϕ[-1, 0]]; 
Print[" ϕ[0,1] = ", ϕ[0, 1]]; 
Print[" ϕ[1,0] = ", ϕ[1, 
          0]]; 
Print[" ϕ[1,0] = ", ϕ[2, 
          0]];

Use Mathematica to make a contour plot of the electrical potential.

ϕ[x_, y_]  =  
        (600/π) Re[
            ArcSin[x  +  I y]]; 
cont  =  Table[-300  +  60  k, {k, 0, 10}]; 
ContourPlot[ϕ[x, y], {x, -4, 4}, {y, -4, 
          4}, 
Contours  ->  cont, 
        PlotPoints -> 50, 
AspectRatio -> 1, 
        ColorFunction -> Hue]; 
Print["Contour plot for"]; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print[""]; 
Print[ϕ[x, y] == c]; 
Print["where  c = ", cont];

The following commands will graph the contours as parametric curves.

Remove[X, Y, z, Z, 
        w]; 
Clear[F, solset]; 
F[z_]  =  (600/π) ArcSin[z]; 
solset  =  Solve[w  ==  F[z], z]; 
Z[u_, v_]  =  
        ReplaceAll[
          (Subscript[solset, [[1, 1, 2]]]), 
          w -> u  +  I v]; 
X[u_, v_]  =  ComplexExpand[Re[Z[u, v]]]; 
Y[u_, v_]  =  ComplexExpand[Im[Z[u, v]]]; 
Print["w = ", F[z]]; 
Print["Z = ", Z[u, v]]; 
Print["Z = ", 
        ComplexExpand[Z[u, v]]]; 
Print[""]; 
Print["X = ", X[u, v]]; 
Print["Y = ", Y[u, v]];

Clear[k]; 
eqns  =  
        N[Table[{X[k 60, v], Y[k 60, v]}, {k, -5, 5}]]; 
ParametricPlot[
        Evaluate[eqns], {v, -500, 
          500}, 
PlotRange -> {{-4, 4}, {-4, 
              4}}, AspectRatio -> 1]; 
Print["ϕ[x,y] = ", ϕ[x, 
          y]]; 
Print[""]; 
Print[ϕ[x, y] == c]; 
Print["where  c = ", cont];

Hence the electrical potential  ϕ[x, y] = ((600 Re[ArcSin[x + i y]])/π)  is harmonic in plane slit along the rays  x < -1, y = 0  and  x > 1, y = 0  and has the desired boundary values.

Example 10.21,  Page 426.  Find the electrical potential  ϕ (x, y)  in the disk  D : |z|<1,   that takes on the boundary values  

\(ϕ \((x, 
            y)\) = \(80\ \ \ for\ \ \ z = e\^\(i\ θ\)\ \ \ and\ \ 0 < θ < π\/2\), 
ϕ \((x, y)\) = \ \ \(0\ \ \ for\ \ \ z = e\^\(i\ θ\)\ \ \ and\ \ π\/2 < θ < 2 \( \(π\)\(.\)\)\)\)

Solution 10.21.

Construct the solution with a conformal mapping.  This example is similar to Example 10.8, Page 317.

Remove[u, U, v, V, w, W, x, X, y, Y, z, 
        Z]; 
Clear[f, F, ϕ]; 
f[
          z_]  =  (((1  -  I) (z  -  I))/(z  -  1)); 
W  =   ComplexExpand[f[x  +  I y]]; 
U  =  
        Together[
          ComplexExpand[Re[W], 
            TargetFunctions -> {Re, Im}]]; 
V  =  
          Together[
            ComplexExpand[Im[W], 
              TargetFunctions -> {Re, Im}]]; ; 
ϕ[x_, y_]  =  (80/π) ArcTan[(V/U)]; 
Print["w = ", f[z]]; 
Print["w = u + \[ImaginaryI] v = ", 
        U  +  I V];  
Print[" "]; 
\(\(Print[\*""ϕ[x,y] = (80/π) ArcTan[(v/u)]""];\)\ \)
Print["ϕ[x,y] = ", ϕ[x, 
          y]];

Use Mathematica to make a contour plot of the temperature.  For computational purposes we use the following formula   \(80\/π\ ArcTan[\(1\  - \ x\^2\  - \ y\^2\)\/\(1\  - \ 2\ x\  + \ x\^2\  - \ 2\ y\  + \ y\^2\)] = \(-80\)\/π\ \((\(-π\)\/2 + 
              ArcTan[1 - x\^2 - y\^2, 1 - 2\ x + x\^2 - 2\ y + y\^2])\)\).

ϕ[x_, 
          y_] = ((-80)/π) (((-π)/2) + 
              ArcTan[1 - (x^2) - (y^2), 
                1 - 2 x + (x^2) - 2 y + 
                  (y^2)]); 
cont  =  Table[10  k, {k, 0, 8}]; 
ContourPlot[ϕ[x, y], {x, -1, 1}, {y, -1, 1}, 
        Contours  ->  cont, 
PlotPoints -> 60, 
        AspectRatio -> 1, 
        ColorFunction -> Hue]; 
Print["Contour plot for"]; 
Print["ϕ[x,y] = ", ϕ[x, y]]; 
Print[""]; 
Print[ϕ[x, y] == c]; 
Print["where  c = ", cont];

The following commands will graph the contours as parametric curves.

Remove[u, v, w, X, Y, z, 
        Z]; 
Clear[F]; 
F[
          z_]  =  (((1  -  I) (z  -  I))/(z  -  1)); 
solset  =  Solve[w  ==  F[z], z]; 
Z[u_, v_]  =  
        ReplaceAll[
          (Subscript[solset, [[1, 1, 2]]]), 
          w -> u  +  I v]; 
\(\(X[u_, v_]\  = \ 
        Together[
          ComplexExpand[
            Re[Z[u, v]], 
\ \ \ \ \ \ \ \ \ \ \ TargetFunctions -> {Re, \ Im}]];\)\ \)
\(\(Y[u_, v_]\  = \ 
        Together[
          ComplexExpand[
            Im[Z[u, v]], 
\ \ \ \ \ \ \ \ \ \ \ TargetFunctions -> {Re, \ Im}]];\)\ \)
Print["w = ", F[z]]; 
\(\(Print["z = ", 
        solset\_\(\([[\)\(1, 1, 2\)\(]]\)\)];\)\ \)
Print["Z = ", Z[u, v]]; 
Print[""]; 
Print["X = ", X[u, v]]; 
Print["Y = ", Y[u, v]];

n  =  8; 
Clear[k]; 
x[r_]  =  
        Simplify[
          X[r Cos[((k π)/n)], r Sin[((k π)/n)]]]; 
y[r_]  =  
        Simplify[
          Y[r Cos[((k π)/n)], r Sin[((k π)/n)]]]; 
eqns  =  N[Table[{x[r], y[r]}, {k, 0, n}]]; 
ParametricPlot[Evaluate[eqns], {r, 0.05, 40}, 
        PlotRange -> {{-1.1, 1.1}, {-1.1, 1.1}}, 
        AspectRatio -> 1];

Hence the electrical potential  \(ϕ[x, y] = 
        80\/π\ ArcTan[\(1 - x\^2 - y\^2\)\/\(1 - 2\ x + x\^2 - 2\ y + y\^2\)]\)  is harmonic in the disk |z| < 1 and has the desired boundary values.

Library Research Experience for Undergraduates

Project I. Write a report on two dimensional electrostatics.

1. Arndat, Richard A. and L. David Roper, (1989), ''Scattering in a Spherical Potential: Motion of Complex-Plane Poles and Zeros,'' Computers in Physics, pp. 65-72.Choudary, A. R. and B. H. Kellett and C. Verzegnassi, (1976), ''A Conformal-Mapping Evaluation of the Pion Form Factor in the Spacelike Region,'' Lettere al Nuovo Cimento, Ser. 2, Vol. 15., No. 16, pp. 579-582.

2. Foster, K. and R. Anderson, (1974), ''Transmission-Line Properties by Conformal Mapping,'' Proceedings IEE, Vol. 121, No. 5, pp. 337-339.

3. Li, Z.-M. and S. P. McAlister, (1989), ''Conformal Mapping and the Current Distributions in (Subscript[Al, x]) (Subscript[Ga, 1 - x])/As Ga As Heterostructures,'' J. App. Phys., V. 66, No. 2, pp. 680-685.

4. Mathis, H.F., (1963), ''Conformal Transformation Charts Used by Electrical Engineers'' Math. Mag., Vol.36, pp. 25-30.

5. Piele, Donald T., Morris W. Firebaugh and Robert Manulik, ''Applications of Conformal Mapping to Potential Theory Through Computer Graphics,'' Am. Math. M., Vol. 84, No. 9. (Nov., 1977), pp. 677-692., Jstor.   

6. Potter, Robert I. and Robert Schmulian and Keith Hartman, (1971), ''Fringe Field Readback Voltage Computations for Finite Pole-Tip Length Recording Heads,'' IEEE Transactions on Magnets,\ Vol. 7, No. 3, pp. 689-695.

7. Wu, Ruey-Beei and Chun Hsiung Chen, (1986), ''A scalar variational conformal mapping technique for weakly guiding dielectric waveguides,'' IEEE J. of Quantum Electronics, V. 22, pp. 603-609.

Section 10.6  Exercises for Two Dimensional Electrostatics
See textbook page 428.

Section 10.7	Two-Dimensional Fluid Flow

Section 10.8	The Joukowski Airfoil

Section 10.9	The Schwarz-Christoffel Transformation

Section 10.10	Image of a Fluid Flow

Section 10.11	Sources and Sinks

GO TO  THE  Table of Contents