stakgold chap 6 raw 1

Stakgold Chapter 6 Raw Notes PhL begun 7.16.09 Chapter Title: Potential Theory 3 6.1 Introduction (88) 4 6.2 Interior Dirichlet Problem for the Unit Circle 4 Separation of Variables [ 2D, r φ ] 4 Verification of the Solution of (6.5). (95) 11 Verification of the Solution of (6.4). 12 6.3 Some properties of harmonic functions 13 A. Mean Value Theorem and the "part problem method" 13 B. The Maximum Principle for Harmonic Functions (101) 16 C. Uniqueness and BC dependence continuity of the solution 16 Exercises: (103) 17 Exercise 6.1. Compute four infinite sums. 17 Exercise 6.2. Evaluate a Poisson Kernel integral with f = 1, then interpret the result. (2D) 17 Exercise 6.3. Solutions outside the unit circle, exterior Dirichlet. (2D) 17 Exercise 6.4. Symmetric polygon problem. (2D) 17 Exercise 6.5. Work things through for Helmholtz in place of Laplace. (2D) (105) 18 Exercise 6.6. The 4 guy 19 Exercise 6.7. Here L = 2 - P where P is a "positive function". (2D) 19 Exercise 6.8. Show uniqueness for Lu = (k u) . 22 Exercise 6.9. Poisson Kernel and the distributional Dirichlet problem (2D) 23 Exercise 6.10. Dirichlet Laplace on a wedge (2D) 23 Exercise 6.11. Dirichlet Laplace on an annulus. (2D) 26 Exercise 6.12. Dirichlet Laplace on a curved quadrilateral. (2D) 29 Exercise 6.13. Dirichlet Laplace on a rectangle. (2D) (108) 33 Exercise 6.14. Complete orthogonal set 35 Exercise 6.15. Laplace solution inside a sphere whose surface has f(θ,φ) (3D) 37 Exercise 6.16. Laplace solution outside a sphere whose surface has f(θ,φ) (3D) 39 Exercise 6.17. Laplace solution for spherical annulus (3D) 39 General Comment: Cauchy, Dirichlet and Neumann BC's 41 Comment on Theory so Far versus Real World Problems 42 6.4 Surface Layers (110) 42 A1. Finding u from a source density function q (n=3) (111, 6.30 resulting u) 43 B1. Suppose the source density function represents an isolated unit dipole. (112, 6.29) 43 A2. Suppose our charge q(x) is spread on a surface σ (simple layer) (112, result 6.30 for u) 45 A3. Proof that u of A2 is finite for x actually on the surface (112) 45 B2. A dipole layer on surface σ (113 result 6.31 for v) 46 B3. Proof that u of B2 is finite for x actually on the surface. (113) 46 C. Detail about x → s where s is on the surface σ which has a δ source. (113-115) 46 D. Extend the Part C analysis from δ sources to layers on the x-y plane (115) 47 E. Now consider the actual curved surface σ. ( 117) 48 Exercises (120) 50 Exercise 6.18 In 2D space, find u and v due to layers on σ = C, a curve. (120) 50 Exercise 6.19 In 2D space, compute the two distributional results for E (121) 51 Exercise 6.20 More 2D results (121) 52 6.5 Integral Equations of Potential Theory (122) 52 A. The double-layer integral equation method for solving Dirichlet (122 3D) 53 Question about the reality of the dipole layer 53 B. A way to approach the combined interior/exterior Dirichlet problem. (123 nD) 54 Example (125) : Simultaneous computation of u in and out for the unit sphere (3D) 56 Question: Why did our integral equation diagonalize to such a simple form? 57 The Neumann Problem in 3D and for general n ≠ 2 (126) 58 The Dirichlet and Neumann Problems in 2D (128) 58 Exercises (130) 59 Exercise 6.21. Interior Dirichlet unit sphere using double-layer integral equation method. 59 Part A. 59 Part B. 60 Part C. 62 Let's review the threaded steps: 65 Exercise 6.22. Interior Dirichlet unit disk using double-layer integral equation method. (130) 66 Exercise 6.23 Neumann on 2D unit circle (130) 74 Review of Chapter 6 up to this point. 75 Comment on Integral Equations 76 6.6 Green's Function for L = -2 (130) 77 Definition and Properties: Theorems 1 thru 6. (130) 77 Solution of the Dirichlet Problem. (135) 79 Eigenvalue Problem for L = – 2 (136) * 82 Two Examples: Rectangle and Disk (139) 82 Example 1: 82 Example 2: 83 Green's Function for Unbounded Regions (142) 83 Exercises (143) -- eleven problems here 84 Exercise 6.24. State Theorem 4 for n = 2 dimensions. 84 Exercise 6.25. Show that G is H-S in n = 2 dimensions. 84 Exercise 6.26. Show that the kernel k = a(x,ξ)/Rm (m < n) generates a CC operator G. 85 Exercise 6.27. Show that g(x|ξ) in 6.90 is symmetric (ie, for the Re external problem). 85 Exercise 6.28. Find the Laplace eigenfunctions for a 2D disk wedge r = 1 86 Exercise 6.29. Find the Laplace eigenfunctions for a 2D annular ring between r = a and r = b. 88 Comments on Neumann's. 89 Exercise 6.30. Eigenfunctions of Laplace for the unit sphere in 3D. 90 Exercise 6.31. Eigenfunctions of Laplace for a cone of the unit sphere in 3D. 91 Exercise 6.32. Eigenfunctions of Laplace for an orange slice of the unit sphere in 3D. 95 Exercise 6.33. Laplace on R where ∂nu = 0 (Neumann) on σ; apply to disk. (145) 96 Exercise 6.34. Repeat last but with "mixed" BC; apply to the disk (145) 97 6.7 Methods for finding the Green's Function (146) 98 (1) The Integral Equation Method 98 Example: the unit circle. (147) 100 Contradiction? 103 Resolution: 104 Plots of the 2D equipotential lines. 106 Comment on the 3D dipole equipotential lines. 106 (2) The Method of Images 110 Examples 110 Example 1: Green's function for the unit sphere (3D) (150) 110 Example 2: Green's function for a half space (3D) 112 Example 3: Neumann's function for a half space (3D) 113 Example 4: Line charge parallel to a plane (3D) 114 (3) The Method of Full Eigenfunction Expansion (153) 115 (4) The Method of Partial Eigenfunction Expansion (154) 115 Long Comment on Example 1 (rectangle in 2D) to come: 115 Example 1. Green's function for the rectangle. (155) 118 Example 2. Green's function for the unit circle. (159) 121 Example 3. Green's function for two horizontal plates with a line charge (2D strip). (161) 122 (5) Complex Variable Method for 2D (164) 123 Exercises (166) 125 Chapter Title: Potential Theory Comment: In the first three sections of this chapter, we look for solutions to Laplace 2u=0 where we specify the solution u on a closed bounding surface. This specification is "the boundary value" of the solution. We do this in Rn for n = 2 and n = 3 mostly. This whole activity is very similar to finding a solution of an ODE where you insist on two boundary values, one at the endpoint of an interval (a,b), which endpoints comprise the bounding closed surface in the case n = 1. In these sections, there are no "sources", and there are no Green's functions used. Only starting in section 6.4 do we start talking about sources and Green's functions. A question remains unanswered as I started section 6.4: why is it that when you have a closed boundary, even with a p = 2 operator L, we need only the value u on the surface? The Cauchy conditions discussed in Chap 5 say we should need the normal derivative as well. For Laplace, which is elliptical, there are no characteristic tangencies to worry about. Stak looks first at the disk problem in 2D and conjectures on page 91 that u on the circle does in fact fully constrain the solution all by itself. He explicitly constructs the solution 6.11 and carefully "verifies" that it solves both posings of our problem. He derives the maximum theorem for harmonic functions and then uses it to show that the solution u is unique: 0 on entire closed boundary => 0 everywhere within. Once you have this solution you can of course compute the normal derivative on the ring and so that is fully determined and it would have been wrong to try to specify it! He pretty much avoids asking my question. But suppose we specified f(ψ) on only part of our circle. How would this change the analysis? For one thing, the constructed solution of 6.11 could no longer be constructed, so we don't have a solution, because f(ψ) is unknown on part of the arc. Looking back at the derivation, we see that the Fourier Series in complex form still has a correct expansion, but you cannot use the projection formula to get the coefficients if you have less than full circle information on u. There are as many different solutions as there are different values of f(ψ) you could supply on the missing arc! So here clearly if we have a bounded surface, knowledge of bounding u pins down the overall solution inside the boundary. But if you are missing some boundary info, then there are many solutions. I suppose in that case you could supply the normal derivative and hope to get a solution that way. One gets the feeling that some theorem like Green's or divergence is operating here and is requiring that closed boundary. Once you know the fact for the disk or the sphere, probably you know the fact for any shape just by slowly morphing from the sphere to that shape. Of course this brings up some obscure topological issues which I don't want to ponder. 6.1 Introduction (88) In physics, we have PDEs with "supplemental data" called "boundary conditions", including at the t=0 boundary in time. The subject here is that of having a "well-posed problem". The solution to the problem should exist (not too few) and should be unique (not too many). The third item is new to me: the solution should be a continuous function of the various functions involved in the PDE and the BC's. If this is missing, you probably have a stability problem with your solution. The Laplace and Poisson equation are defined, Poisson with a minus sign. The various integral theorems are stated: divergence and the two Green's identities. All old hat to the new Phil. 6.2 Interior Dirichlet Problem for the Unit Circle The phrase interior Dirichlet implies that data is specified on a closed boundary and you want to find the solution of your PDE inside that boundary. Here we are in R2 with a unit circle as boundary, and with Laplace as our PDE. The boundary values are specified by f(φ). When u satisfies 2u = 0 in some region, u is said to be harmonic in that region, a new term for me. [ I must have known this once...] At this point Stak "poses" the problem in two distinct ways. In 6.4 he wants the solution to be continuous on the closed disk. This implies that if you pick a small open disk around a boundary point, if you approach the boundary point along any path in that disk, things should be continuous. If f(φ) is itself discontinuous, then obviously the problem is not posed correctly. In this case, he says things only need be continuous if we approach a boundary point radially! That is what 6.5 says. Below that he requires that averaged over φ our u should approach f as we go to the rim. Then we avoid saying u(1,φ) = f(φ). Point of confusion: the norm defined in the text here is really "the L2 norm", but he puts a 1 subscript on it. The L1 norm involves the abs value not squared inside the sum or integral. The 1 subscript here means this is a 1D norm involving only the variable φ. Later he will have a different L2 norm that involves integration over both r and φ, and this will have a 2 subscript. Separation of Variables [ 2D, r φ ] Here we do "the usual thing" in r and φ and the separation constant is called λ. Solution is u(r,φ) = R(r) Φ(φ) // which becomes Rn(r) Φn(φ) If we want Φ to be C and C1 at π, we add page 92 A. The Φ equation just has sin,cos or expo type solutions where the eigenvalues are λ = n2, for example sin(nφ) or cos(nφ) or einφ. The eigenvalue quantization is caused by that requirement that Φ match at its π boundary. The solution to the R equation is easily found by trying a power series, and the resulting series is quite short. Solutions are R = rn and r-n . If n = 0, however, the two solutions are 1 and ln(r). Aside: Recall our Frobenius theory of ODE's. When the two roots of the indicial equation are the same, as they are here since - n = n when n = 0, the ln(r) solution appears. ___________________________ notes added later concerning general forms___________________ In the usual manner, we try u(r,φ) = R(r) Φ(φ) and identify a separation constant λ. The usual azimuthal conditions quantize this to be λn = n2 with n = 0,1,2... General solution is then lincomb of Rn(r) Φn(φ) where Φn(φ) = an eiφn + bn e-iφn and Rn(r) = (1-δn,0)) [ Anrn + Bnr-n] + δn,0 [ C + D ln(r) ]. In this section, D = 0 (ignore the log) because we are only interested in solutions inside the unit circle finite at r=0 (Theorem 1) or outside with the solution bounded at r=∞ (Theorem 2), in which case Rn(r) = Anrn + Bnr-n where C = A0+B0. So general solution is u(r,φ) = Σn=0∞[ an eiφn + bn e-iφn][ Anrn + Bnr-n ] In order to do simple Fourier series work, it is nicer if this series has the form Σn=-∞+∞. Let's try to find such a series that equals the original series shown above. Write the new series this way (just prime all coefficients) u(r,φ) = Σn=-∞+∞ [ a'n eiφn + b'n e-iφn][ A'nrn + B'nr-n ] Note inserted: Assuming we justify the above, change sign of n in the [....]Bn' term to get u(r,φ) = Σn=-∞+∞{ [ a'n eiφn + b'n e-iφn][ A'nrn] + [ a'-n e-iφn + b'-n e+iφn][ B'-nrn] } = Σn=-∞+∞ { [a'nA'n + b'-nB'-n] eiφn + [b'nA'n + a'-nB'-n] eiφn } rn = Σn=-∞+∞ {αneiφn + βn eiφn } rn // a form appearing in 6.143 later on page 175 as our trial, where the primed coefficients are defined also for negative n. If we reflect the negative side, we get this result: = [ a'0 + b'0][ A'0 + B'0 ] + Σn=1∞[ a"n eiφn + b"n e-iφn][ A"nrn + B"nr-n ] where a"n = a'n + b'-n and b"n = b'n + a'-n and similarly for A and B, all for n ≠ 0. Here we identify things term by term! This will then be the same as our original series (n ≥0, unprimed coefficients) if we can arrange to have [ a'0 + b'0][ A'0 + B'0 ] = [ a0 + b0][ A0 + B0 ] a"n = an , b"n = bn, A"n = An, B"n = Bn // for the n > 0 terms So assume we are given the second series Σn=-∞+∞ with primed coefficients. Here is one way we could construct the "exactly equal" first series Σn=0 with unprimed coefficients: a0 = a'0 b0 = b'0 A0 = A'0 B0 = B'0 an = a'n + b'-n bn = b'n + a'-n An = A'n + B'-n Bn = B'n + A'-n n>0 Suppose we had tried our second-series form with bn' = 0 for all n. Then we could say that our two forms were equal provided that a0 = a'0 b0 = 0 A0 = A'0 B0 = B'0 (*) an = a'n bn = a'-n An = A'n + B'-n Bn = B'n + A'-n n>0 So, if we select our unprimed coefficients in this manner, we have shown that Σn=0∞[ an eiφn + bn e-iφn][ Anrn + Bnr-n ] = Σn=-∞+∞ a'n eiφn [ A'nrn + B'nr-n ] So now we have two questions: (1) suppose we are handed the RHS series with all coefficients specified. How would we construct the an equivalent LHS series? Answer: follow the results stated above in (*), then we get LHS = Σn=0∞[ a'n eiφn + a'-n e-iφn][ (A'n + B'-n)rn + (B'n + A'-n)r-n ] Notice that all four sets of coefficients are completely independent and arbitrary!! (2) suppose we are handed the LHS series with all coefficients specified. How would we construct an equivalent RHS series? One way is to set all B'n coefficients to zero. Then we have a0 = a'0 b0 = 0 A0 = A'0 B0 = 0 (*) an = a'n bn = a'-n An = A'n Bn = A'-n n>0 which we can invert to get (everything here is for n> 0 only) a'n = an A'n = An B'n = 0 a'-n = bn A'-n = Bn B'-n = 0 a0' = a0 A0' = B0 B0' = 0 Let's define a "integer-argument Heaviside function" as follows: θ(n) = 1 for n = 0,1,2... but = 0 for n = -1,-2,-3... Then we can rewrite our conditions above as follows, where now n = arbitrary integer, a'n = anθ(n) + b-nθ(-n) A'n = Anθ(n) + B-nθ(-n) B'n = 0 so our RHS series becomes RHS = Σn=-∞+∞ (anθ(n) + b-nθ(-n)) eiφn [(Anθ(n) + B-nθ(-n)) rn Conclusion: the following two forms are each completely general solutions for u(r,φ) u(r,φ) = Σn=0∞[ an eiφn + bn e-iφn][ Anrn + Bnr-n ] u(r,φ) = Σn=-∞+∞ a'n eiφn [ A'nrn + B'nr-n ] Obviously, we can still have "most general solutions" in the second case if we set a'n= 1, and in the first case if we set any of the four constant sets to 1. So here are two "most general forms" u(r,φ) = Σn=0∞[ eiφn + bn e-iφn][ Anrn + Bnr-n ] u(r,φ) = Σn=-∞+∞ eiφn [ A'nrn + B'nr-n ] where we are always ignoring the ln(r) term for the reason given above. W could always add it back in. Now suppose we want an interior-of-disk solution. In the first form we set Bn = 0 for n>0. In the second case, we need to restrict it as follows u(r,φ) = Σn=-∞+∞ eiφn [ A'n θ(n)rn + B'n θ(-n)r-n ] In the second term, -n > 0 for all nonzero contributions, so then -n = |n|. And in the first term we can write that n = |n| for the same reason. So we can rewrite this as u(r,φ) = Σn=-∞+∞ eiφn [ A'n θ(n)+ B'n θ(-n)] r|n| We can then define Cn' = [ A'n θ(n)+ B'n θ(-n)] for all n to get u(r,φ) = Σn=-∞+∞ eiφn Cn' r|n| = most general form inside the disk. Now suppose we want an exterior-of-disk solution. In the first form we set An = 0 for n>0. In the second case, we need to restrict it as follows u(r,φ) = Σn=-∞+∞ eiφn [ A'n θ(-n)rn + B'n θ(n)r-n ] In the second term, n > 0 for all nonzero contributions, so then n = |n|. And in the first term we can write that n = - |n| for the same reason. So we can rewrite this as u(r,φ) = Σn=-∞+∞ eiφn [ A'n θ(-n)+ B'n θ(n)] r-|n| We can then define Dn' = [ A'n θ(-n)+ B'n θ(n)] for all n to get u(r,φ) = Σn=-∞+∞ eiφn Dn' r-|n| = most general form inside the disk. Conclusion: our most general form for u(r,φ) for the unit disk is this: u(r,φ) = Σn=-∞+∞ fn± eiφn r±|n| + for interior, – for exterior You see that both eiφn and e-iφn are in there, each with a different coefficient for a given sign, but both with the exact same power of r, for example r|n| for the interior. I had been wondering how the two independent angle functions got in there, this is how. __________________________________________________________________________ So our eigenvalue n solution to our Laplace equation in polar coordinates inside the unit circle is then simply un = an r|n|eiφn . You can then superpose these general solutions in an effort to match the boundary condition f(φ) on the ring, and so we arrive at general solution 6.8. Theorem 1 says that the solution we found is in fact a solution. He does a Fourier Series transform from φ to n, so that u(r,φ) becomes un(r) with inversion as shown. He takes our Laplace equation and processes it as in p 93 A, and after processing he shows that un(r) satisfies our separated r-only Laplace equation. I think a key idea here is that if r < 1, then the sum shown in p 93 B converges. Theorem 2 I think shows why he did theorem 1. Here we have r-|n| in the sum, and that makes us converge for r > 1. [ Both these theorems should be talking about r > 1 and r < 1, not "a", this is I think a Stak booboo. He probably had a everywhere, then changed his plan to have a = 1 until the end, but he failed to update these sections.] I accept his interior and exterior expansion representations as given in the two theorems! A nice thing is that we then have u±(1) = Σn=-∞∞ (An einφ) which must be our f(φ) BC function as p 94 A, and then we know the An as shown in B, so the coefficients An are just the fn Fourier coefficients. So we now have a solution 6.9 with coefficients 6.10. the first is an infinite sum, the second an integral. If we plug the second into the first and fiddle (I did it all) we get 6.11 which involves just one integral and no sum, so it is a simpler form. It says u(r,φ) = ∫-ππ pk(r,φ,φ')f(φ')dφ' where pk is the Poisson kernel shown. Now suppose we have f(φ) = δ(φ-φ1) as our BC. [ Don't confuse this with a delta function source, which is completely different! ] The solution is u(r,φ) = pk(r,φ,φ1). So here we have a distribution for our BC, and the solution to 2u = 0 is a regular function which happens to diverge at the one point φ = φ1 and r = 1. We then superpose the solutions with amplitude f(φ) on the BC to get our full answer. So this is an unusual superposition situation. We superpose solutions to a homo equation such that our BC is met. I could not resist plotting the Poisson kernel in Maple. I converted it to x,y and created a little marker function "f" to show where the unit circle lines (plot3D does not plot multiple functions). f := piecewise( (x^2+y^2 < .9),0,(x^2+y^2 > 1.1),0, (.9 < x^2+y^2) and ( x^2+y^2< 1.1),1); plot3d((1-x^2-y^2)/(1 + x^2+y^2 - 2*x)+f(x,y),x=-1..1,y=-1..1,numpoints=4000); A picture is worth a thousand words. Continuation Confusion: If you expand the x,y range, you get an interesting result: In the derivation, the series on the bottom of page 94 converge only for r < 1, but once you have the result, you can continue it out of this range, and thus you are continuing u(r,φ) and getting a viable solution? I thought when I first did these notes that the above really was the Laplace solution in and out and I got it by just "continuing" the formula 6.11. But it now seems pretty obvious that this is not the case since as we approach our δ BC value from r>1, we go negative instead of positive. The above is in fact a correct plot of the continuation of the kernel of 6.11, but that continuation is not in fact the solution of Laplace outside the disk! As I learned later (Exercise 6.2 below), the solution for r > 1 is the formula 6.11 with an overall minus sign. Stak verifies this explicitly in Exercise 6.4 following. That true plot is then this: Remember that we are forcing u = 0 all around the ring except at the δ location. So this is the true solution to our problem where we show both inside and out with a δ angular source. Note that this is not the same as having a δ source, this is a BV δ. Without my marker ring, you would see that the continuation and the correct solution are continuous at r = 1. But there is a dimple at the ring in the correct solution, which is a little more obvious if we cut off the delta peak , plot3d(abs(1-x^2-y^2)/(1 + x^2+y^2 - 2*x)+ 0*f(x,y),x=-2..2,y=-2..2,numpoints=4000, view=-2..4); And if we just plot our Poisson kernel radially as a function of r at some fixed θ angle we get plot( abs(1-r*r)/(1 + r*r - 2*r*.4),r=0..4); The reader at first asks: why would the potential have a discontinuous derivative like this away from the δ peak? The answer is that in effect there is "charge" all around the ring, an effective source distribution, which causes a discontinuity in the E field. But at this point in the book, we are not ready yet for this interpretation. Note also that the two "stacks" appearing above merge into one peak at r=1. Notice that the "flat area" in these pictures is u = -1, not u = 0. You can see that for large r, the PK is -1 from its formula on page 94, This is really quite an amazing formula I think. If we go to r = .9 and plot around the circle, we get this: plot((1-.95^2)/(1+.95^2 - 2*.95*cos(theta)), theta=0..2*Pi) ; showing how "flat" this is around the ring, only bumping up at φ = 0. (called it theta accidentally) We know of course that solution u at r=1 would give a δ(φ) for the above plot. Notice by the way that at any point, a "harmonic function" must have equal and opposite curvature in the two directions, Stak never mentions this visual fact. So plot is either gently curved in both directions, or highly curved in both directions, and in all cases the curvature sign is opposite. Verification of the Solution of (6.5). (95) Recall that we had two "posings" of our Laplace on unit circle problem with f(φ) on the ring. Here, Stak wants to show that our solution meets the requirements of the second posing 6.5. The first section proves the "1st condition" shown in (6.5) concerning what happens when you approach the ring in a radial direction, summed over all approaches at φ, using the L2 norm. He starts by forming a Cauchy sequence where n and m are partial sum limits. If m and n are really large, the norm shown in p 95 A becomes really small, though he does not point that out. This Cauchy sequence converges uniformly over r ≤ 1, I agree. But I don't see how this result is used. Result B is more to the point, since r = 1 means on the ring. Result C looks useful since it is close to our desired "1st condition". In fact, it looks to me that it IS our first condition, so what more is there to say? Well, C contains u(1,φ) and not f(φ). The thing we want is on the LHS of D which I suppose is trivial. Then make each term < ε/3 and OK we end up with that "1st condition". This is a lot of technical detail, but he is determined to "do it once". The second section shows that the solution u is unique. He assumes some other solution, v is the difference, then he shows v = 0, the usual method of doing this kind of thing. The third section shows something not mentioned in 6.5, but mentioned at the start of the chapter, and it has to do with stability. You want to show that the solution u has a "continuous dependence on the boundary condition f ". He imagines we have some linear operator A such that u = Af. That is what 6.11 is, for example, and A is an integral operator. To show that the mapping u = Af is continuous, we need only show it is bounded (Chap 2 showed bounded = continuous for functionals and operators in Rn). He shows this boundedness and concludes therefore that if f changes by a small amount, then u will change by a small amount, that is to say, Δu ≤ constant * Δf . If we want Δu < ε, we can find δ such that Δf < δ will achieve that goal. Verification of the Solution of (6.4). Now we shall assume f(φ) is continuous and in this case our posing is 6.4. The main goal of this section is to show that u is continuous at the boundary ring. We know it is continuous inside. The combination of these two items is the "1st condition" of (6.4). His method is to first proof a sort of lemma ("first goal") and then on the top of page 99 he shows ring continuity. His first goal here is to show that u(r,φ) → f(φ) as r → 1- uniformly for all φ. This consumes all of page 98. He does some strange shifting of the angle variable, and defines the Poisson kernel now in θ bottom of page 97. As my 1D graph above shows, he wants to show that along the ring, this PK thing approaches δ(θ) as r → 1. That would then show u(r,φ) → f(φ) just from 6.11. He wanders around with lots of tech detail on page 98 and at the top of 99 he has shown that we do have our result u(r,φ) → f(φ) as r → 1- . So the first goal is met. Now let's focus on the logic top of page 99. Result B is trivial, add and subtract and CSI. Our uniform convergence implies C. D is the continuity of f. So combine B,C and D to get A. This is all very fiddly. We are just showing that as we get radially close to the ring, if φ is close to φo, then u is close to f(φo). This is the issue of "continuity at the boundary". It is a 2D continuity thing. If φ-φ0 < δ, and if R < r < 1, then u(r,φ)-f(φo) < ε. Given any ε, we and any R, we can find δ. Given and ε and any δ, we can find R. The next tiny section on page 99 just says he is not going to worry about "continuous dependence on the BC" and "uniqueness" for this posing. The last section rescales things from a = 1 to a = a, the radius of the circle. Just replace r with r/a everywhere. The new PK then appears in 6.17 6.3 Some properties of harmonic functions A. Mean Value Theorem and the "part problem method" We have a formula for u(r,φ) for our ring BC Laplace problem. If we set r = 0, the Poisson Kernel collapses into the number 1, and the result is 6.18 which is "the mean value theorem": the solution u at the center is the average of the value on the ring. Very simple. On page 100 Stakgold gives a symmetry-based derivation of this result. I was pretty confused about this at first, but finally figured it out. Consider two "part problems" where we have f = 1 on a finite piece of the circle, and f = 0 on the rest of the circle. If we were doing a problem with "sources" (as if the arcs were regions of charge), we would be happy to say that the solution to the total problem was the sum of the solutions to the part 1 and part 2 problems. But we don't have sources, we have only BC's. We know that part 1 and part 2 have solutions for u. Explicitly, you could find these solutions by using formula (6.11). These are non-trivial solutions by the way, not just u = constant. Here is my attempt to plot a typical such solution: [ see Visio notes on plotting for more explanation of things ] readlib(addcoords)(z_cylindrical,[z,r,theta],[r*cos(theta),r*sin(theta),z]); f := int((1/(2*Pi))*(1-r^2)/(1 + r^2 - 2*r*cos(phi -psi)), psi=0..Pi/8); plot3d(f,r=0..1,phi=-Pi..Pi,coords=z_cylindrical, numpoints=4000); The first line computes the integral of 6.11 for the part problem where f = 1 on the part of the arc from 0 to Pi/8 degrees. Here is the result of that integral according to Maple And here is the plot The "wedge" on the back side is an artifact I am pretty sure which arises from the branches of the arctan function that Maple uses. The real solution is flat through that back wedge area. I have had this kind of problem before where something is pushed down from where it should be. The important part of this plot is the little "awning: area on the left side. The top of this awning is where f = 1, and it quickly drops down to f = 0 as you move off the (0,π/8) arc. In the real plot those edges are infinitely steep on the unit circle. The solution u(r,φ) is obviously non-trivial in the interior of the circle. I just wanted to plot this in the course of this discussion. Now back to our "parts" approach and the above picture which I repeat Suppose you find the solutions to the part 1 and part 2 problems. Each is a rotated version of the Maple plot I show above. If you add these two solutions, you do in fact get a solution to the "total" problem! Why? Because that sum solution meets the boundary conditions of the "total" problem and we know that the solution is unique, so we have found it! You are always allowed to superpose two solutions of the Laplace equation, but the key idea here is that when you do that here, you meet the new BC! What happens to the value at r = 0 in the above situation? Since we really are adding the two solutions on the left, we conclude that utotal(r=0) = upart1(0) + upart2(0) !! But by rotational invariance of the problem, if part 1 and part 2 have equal arc lengths, we know upart2(0) = upart1(0). Therefore we really can say utotal(r=0) = 2 upart1(0). Now put N equal arcs around the circle so that the total has f = 1 on the entire circle. Then we have utotal(r=0) = N upart1(0) But in this case, we know that utotal(r=0) = 1 because the only solution to Laplace for a constant f = 1 on the closed circle boundary is the trivial solution f = 1. So we have utotal(r=0) = N upart1(0) = 1 => uparti(0) = (1/N) i = 1...N Now suppose we "scale" the value of f separately on each section of arc, so they have values fi. We then know that uparti(0) = (1/N) (fi/1) This arises since things are linear. If we double f on this part, we double u everywhere, including at the center. When we do this, the superposition of solutions idea still works. In our example above, suppose we had f1 and f2 instead of 1 and 1 for our arc values. The sum solution then has the right BC and must be the correct solution. Therefore, the correct solution for N arcs with fi is this utotal(r=0) = Σi(1/N) (fi/1) = (1/N) Σifi and this is the essential result that u in the middle is the average of the fi values around the rim. We can now think of N = 2π/Δφ so we have u(0) = (1/2π) ΣifiΔφ → (1/2π) ∫dφ f(φ) and there is our mean value theorem, 6.18. So the MVT is true for a circular boundary, and also it is true for a regular polygonal boundary with potential a constant on each edge. The idea can be repeated in any number of dimensions, but n = 3 is good to think about. The arc is then a patch. We cover the sphere with N equal size patches. The solution additivity idea is still the same. Each step is really the same and we end up with utotal(r=0) = (1/N) Σifi where fi is the constant value on each patch of the sphere. But now we identify N = (4π)/ΔΩ so utotal(r=0) = (1/4π) Σifi ΔΩ and this then becomes (1/4π) ∫f(θ,φ) dΩ . In n dimensions, presume the solution is (1/Sn(1)) ∫f dΩn . B. The Maximum Principle for Harmonic Functions (101) Unless a Laplace solution is constant on a region R, the maximum of u over R must occur on the boundary. This is proved by contrapositive. Suppose x0 out in the middle is the maximum. Draw a circle around it. If there is even one point on the boundary of that circle x with u(x) < u(x0), then the average of u on that circle is less than u(x0) since no points on the boundary can have u larger than u(x0). This contradicts the MVT which says the average of u on the circle equals the value at xo. Therefore, there are no points on the boundary of C0 with u(x) < u(x0). We know there are no points with > u(x0) as well, so all points on the boundary must then be u(x0). Since this is true for any concentric circle smaller within C0 as well, we conclude that we have u(x0) at all points inside the circle C0 . In other works, u is constant on the disk C0 with this supposed maximum value. If we go to an overlapping circle with center at x1, we get to start all over with that circle and we conclude that u is constant on that entire circle as well. But we can cover all of the region R by a set of circles, so in this way we end up with u being constant on all of R. That is the basic idea, there are issues of R being closed and including its boundary (in my notation). THEREFORE, with the assumption that we have a max value somewhere inside, we conclude that u must be constant over R with that max value. That is fine, but if we have a NON constant solution u in R, then it cannot have a max at some interior point like x0 . QED. This same argument also applies to the minimum value!!! C. Uniqueness and BC dependence continuity of the solution Lemma: if u = 0 on a boundary of R, then u = 0 in all of R. Proof #1: If there were a non-trial solution with some values u ≠ 0 in the interior, then we would have a max or min somewhere inside R, but a non-trivial solution cannot do that, QED. Proof #2: By Green #1 which is 6.1A, he has another proof, but it requires that u exist and it requires that things are not singular on the boundary. This is a weaker proof, but he claims it is valid if you require that u be C1 in both variables. Proof of uniqueness: if there are two different non-trivial solutions as shown top p 103, then difference must be a non-trivial solution with u = 0 on the boundary. But our Lemma says that the only solution to this is u = 0, and therefore u1 = u2 and the solution (assuming a solution exists) is unique. Proof of continuous dependence on the boundary data: This is pretty clever. It is the usual ε and δ thing, but the inequality required comes from the max principle! So if you require || u2-u1|| < ε, you can find a δ which makes this true. In this case, δ = ε is what happens to work. As you get close in the range, you get close in the domain. Exercises: (103) // there are a lot of them, and I need to do at least some of these! Exercise 6.1. Compute four infinite sums. Did this on scratch paper, here is a summary. (a) Set z = α eiθ, write out the hint formula, rationalize right side so denominator is real. On the LHS write out αneinθ as real and imag parts, then compare LHS and RHS and you get the two formulas given here. (b) Integrate the hint formula, change summation index to be m = n+1. Verify that constant is zero since both sides vanish when z = 0. (c) Integrate the right formula in 6.21 wrt θ. Set x = cosθ so dx = -sinθdθ, integral is trivial. Verify no extra constant because set α = 0 and both sides give 0. As usual, log means ln in Stak. Exercise 6.2. Evaluate a Poisson Kernel integral with f = 1, then interpret the result. (2D) Just set z = eiψ and do the math and the integral shown follows after a few scratch lines. If r < 1, residue is 2πi * 1/(r - 1/r) and get J = +2π as claimed. If r > 1, it is the 1/r pole that is inside the same CCW contour, residue is 2πi * (1/r- r) and of course this then gives - 2π. Interpretation: This thing is the solution u(r,φ=0) of Laplace where f(ψ) = 1 on the entire ring. Looking at 6.11 with its Poisson kernel, if J = 2π we then get u(r,φ=0) = +1 for any r inside the circle, which of course is our known constant solution. Formula 6.11 does not apply to r > 1. But it would be easy to make a new version of it that did apply to r > 1. You would use z = 1/r instead of z = r in its derivation. The upshot is that you get the 6.11 result with r → 1/r. That gives the exact same formula but with an overall minus sign. We then get u = +1 for r>1 which again is the right constant answer. Sort of an interesting problem, much is left to the reader. Exercise 6.3. Solutions outside the unit circle, exterior Dirichlet. (2D) (a) Show u(1/r,φ) is harmonic. Well p 92 C we see the solution to the radial equation, and r → 1/r gives the exact same two solution terms, with different constants (C and D are swapped). For n = 0, you get that ln(1/r) = ln(1) - ln(r) = -ln(r) so get the same answer there with D → -D. So for r > 1, all the math on page 93 will have r → 1/r to maintain convergence of sums, we just pick "the other term" to play with. (b) Here we get the formula I anticipated in the previous exercise, showing the overall minus sign relative to 6.11. We have here the exterior Dirichlet problem with BC on the unit circle. You can think of the whole thing as a mapping of the r,φ plane where the inside of the circle becomes the outside, and the center of the circle becomes the point at ∞. With this mapping in mind, we are not surprised to find the mean value theorem quoted applying at ∞ instead of at the center of the circle. And the former limit of u to the ring from r < 1 is now the limit of v to the ring from r > 1. So don't associate the word "interior" with Dirichlet. It could be exterior. Exercise 6.4. Symmetric polygon problem. (2D) Think of this as k "part problems". We can solve each of these, and then the sum of these k solutions will match the "whole problem" BC. Here the whole problem is a regular polygon with n edges, and u = 1 on k of those edges, and u = 0 on n-k edges. Each of the k part problems has u = 1 on its edge and u=0 on all the other edges. This solution has some ue(0) at the center. If we add our k solutions, we get k ue(0) at the center. But if we had 1 on every edge, we would have u(0) = n ue(0) = 1 since this has the constant solution, so ue(0)= 1/n. Therefore in our whole problem here with the k edges, we have u(0) = k ue(0) = k * 1/n = k/n, as they claim. Exercise 6.5. Work things through for Helmholtz in place of Laplace. (2D) (105) If we go back to the ODE and separate variables, we still get 6.6 and 6.7 is replaced with this: r ∂r(r ∂rR) + k2r2R = λ R R = n2 so you see the extra k2 term sitting there. Write this out as r2 ∂2r R + r ∂r R + (k2r2- n2) R = 0 Now rescale to variable x = kr so that ∂r = k ∂x and we get (x/k)2 k2 ∂2x R + (x/k) k ∂x R + (x2- n2) R = 0 x2 ∂2x R + x ∂x R + (x2- n2) R = 0 and of course this is Bessel's equation with solutions of the type Jn(kr). Bessel is "what you get" when you worry about circular things like unit circles and drum heads and polar coordinates. (a) If the circle BC is f = 1, there is no φ dependence, so must have n = 0, so solution is C J0(kr). To normalize this to 1 at r = a, we get u = J0(kr)/ J0(ka) as claimed. So the solution here is NOT a constant as it was in Laplace world. The "symmetry argument" is to make N part problems, each makes ue(0) at the center (e = arc?) . If we had 1 on all N arcs, superposition tells us u(0) = N ue(0) . But this must be J0(0)/J0(ka) = 1/J0(ka). So we find that N ue(0) = 1/J0(ka) so that ue(0) = 1/[N J0(ka)]. In our linear part problems, we weight each arc with fi and then superpose those to get u(0) = Σi fi ue(0) = 1/[N J0(ka)] Σi fi . We then set N = 2π/Δφ so that u(0) = 1/[2πJ0(ka)] * Σi fi Δφ and this becomes u(0) = 1/[2πJ0(ka)] * ∫dφf(φ), exactly as claimed in 6.24. (b) The normalized Bessel solutions are Jn(kr)/ Jn(ka) and these play the role of (r/a)n in the math. Looking then at p 93 B with r → r/a, we get u = Σn fn Jn(kr)/ Jn(ka)einφ with full n sum. The Bessels should each have |n| as the order, but going n to -n creates (-1)n for each of the J's so these cancel, so we can then omit the abs value. The former coefficients an are just named fn here. QED. (b') We are now supposed to jam this sum into the MVT. Manana. What about drumheads? It is true that f = 1 on the circle forces you to have only the u = J0(kr)/ J0(ka) solution. But what about the case f = 0 on the circle? In this case, you get the more general solution, Σm fm Jm(kr) eimφ , but now you need km,n such that Jm(km,n a) = 0, where the n index is listing off the zeros of the function Jm(ka), So the drumhead solution (wave amplitude = 0 on the ring) has the form Σm eimφ Σn fm,n Jn(km,n r). Here are some of the drumhead modes: http://oak.ucc.nau.edu/jws8/dpgraph/drumheads.html Exercise 6.6. The 4 guy is a p = 4 operator, so Cauchy would require u and normal u', u" and u'''. But I think for a closed boundary, we only need u and u' on the boundary and that is what is specified here as f and g. (a) It is very easy to show that these particular solutions are both 4-harmonic and satisfy the conditions shown. Again a ring of f=0 and f'=1 results in the little parabolic shape shown. (b) MVT: Contribution of one of N arcs with f=1 and g=0 is ue1(0), and we assume f = g = 0 everywhere on the ring off this arc. Add up a full set of arcs and find that N ue1(0) = 1. Thus ue1(0) = 1/N. Contribution of one of N arcs with f=0 and g=1 is ue2(0), and we assume f = g = 0 everywhere on the ring off this arc. Add up a full set of arcs and find that N ue2(0) = -a/2. Thus ue2(0) = -a/2N . Now weight each arc with fi and gi and get u(0) = Σi fi ue1(0) + Σi gi ue2(0). Thus, u(0) = Σi fi ue1(0) + Σi gi ue2(0) = 1/N Σi fi -a/2N Σi gi Replace N = 2π/Δφ as usual so 1/N = Δφ/2π and get u(0) = (1/2π) { ∫dφ f(φ) - (a/2) ∫dφ g(φ) }. So the new feature here is that our "basis" of arc contributors is doubled from p = 2 examples. We have here 2N basis elements which we can weight with fi and gi to contribute to u(0). Exercise 6.7. Here L = 2 - P where P is a "positive function". (2D) What is the MVT in this situation? Since P is likely non-symmetric, we can't use our usual superposition argument. So I don't know the mean value in the center of a circle. Stak does not tell us to work in n = 2 here, so I suspect the result is valid in any number of dimensions. I will work in n=2 at least to start. Suppose we try a delta source for P and then superpose later. We know the solution for u due to this point source from chapter 5 for n = 2: u(r) = -(P1/4π) ln |r - r1| + B1. I suppose if we integrate over our region, we get u(r) = –∫dV1 { P(r1)/4π* ln|r - r1| } + B1 ∫dV1 P(r1) = –∫dV1 { P(r1)/4π* ln|r - r1| } + C where C would have either sign. So this is a Green's like solution to our problem, so at least we have something! To this we have to add an arbitrary homo solution, so we could say u(r) = –∫dV1 { P(r1)/4π* ln|r - r1| } + uh(r) where uh is an arbitrary harmonic solution and this then absorbs the constant C. Using this fact, what do we then know? If we put r on the boundary, we know this: u(rσ) = –∫dV1 {P(r1)/4π* ln|rσ - r1| } + uh(rσ) = f(rσ) An interesting way to write this is as follows: uh(rσ) = f(rσ) + ∫dV1 {P(r1)/4π* ln|rσ - r1| } ≡ F(rσ) Here, P is specified, and f is specified, so we then have a harmonic problem with the boundary value specified by the above F. Can we solve for uh(r) inside the boundary? This is the interior Dirichlet problem, and I think in theory we can solve it. So far in our reading in Chap 6, we have only learned how to do this if the boundary is a circle (sphere) and then the solution is 6.11. But even without knowing the solution, we know that uh(r) has its max and min on the boundary. Now let's consider a point rσ2 on the boundary where we know that f is positive. Then we have this fact: –∫dV1 {P(r1)/4π* ln|rσ - r1| + uh(rσ2) > 0 at the one point rσ2 It is not clear to me what I learn from this fact. For example, I see no reason why solution u(rσ) should have its max value at this same point. But this does tell us that there is at least one point on the boundary where uh(rσ2) > 0 since the integral is positive. Can the previous formula give us a MVT for this problem? At least for a disk? Let's look at the geometry. I just showed on scratch this well known fact: (n=2) |r1 - r2|2 = r12 + r22 - 2r1r2cos(θ1- θ2) So we can write our general solution a bit more explicitly like this: u(r,θ) = –∫dθ1∫r1dr1 P(r1, θ1) (1/4π) ln( [ r12 + r2 - 2r1rcos(θ1- θ)]1/2/a) + uh(r,θ) where we have added 1/a inside the ln, absorbing the constant into uh again. If we are interested in the value at r = 0, we get a simpler result: u(0,θ) = –∫dθ1∫r1dr1 P(r1, θ1) (1/4π) ln(r1/a) + uh(0,θ) = – (1/4π) ∫dθ1∫ dr1 r1 ln(r1/a) P(r1, θ1) + uh(0,θ) Let's now restrict our interest to a disk around the origin of radius a. Assume that on that disk we have some boundary condition f(θ). Then u(a,θ) = –∫dθ1∫r1dr1 P(r1, θ1) (1/4π) ln( [ r12 + a2 - 2r1acos(θ1- θ)]1/2/a) + uh(a,θ) = f(θ) The first integral is entirely determined by function P, so call this integral fP(θ) (without sign). Thus u(a,θ) = –fP(θ) + uh(a,θ) = f(θ) which tells us that uh(a,θ) = f(θ) + fP(θ) ≡ F(θ) So we are in effect GIVEN a specific boundary value for our harmonic solution, it is F(θ). The MVT for harmonic functions then tells us this: uh(0,θ) = (1/2π) ∫dψ F(ψ) and then we know the following u(0,θ) = – (1/4π) ∫ dr1 r1 ln(r1/a) ∫dθ1 P(r1, θ1) + (1/2π) ∫dψ F(ψ) So this then is our "mean value theorem" for this problem. We can express the solution at the center in terms of the boundary value f(θ) and two integrals of the given function P. Here is this result in full: u(0,θ) = (1/2π) ∫dψ f(ψ) – (1/4π) ∫dθ1∫dr1 r1 ln( [ r12 + a2 - 2r1acos(θ1- θ)]1/2/a) P(r1, θ1) – (1/4π) ∫dθ1 ∫ dr1 r1 ln(r1/a) P(r1, θ1) Again, this is our MVT for the solution to this problem, for what it's worth. Since P is positive, both integrals are negative numbers. We could just write the above as u(0,θ) = (1/2π) ∫dψ f(ψ) - A where A is some positive number that is a function of P. It took a lot of hot air to arrive at this simple result! We can now define a new function w(r,θ) = u(r,θ) + A. Then we have w(0,θ) = (1/2π) ∫dψ f(ψ) and this says that our function w has the same MVT formula as does a harmonic function, although we know that w and u are not harmonic. I think we can then apply our page 101 argument to w(r,θ) to conclude that w(r,θ) must have both its max and min values on the boundary σ . But then this must be true for u(r,θ) as well. So I arrive at a conclusion stronger than the one claimed in the text. I claim that the min and max are both on the boundary without regard to the sign of f(ψ) at any boundary point. So I have spent enough time on this one exercise, and I did not quite arrive that the correct answer, so it goes, can't do everything right, and there is no place to look (yet) for what I did wrong. One could obsess on this, but then one would not move forward in the book. What about uniqueness of the general Dirichlet problem with a source? In the above, I show that u is a particular solution plus a general homo solution. The question is then whether the homo solution part is unique or not. But it has boundary conditions F(θ) as noted above on the ring. If there were two homo solutions, the difference function would have 0 on the ring, and that means the difference function is 0 everywhere, so there can only be one solution. Exercise 6.8. Show uniqueness for Lu = (k u) . [ This L operator appears in the static heat conduction situation, see Vol I p 328 A.10. ] On page 102 we got two different proofs for the case that k = const. The first used the min/max for 2 and the second used Green. I don't know offhand whether there is a min/max for this new operator, I would have to compute its MVT first. We don't know k is symmetric. Let's try the Green's method first: Lu = 0 (k u) = 0 = k 2u + k u Try to imitate page 102. Multiply by u then integrate: 0 = ∫dV uk 2u + ∫dV u k u (*) where now we have a new second term. Consider Green #1 on the first term ∫dV (uk) 2u = ∫dS (uk) (∂k/∂n) – ∫dV (uk) u = – ∫dV (uk) u = – ∫dV [ uk + ku] u = – ∫dV u k u – ∫dV k u u Therefore, equation (*) becomes 0 = – ∫dV u k u – ∫dV k u u + ∫dV u k u = – ∫dV k u u = – ∫dV k |u|2 which is the result we want. Since the integrand is positive definite and the integral is 0, the integrand must be 0 which in this case means we must have u = 0 which means u = constant. But since u = 0 on the boundary, we must have u = 0 everywhere. So this is the new version of "lemma" on page 102. The uniqueness proof is then the same as with just 2. If two solutions, difference solution = 0, etc. QED. As for a MVT, I will let that ride for now, we were not asked for it here. Exercise 6.9. Poisson Kernel and the distributional Dirichlet problem (2D) (a) This section p 106 is just a discussion of what I already know, that we can think of the Poisson kernel as the Green's function so to speak with an f(ψ) that is a delta function in the sense of 6.11. The kernel is the solution to the Laplace equation with a delta BV on the ring where approach is from the inside. He refers to this Green's function problem as a "distributional Dirichlet problem" in the sense that f(ψ) is not a classical function but is a distribution type function, etc. (b) Now first a theorem: the φ derivative of a harmonic function is harmonic. This seems clearest from 6.8 where said derivative just provides a new set of coefficients. Or look at 6.6 and let Ψ = ∂φΦ and it is easy to show that 6.6 is then true for Ψ if it is true for Φ, so Ψ. So our Poisson kernel is harmonic, then its derivative as shown is also. Just as limr→1 sr(θ) = δ(θ), we find that limr→1 sr(θ)' = δ'(θ), where r approaches from the inside. I think this is obvious from a Chap 5 theorem that derivatives of sequences converge, etc. But here he just does out a few details. I don't see a lot of profit from this exercise. Exercise 6.10. Dirichlet Laplace on a wedge (2D) Stakgold really should have started with problems on rectangular areas, then moved to things like wedges. The rectangle appears below in a later exercise. The general idea is that in some coordinate system, we are going to have 0 at two limit curves or lines, and this is going to cause eigenvalues for THAT coordinate, and those will in turn determine the nature of the eigensolutions in both dimensions. Stak really has not said much about this in the text. Also, there is that idea I have from some unknown place that says solutions will be wavelike in one direction and expo like in the other dimension. [ see below ] So here we have a wedge lined up with φ = 0 to make it easy. We always need to go back to our separation of variables discussion on top of page 92, which I write here as -∂φ2Φ = σ2Φ r∂r(r∂rR) = σ2R σ2 = λ The eigenvalue problem here is going to be in the φ direction, so we start with that part of the problem. Our general solution is u = R(r)Φ(φ) and Φ(φ) = Aσ sin(σφ) + Bσcos(σφ), but u = 0 on φ = 0 eliminates the cos, so we have this simple form Φ(φ) = Aσ sin(σφ). We then require u = 0 at φ = α which tells us that sin(σα) = 0 so σα = πn where n is any integer. Thus, we arrive at our eigenvalues: σn = (π/α)n = πn/α. Now that we know the legal values for σ, we THEN examine the other direction r. The solution of the R problem in general is Cσrσ + Dσr-σ. However, if σ = 0 is a legal eigenvalue, this solution expression does not represent two independent solutions, our indicial equation has identical roots, we know that C0 + D0 ln(r) is then the general solution. So we can write our general R form as R = Σ'n=-∞∞ (Cnrπn/α + Dnr-πn/α) + C0 + D0 ln(r) where the sum Σ' is missing the n = 0 element. There is no loss of generality if we replace this with Σn=1∞ since the summand basis functions are symmetric under n → -n . So R(r) = Σn=1∞(Cnrπn/α + Dnr-πn/α) + C0 + D0 ln(r) Our wedge edges meet at r = 0 and at this point we require that R = 0. That knocks out the negative powers but let's keep the ln(r) term for a moment longer, so we then must have R(r) = Σn=1∞Cnrπn/α + [C0 + D0 ln(r) ] What we really mean by this is that we have a viable set of "basis functions in the r direction" Rn(r) = { Cnrπn/α for n = 1,2,3... , C0 + D0 ln(r) for n = 0 } The ln(r) term blows up at r = 0, but we need to first look at the full solution which includes Φ' Φn(φ) = An sin(πnφ/α) Our general solution is u = Σn Rn(r) Φn(r) = Σn=1∞ Cn rπn/α sin(πnφ/α) + sin(π0φ/α) [C0 + D0 ln(r) ] = Σn=1∞ Cn rπn/α sin(πnφ/α) so the constants C0 and D0 play no role in the solution. Now our only remaining task is to require that u = f on our boundary at r = b I shall call it, anticipating the quadrilateral problem yet to come where I want a to be the "inner" curve radius. So u(b,φ) = f(φ) = Σn=1∞ Cn bπn/α sin(πnφ/α) ≡ Σn=1∞ bn sin(πnφ/α) Now we want to cast this into the form of a Fourier Series as on Schaum p 131, but care is needed. The unusual situation here is that f(φ) is only "defined" on (0,α), but because it must have the form shown as the above series, it turns out f(φ) is also "defined" on (α,2α) such that the odd coefficients have opposite sign. [ See wedge meets fourier series.doc ]. The periodicity of f(φ) is in fact 2α, so we have to use the Fourier Series with L = α. When we do this, we get the above coefficient to be bn = Cn bπn/α = (1/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) n = 1,2,3... Only after some work can one show that the (α,2α) portion of this integral replicates the bottom portion in value, so we then have bn = Cn bπn/α = (2/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) n = 1,2,3... Thus, the final solution to our wedge problem can be written in this way: u(r,φ) = Σn=1∞ Cn rπn/α sin(πnφ/α) = Σn=1∞ Cn bπn/α(r/b) πn/α sin(πnφ/α) = Σn=1∞ bn (r/b) πn/α sin(πnφ/α) so our final result is this for the wedge solution for wedge of radius "b" is this u(r,φ) = Σn=1∞ bn (r/b) πn/α sin(πnφ/α) bn = (2/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) Stak's solution on page 107 uses b = 1 and then my results here agree with page 107 A and B. Comments: There is another wedge problem that we sometimes hear about: where u = 0 on the entire shaded "piece of metal" which is infinite, and we try to impose some f(φ) on a circular arc inside the wedge. I now think that the solution found above would in fact be correct for the region inside the wedge, but obviously not outside the wedge. Special case: try α = 2π so we have a full disk with a u=0 slit. Try f(φ) = sin(φ/2) which vanishes at the slit, so is consistent with the slit. Then we have u(r,φ) = Σn=1∞ bn sin (nφ/2) = f(φ) = sin(φ/2) => b1 = 1 and others 0 Solution is then u(r,φ) = rπ/α sin(φ/2) = sin(φ/2). Exercise 6.11. Dirichlet Laplace on an annulus. (2D) In this problem, the eigenvalues are in the angular direction and will be the same as they were in the full disk problem. Our general angular function is of the form Φn(φ) = Cncos(φn) + Dnsin(φn) probably general enough restricting to positive n. and n = 0 gives C0 which is fine. The radial functional form is this Rn(r) = { Anrn + Bnr-n for n = 1,2,3...; A0 + B0 ln(r) for n = 0 } So our general form solution might look like this: u(r,φ) = Σn=1∞ (Anrn + Bnr-n )( Cncos(φn) + Dnsin(φn) ) + (A0 + B0 ln(r) ) Stak does not tell us "what to do" so I will impose f(φ) on the outer boundary at r = b, and 0 on the inner boundary at r = a. The r=a one gives this: (I set C0= 1) u(a,φ) = Σn=1∞ (Anan + Bna-n )( Cncos(φn) + Dnsin(φn) ) + (A0 + B0 ln(a) ) = 0 This implies the following: Anan + Bna-n = 0 A0 + B0 ln(a) = 0 The other condition then is this: u(b,φ) = Σn=1∞ (Anbn + Bnb-n )( Cncos(φn) + Dnsin(φn) ) + (A0 + B0 ln(b) ) = f(φ) If I set L = π, this fits exactly into the Fourier Series form of Schaum p 131, 2L = 2π = periodicity. Rewrite the above like so: f(φ) = (A0 + B0 ln(b) ) + Σn=1∞ (Anbn + Bnb-n )Cn cos(φn) + Σn=1∞ (Anbn + Bnb-n )Dn sin(φn) = a0/2 + Σn=1∞ an cos(φn) + Σn=1∞ bn sin(φn) where we make these identifications a0/2 = A0 + B0 ln(b) an = (Anbn + Bnb-n )Cn bn = (Anbn + Bnb-n )Dn We know at once what these guys are in terms of integrals an = (1/π) !Syntax Error, Idφ f(φ) cos(nφ) valid also for n = 0,1,2... bn = (1/π) !Syntax Error, Idφ f(φ) sin(nφ) valid also for n = 1,2... We are now done! But there is certain simplification to be had. The first task is to use the u = 0 conditions to eliminate the Bn constants: Anan + Bna-n = 0 Bn = – An a2n A0 + B0 ln(a) = 0 B0 = -A0/ln(a) then we have a0/2 = A0 + B0 ln(b) = A01 + [-A0/ln(a) ] ln(b) = A0[ 1 - ln(b)/ln(a) ] an = (Anbn + Bnb-n )Cn = (Anbn + [– An a2n]b-n )Cn = ( bn - a2nb-n) Cn bn = (Anbn + Bnb-n )Dn = An( bn - a2nb-n) Dn which we summarize now as a0/2 = A0[ 1 - ln(b)/ln(a) ] an = An Cn( bn - a2nb-n) = bn = An Dn ( bn - a2nb-n) Now let's go back to our general form which was this u(r,φ) = Σn=1∞ (Anrn + Bnr-n )( Cncos(φn) + Dnsin(φn) ) + (A0 + B0 ln(r) ) and let's get rid of the Bn and B0 as we did above u(r,φ) = Σn=1∞ (Anrn – An a2n r-n )( Cncos(φn) + Dnsin(φn) ) + (A0 – A0 ln(r)/ln(a) ) = Σn=1∞ (rn – a2n r-n ) [ AnCncos(nφ) + AnDnsin(nφ) ] + A0 (1 - ln(r)/ln(a) ) We could at this point replace in terms of our an and bn etc to get u(r,φ) = Σn=1∞ (rn – a2n r-n )/ ( bn - a2nb-n) [ ancos(nφ) + bnsin(nφ) ] + (a0/2) [1 - ln(r)/ln(a) ] / [ 1 - ln(b)/ln(a) ] We should now check our boundary values: u(a,φ) = Σn=1∞ (an – a2n a-n )/ ( bn - a2nb-n) [ ancos(nφ) + bnsin(nφ) ] + (a0/2) [1 - ln(a)/ln(a) ] / [ 1 - ln(b)/ln(a) ] Every term now has a factor of 0, so we correctly recover u(a,φ)= 0. Out at r = b we get u(b,φ) = Σn=1∞ (bn – a2n b-n )/ ( bn - a2nb-n) [ ancos(nφ) + bnsin(nφ) ] + (a0/2) [1 - ln(b)/ln(a) ] / [ 1 - ln(b)/ln(a) ] = Σn=1∞[ ancos(nφ) + bnsin(nφ) ] + (a0/2) which is exactly our Fourier Series for f(φ), so that looks good too. Can we simplify our answer? (rn – a2n r-n )/ ( bn - a2nb-n) = rn (1 - (a/r)2n) / bn(1 - (a/b)2n) = (r/b)n (1 - (a/r)2n)/ (1 - (a/b)2n) Well, does not do much for me. Here then is the complete answer to this problem of the annulus with 0 on the inner radius a, and f(φ) on the outer radius b: u(r,φ) = Σn=1∞ (rn – a2n r-n )/ ( bn - a2nb-n) [ ancos(nφ) + bnsin(nφ) ] + (a0/2) [1 - ln(r)/ln(a) ] / [ 1 - ln(b)/ln(a) ] an = (1/π) !Syntax Error, Idφ f(φ) cos(nφ) valid for n = 0,1,2... bn = (1/π) !Syntax Error, Idφ f(φ) sin(nφ) valid for n = 1,2... What happens if we put f(φ) = 1? I think only a0 survives as a0 = 2 and our solution is: u(r,φ) = [1 - ln(r)/ln(a) ] / [ 1 - ln(b)/ln(a) ] At r = b this is indeed 1. At r = a, it is indeed 0. So it works! What happens if a > b ? Nowhere in the above work did I use the fact that I was assuming b > a (which I was assuming). So the answer the other way would be the above with a ↔ b. Ie, then we would have f(φ) specified at radius r = a instead of r = b. What happens if I want to put fa(φ) on the inner circle at a and fb(φ) on the outer circle at b? The solution will be the SUM of the two problems above, so here it is. First, define: ana = (1/π) !Syntax Error, Idφ fa(φ) cos(nφ) valid for n = 0,1,2... bna = (1/π) !Syntax Error, Idφ fa(φ) sin(nφ) valid for n = 1,2... anb = (1/π) !Syntax Error, Idφ fb(φ) cos(nφ) valid for n = 0,1,2... bnb = (1/π) !Syntax Error, Idφ fb(φ) sin(nφ) valid for n = 1,2... Then here are the two solutions we will want to add: ua(r,φ) = Σn=1∞ (rn – a2n r-n )/ ( bn - a2nb-n) [ anacos(nφ) + bnasin(nφ) ] + (a0a/2) [1 - ln(r)/ln(a) ] / [ 1 - ln(b)/ln(a) ] ub(r,φ) = Σn=1∞ (rn – b2n r-n )/ ( an - b2na-n) [ anbcos(nφ) + bnbsin(nφ) ] + (a0b/2) [1 - ln(r)/ln(b) ] / [ 1 - ln(a)/ln(b) ] I found something resembling this on the web, but it is not taken to closed form as I have done. And there is a massive use of symbols with meanings different from mine, but here it is for the record. http://books.google.com/books?id=QGw2EkI3WNIC&pg=PA7&lpg=PA7&dq=laplace+annulus&source=bl&ots=rvzNcD0nPJ&sig=tH6U4g9r57j6GcunK2_Ve9FZ50I&hl=en&ei=VhNlSqq2N4ygsgOQpIHqDg&sa=X&oi=book_result&ct=result&resnum=9 Exercise 6.12. Dirichlet Laplace on a curved quadrilateral. (2D) So here we have a combination of the wedge problem and the annulus problem, but two different BC sets are provided. In (a) we have u=0 at angles 0 and α, so this will generate the same eigenvalues as the wedge problem 6.10. In (b), however, we have u = 0 at r = a and r = b, so this is going to produce some brand new eigenvalues, this is the first problem with u = 0 at two radii. Part (a) Here we have u = 0 on the sides, as in the wedge problem. We can then crib some of our results from that problem above 6.10 R(r) = Σn=1∞(Cnrπn/α + Dnr-πn/α) + C0 + D0 ln(r) Φn(φ) = An sin(πnφ/α) // no cos due to BC u=0 at φ = 0, same in (b) to come u(r,φ) = Σn=1∞(Cnrπn/α + Dnr-πn/α) An sin(πnφ/α) + (C0 + D0 ln(r)) A0 sin(π0φ/α) Obviously we can throw out the n=0 second term due to our φ = 0 BC, combine coefficients to get u(r,φ) = Σn=1∞(Cnrπn/α + Dnr-πn/α) sin(πnφ/α) This is similar to the wedge problem, but now we have Dn coefficients that must be kept since r = 0 is not part of the "volume" of interest. Here are our two radial BC conditions u(b,φ) = Σn=1∞(Cnbπn/α + Dnb-πn/α) sin(πnφ/α) = f(φ) u(a,φ) = Σn=1∞(Cnaπn/α + Dna-πn/α) sin(πnφ/α) = 0 From the second, we conclude that Cnaπn/α + Dna-πn/α => Dn = -Cn a2πn/α which will let us dump the Dn unknowns. The other condition is this: f(φ) = Σn=1∞ bn sin(πnφ/α) where bn = (Cnbπn/α + Dnb-πn/α) Periodicity 2L of f(φ) is set by the n=1 term and we find that sin(πn[φ + 2L]/α) = sin(πnφ/α) => π2L/α = 2π L = α As in the wedge problem, our periodicity is 2L = 2α, even though data provided on (0,α). And as before, our solution is given by bn = (1/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) = (2/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) n = 1,2,3... But this time we have bn = (Cnbπn/α + Dnb-πn/α) Dn = -Cn a2πn/α so we then get bn = (Cnbπn/α + [-Cn a2πn/α]b-πn/α) = Cn(bπn/α – a2πn/α b-πn/α ) Therefore, our final solution is this: u(r,φ) = Σn=1∞(Cnrπn/α + Dnr-πn/α) sin(πnφ/α) = Σn=1∞( Cn rπn/α + [-Cn a2πn/α]r-πn/α) sin(πnφ/α) = Σn=1∞Cn (rπn/α – a2πn/α r-πn/α) sin(πnφ/α) = Σn=1∞bn (rπn/α – a2πn/α r-πn/α)/ (bπn/α – a2πn/α b-πn/α) sin(πnφ/α) As usual, we want to check our BC's. At r = 1, first factor is 0, good. At r = b, ratio is 1 and we get the expected FS of f(φ), good. Summary: u(r,φ) = Σn=1∞bn (rπn/α – a2πn/α r-πn/α)/ (bπn/α – a2πn/α b-πn/α) sin(πnφ/α) where bn = (2/α) !Syntax Error, I dφ f(φ) sin(πnφ/α) Part (b) This is different because we have f set on one of the radial edges and all else is 0. Our whole framework of part (a) does not apply because we don't have that same u(φ=0) = 0 and u(φ=α) BC situation which forced our eigenvalues σ. So maybe go back to this: u(r,φ) = Σσ { Aσ rσ + Bσ r-σ ) cos(σφ) + Σσ { Cσ rσ + Dσ r-σ )sin(σφ) + A0 + B0 ln(r) where we don't yet know the eigenvalues (I have assumed σ = 0 is one, hence the last two terms) and λ = σ2. Before continuing, we know that, regardless of the eigenvalues, only the sin terms above contribute due to our BC which says u = 0 at φ = 0. So we can simplify right away to this: u(r,φ) = Σσ { Cσ rσ + Dσ r-σ ) sin(σφ) Our eigenvalues are going to come from the radial equation, not the azimuthal equation. We have R(r) = Cσ rσ + Dσ r-σ R(a) = Cσ aσ + Dσ a-σ = 0 => Cσ bσaσ + Dσ bσa-σ = 0 R(b) = Cσ bσ + Dσ b-σ = 0 => Cσ bσaσ + Dσ aσb-σ = 0 => bσa-σ = aσb-σ => b2σ = a2σ => (b/a)2σ = 1 => e2σln(b/a) = 1 2σ ln(b/a) = 2nπ i σ = n iπ/ ln(b/a) n = 0, ± 1, ± 2, ... Let's assume then that σ takes these new eigenvalues. What can we then say about Aσ and Bσ ? Both our conditions are the same, and we have Aσ aσ = – Bσ a-σ => Bσ = – Aσ a2σ Aσ bσ = - Bσ b-σ => => Bσ = – Aσ b2σ both are consistent with our σ set So now that we know the values of σ, we have rσ = r n iπ/ ln(b/a) = eiπn ln(r)/ln(b/a) sin(σφ) = sin(n iπφ/ ln(b/a)) = i sh[n πφ/ ln(b/a)] It seems useful to re-express the rσ by including a constant (r/a)σ = (r/a) n iπ/ ln(b/a) = eiπn ln(r/a)/ln(b/a) We can of course replace the two r solutions with sin and cos and write this new version of our general form of the solution u(r,φ) = Σσ { Eσ cos[ nπln(r/a)/ln(b/a) ] + Fσ sin[ nπln(r/a)/ln(b/a) ] } sh[n πφ/ ln(b/a)] and now we see more clearly the idea that the two u = 0 BC's in radius give us sine like solutions in the r direction, and then the φ direction "goes exponential". The inner radius condition says u(a,φ) = Σn=1∞ { En cos[ nπln(a/a)/ln(b/a) ] + Fn sin[ nπln(a/a)/ln(b/a) ] } sh[n πφ/ ln(b/a)] = 0 But ln(a/a) = ln(1) = 0 so this says that Σn=1∞ En sh[n πφ/ ln(b/a)] = 0 hence En = 0. So we are left with this for our general form: u(r,φ) = Σn=1∞ Fn sin[ nπln(r/a)/ln(b/a) ] sh[ n πφ/ ln(b/a) ] Let's now define these new objects: x = ln(r/a) ln(b/a) = k so our general form becomes (this time, write the sh first) u(r,φ) = Σn=1∞ Fn sh(nπφ/k) sin(nπx/k) and our BC at r = a becomes u(a,φ) = Σn=1∞ Fn sh(nπφ/k) sin(nπ0/k) = 0 while our BC at r = b becomes u(b,φ) = Σn=1∞ Fn sh(nπφ/k) sin(nπ) = 0 so these two are verified. We can see that φ = 0 is still OK, so our final BC is this: u(r,α) = Σn=1∞ Fn sh(nπα/k) sin(nπx/k) = f(r) // x = ln(r/a) , (r/a) = ex , r = aex We want as usual to interpret this as a Fourier Series Schaum p 131. First question: what is the periodicity of this function of x? Δx π/k = 2π so that Δx = 2L = 2k = periodicity and k = L. So we now write the above as F(x) ≡ f(r = aex) = Σn=1∞ bn sin(nπx/k) bn = Fn sh(nπα/k) bn = (1/k) !Syntax Error, Idx F(x) sin(nπx/k) = (2/k) !Syntax Error, Idx F(x) sin(nπx/k) This situation is the same as that encountered in the wedge exercise 6.10 where α → k, and this I think then justifies the second form shown on the right for bn. So here then is the solution to our part (b) problem: u(r,φ) = Σn=1∞ bn [sh(nπφ/k)/ sh(nπα/k) ]sin(nπ ln(r/a)/k) k = ln(b/a) bn = (2/k) !Syntax Error, Idx f(aex) sin(nπx/k) One more set of checks: (1) at φ = 0 get u = 0, check (2) at r = a get u = 0, check (3) at r = b get u = 0, check (4) at φ = α get u = Σn=1∞ bn sin(nπx/k) = F(x), check. I think I have done this correctly, but I could have made some error in the "general form". I don't see how to do part (b) in some other way as Stak hints page 108 top. Notice by the say that all our solutions have a certain structural form. Notice that the above solution is extremely similar to that of problem 6.13 yet to come which is in Cartesians. I am sure I am on the right track here, but cannot find a annular wedge solution on the web. Exercise 6.13. Dirichlet Laplace on a rectangle. (2D) (108) The Cartesian story is this: (∂12 + ∂22)u(1,2) = 0 (∂12 + ∂22)f(1)g(2) = 0 ∂12 f(1)g(2) = – ∂22f(1)g(2) g(2) ∂12 f(1) = - f(1) ∂22 g(2) ∂12 f(1) / f(1) = – ∂22 g(2) / g(2) = – λ = – σ2 the separation constant The separate ODE's are then ∂12 f(1) = – σ2 f(1) ∂22 g(2) = + σ2 g(2) which shows the famous "minus sign between the directions". Now, in our problem, we have u = 0 at x1 = 0 and a so we want sine in the 1 direction so that σ2 > 0 and σ = real. Our solution is then f(x1) = Aσ sin(σx1) + Bσ cos(σx1) but if σ = 0 get f(x1) = A0 + B0 x1 If we want f(0) = 0, we keep only sine terms and we can set A0 = 0 as well, so f(x1) = Aσ sin(σx1) , B0 x1 for σ = 0 If we want f(a) = 0, we learn the eigenvalues to be σ = nπ/a and we learn that B0= 0 fn(x1) = Ansin(nπx1/a) , 0 for n = 0 Meanwhile, in the other direction we know we have g(x2) = Cσ sinh(σx2) + Dσ cosh(σx2) but if σ = 0 get g(x2) = C0 + D0 x2 But requiring that g(0) = 0 gets rid of the C0 so we get g(x2) = Cσ sinh(σx2) , D0x2 for n=0 Our most general form is then this ( since A0 + B0 x1 = 0 ) u(x1, x2) = Σn=1∞ Ansin(nπx1/a) Cn sinh(nπx2/a) + (0) D0x2 So let's cut down on the constants a bit and say u(x1, x2) = Σn=1∞ An sinh(nπx2/a) sin(nπx1/a) where I have put the sine second in the formula. Now at x2 = b we get u(x1, b) = Σn=1∞ An sinh(nπb/a) sin(nπx1/a) = f(x1) This has the usual Fourier Series form in x1. The periodicity is πΔx1/a = 2π so Δx1 = 2L = 2a and L = a. Write the above as u(x1, b) = Σn=1∞ bn sin(nπx1/a) bn = An sinh(nπb/a) bn = (1/a) !Syntax Error, Idx1 f(x1) sin(nπx1/a) Our general result is this: u(x1, x2) = Σn=1∞ bn [ sinh(nπx2/a)/ sinh(nπb/a)] sin(nπx1/a) bn = (1/a) !Syntax Error, Idx1 f(x1) sin(nπx1/a) = (2/a) !Syntax Error, Idx1 f(x1) sin(nπx1/a) As we learned back in problem 6.10, the upper half of the integral replicates the lower. He then asks us to examine how the solution u(x1, x2) approaches f(x1) along the top of the rectangle. Basically x1 is just a parameter in the limit we look at, x2 → b. I think then this is just standard Fourier Series convergence. We have u(x1, x2) - f(x1) = Σn=1∞ bn [ sinh(nπx2/a)/ sinh(nπb/a)] sin(nπx1/a) - Σn=1∞ bn sin(nπx1/a) = Σn=1∞ bn sin(nπx1/a) { sinh(nπx2/a)/ sinh(nπb/a) - 1 } = Σn=1∞ bn sin(nπx1/a)/ sinh(nπb/a) { sinh(nπx2/a) - sinh(nπb/a) } = Σn=1∞ kn { sinh(nπx2/a) - sinh(nπb/a) } Then || u(x1, x2) - f(x1)|| ≤ Σn=1∞ | kn| | sinh(nπx2/a) - sinh(nπb/a)| → 0 Exercise 6.14. Complete orthogonal set Here we basically have expansion of a function f(x) on orthonormal set φn(x), f(x) = Σnfnφn(x) fn = !Syntax Error, Idx φn(x)* f(x) with a little twist added. For each term we add as a factor a function of another variable y which we want to consider as y → 1, and we have vn(y) → 1 in this limit. Then we want to show this: f(x) = limy→1 {Σnfn vn(y)φn(x)} in the L2 norm sense of convergence This is one of those little "order interchange" theorems. One limitation is that the vn(y) have to be smooth and bounded in the neighborhood of y = 1 with the same bound for all n (uniform over n). I won't prove this, but I know I could. No doubt he will use this little theorem soon. ___________________________________________________________________________ Comment: I have been preferring to do things in sin and cos rather than einφ because it makes things that are real look real. Looking at p 94 we know that, from B, an = projection onto e-inφ an* = projection onto e+inφ an* = a-n Then in the sum of A we have aneinφ + a-ne-inφ = aneinφ + an*e-inφ = 2 Re [aneinφ] = real If you go the sine and cosine route for the disk, what happens in angle world? You have both Φ(φ) = Σn=1∞ ( Ancos(nφ) + Bnsin(nφ) ) and nothing kills either of these terms, as in the wedge case. So our general solution here is: u(r,φ) = Σn=1∞ ( Ancos(nφ) + Bnsin(nφ) ) rn + A0 (C0+ D0 ln(r)) and we can toss the ln(r) term I suppose if we want a non-singular u at r = 0. Then we have u(r,φ) = Σn=1∞ ( Ancos(nφ) + Bnsin(nφ) ) rn + K Now we impose our condition that u(a,φ) = Σn=1∞ ( Ancos(nφ) + Bnsin(nφ) ) an + K = f(φ) and we have a match to our Fourier Series in the real form. Periodicity 2L = 2π and L = π so An an = an = (1/π) !Syntax Error, Idφ f(φ) cos(nφ) n = 0,1,2,3 Bn an = bn = (1/π) !Syntax Error, Idφ f(φ) sin(nφ) n = 1,2,3 K = a0/2 where a0 = (1/π) !Syntax Error, Idφ f(φ) K = (1/2π) !Syntax Error, Idφ f(φ) = mean value So here we recover the fact that u(r=0) is the mean value, and here is our complete solution u(r,φ) = Σn=1∞ ( ancos(nφ) + bnsin(nφ) ) (r/a)n + a0/2 an = (1/π) !Syntax Error, Idφ f(φ) cos(nφ) a0/2 = (1/2π) !Syntax Error, Idφ f(φ) bn = (1/π) !Syntax Error, Idφ f(φ) sin(nφ) Schaum comments on p 131 how the real and complex form of the Fourier Series Transform are related, as I have commented above. So in event, for the full disk, both the sin and cos contribute because there is nothing to knock out one or the other as on the wedge with a 0 edge. ___________________________________________________________________________ Exercise 6.15. Laplace solution inside a sphere whose surface has f(θ,φ) (3D) Now we are ready to do this thing in 3D. We are going to impose f(θ,φ) on the outside. When we separate in 3D we get 6.25 which is well known to me, with coefficients as in p 109 A. He then does a pretty obscure series of steps! First, he sets θ = 0 so that z = cosθ = 1, and quotes some assoc. Legendre properties which say that on the polar axis, only m = 0 terms contribute to the sum which is pretty reasonable. So if we keep only m = 0 in 6.25 and install a0m looking up N0n, we do indeed get result C. We can then do the sum on n using a pretty obscure fact shown in A.14 which I will do later below. Now the question is how do you get from D to E ? This has caused me confusion now and in the past, and you can be sure that every Stak reader winced at this point on page 109. There are several ideas involved here, each with a point or two of potential confusion: 1. If we have a full 4π integration region, we know that an integral I = ∫ dΩ f(Ω) will have the same value no matter what coordinate system we use for its evaluation. Suppose Ω' = RΩ. Then the integral is I = ∫ dΩ' f(Ω') . Notice that this would not be true if we had some specific solid angle region of integration! That is because the ∫ 's would not have the same endpoints. Here is a 1D example of this idea: I = ∫dx f(x) = ∫dx'f(x') x = x' + 5 These integrals are only equal if they are over the full real line, the entire "group space" so to speak. The endpoints in both integrals will be the same. The range of both integrations is the same. 2. Imagine that we are given a general equation of this form, g(Ω1) = ∫dΩ' F(Ω', Ω1) where we have a full 4π integration and where Ω1 is a "variable" θ1,φ1 . If we look at the above equation in a "rotated coordinate system" we will see this: g(RΩ1) = ∫dΩ" F(Ω", RΩ1) where Ω" is the rotated version of Ω', so we could say that Ω" = RΩ'. But since the integration is full 4π range, we could replace it with the original Ω' integration and say g(RΩ1) = ∫dΩ' F(Ω', RΩ1) Since we can get to any direction Ω2 by a rotation from Ω1 so that Ω2 = RΩ1, we have shown that g(Ω2) = ∫dΩ' F(Ω', Ω2) Therefore, we have shown that this equation must be true for ANY direction Ω2 if it is true for some one particular direction Ω1. As a specific example, let Ω = R Ω0 where Ω0 = (0,-) and Ω = (θ,φ) so we in fact know exactly what R would be, but no need to write it down. Then these two facts would be true: g(Ω0) = ∫dΩ' F(Ω', Ω0) Ω0 = (0,-) g(Ω) = ∫dΩ' F(Ω', Ω) Ω = (θ,φ) 3. Now suppose we look at the special case of the above where F(Ω', Ω2) = f(Ω') h(' 2) Then our two equations above would become: (where r0 =(r,0,-) ) g(Ω0) = ∫dΩ' f(Ω') h(' 0) r' = (r,θ',φ') ' 0 = cosθ' Ω0 = (0,-) g(Ω) = ∫dΩ' f(Ω') h(' ) r = (r,θ,φ) ' ≡ cosγ Ω = (θ,φ) This then explains how we get from D to E! I like "my words" much better than "his words". Then we write: = (sinθcosφ, sinθsinφ, cosθ) ' = (sinθ'cosφ', sinθ'sinφ',cosθ') ' = sinθcosφ sinθ'cosφ' + sinθsinφ sinθ'sinφ' + cosθ cosθ' = cosγ = cosθ cosθ' + sinθ sinθ' { cosφ cosφ' + sinφ sinφ' } = cosθ cosθ' + sinθ sinθ' cos(φ-φ') // as claimed in p 110 A Sphere of radius a. If we change the radius of the sphere from 1 to a, then we modify 6.25 to have (r/a)n, and then equation p 109 A is the same. Exercise 6.16. Laplace solution outside a sphere whose surface has f(θ,φ) (3D) Let's just rework page 109 for r > 1. Then 6.25 would have r-n-1 because these are all the solutions that are 0 at infinity, one of our added BC's. Thus, we have already confirmed p 110 B. Result p 110 C then applies at once because the power of r plays no role in this. But what happens to 6.26? Let's go through the p 109 work. We would still find that only m = 0 contributes so C would be true, but with r-n-1. At this point, we would need the r > 1 version of A.14 and my conjecture is that you replace rn by –r-n-1 (see notes above), so we then get D with an overall minus sign. Then E is still true but with a minus sign, and then 6.26 is true as well for r > 1 but with a minus sign. We could apply our same rule to 6.26 of r → 1/r to reach this same conclusion. And this is the same conclusion I reached for 6.11 in the 2D case. Exercise 6.17. Laplace solution for spherical annulus (3D) Solution general form is now this: u(r,θ,φ) = Σn=0 Σall m (Anm rn + Bnmr-n-1) Ynm(θ,φ) The two boundary conditions are now u(b,θ,φ) = Σn=0 Σall m (Anm bn + Bnmb-n-1) Ynm(θ,φ) = fb(θ,φ) u(a,θ,φ) = Σn=0 Σall m (Anm an + Bnma-n-1) Ynm(θ,φ) = fa(θ,φ) Rewrite these as u(b,θ,φ) = Σn=0 Σall m bnm Ynm(θ,φ) = fb(θ,φ) bnm = Anm bn + Bnmb-n-1 u(a,θ,φ) = Σn=0 Σall m anm Ynm(θ,φ) = fa(θ,φ) anm = Anm an + Bnma-n-1 Apply page 395 E orthogonality with Nnm to get bnm = (1/Nnm)∫dΩ Ynm(Ω)* fb(θ,φ) anm = (1/Nnm)∫dΩ Ynm(Ω)* fa(θ,φ) Then reverse solve for Amn and Bmn and you are done. bnm = Anm bn + Bnmb-n-1 anm = Anm an + Bnma-n-1 anbnm = Anm anbn + Bnm an b-n-1 bnanm = Anm anbn + Bnm a-n-1bn anbnm – bnanm = Bnm(an b-n-1– a-n-1bn) Bnm = (anbnm – bnanm) / (an b-n-1– a-n-1bn) bnm = Anm bn + Bnmb-n-1 anm = Anm an + Bnma-n-1 a-n-1bnm = Anm a-n-1bn + Bnm a-n-1b-n-1 b-n-1anm = Anm an b-n-1 + Bnma-n-1 b-n-1 a-n-1bnm - b-n-1anm = Anm(a-n-1bn - an b-n-1) Anm = (a-n-1bnm - b-n-1anm) / (a-n-1bn - an b-n-1) Bnm = (anbnm – bnanm) / (an b-n-1– a-n-1bn) // from above = – (anbnm – bnanm) / (a-n-1bn - an b-n-1) So we end up with Anm = + (a-n-1bnm – b-n-1anm) / (a-n-1bn - an b-n-1) Bnm = – (anbnm – bnanm) / (a-n-1bn - an b-n-1) Then our solution is this: u(r,θ,φ) = Σn=0 Σall m [ (a-n-1bnm – b-n-1anm) rn – (anbnm – bnanm)r-n-1 ] / (a-n-1bn - an b-n-1) Ynm(θ,φ) bnm = (1/Nnm)∫dΩ Ynm(Ω)* fb(θ,φ) anm = (1/Nnm)∫dΩ Ynm(Ω)* fa(θ,φ) Limits? The limit a → 0 is ill-defined because in the limit you will have u(0) = MVT and it cannot just be specified by some fa. In this limit, u(0) for sphere is function of fb, but in our limit here it would instead be a function of fa averaged. So the limit of this problem does not give the sphere situation! But maybe we can take the a→ 0 limit in this sense: assume that fa and fb are the same order of magnitude so same is then true of bnm and anm. Then for some large r away from the center, we approximate our solution above in this way: u(r,θ,φ) = Σn=0 Σall m [ (a-n-1bnm) rn / (a-n-1bn) Ynm(θ,φ) = Σn=0 Σall m bnm (r/b)n Ynm(θ,φ) bnm = (1/Nnm)∫dΩ Ynm(Ω)* fb(θ,φ) This then replicates the page 109 top results for just a sphere of radius b. So our spherical annulus meets this limit OK. I think I did it right. Another comment: I don't see this as a superposition problem say of exterior of sphere a and interior of sphere b. The simple reason is that for interior of sphere b, for example, we don't have u = 0 on the sphere at r = a, which would be necessary to superpose the two BC's. General Comment: Cauchy, Dirichlet and Neumann BC's Recall the general Cauchy BC of the previous chapter, and apply to Laplace 2D. We have an elliptic equation, and we can pick a curve C which avoids "problems" (ie, avoids being tangent to either of the two characteristic curve families). If we specify u and all normal derivatives up to p-1 = 1 on this curve, the solution is very likely fully determined. So in our case here we have u and ∂u/∂n as our Cauchy data. The problem is determined if C is any finite curve! It does not have to be a closed boundary of your volume of interest. Now here we have been specifying only u on our boundary, not ∂u/∂n. It appears that if we specify u on a closed boundary like the unit circle, then the solution is fully determined and no Cauchy data need be (nor can be) supplied for ∂u/∂n. Probably a property of a closed boundary is that it will be tangent to characteristic curves and some set of points (a circle would have two points for each curve family) and that somehow disqualifies the Cauchy BC approach. That would explain why trying to specify ∂u/∂n would overspecify things and lead to contradictions for ∂u/∂n. [ but elliptic has no char curves! ] The word Dirichlet is used when you specify just u on some closed boundary, like our f(ψ) for the unit circle problem. We will probably learn later that you can instead specify only ∂u/∂n and that is called the Neumann BC case. And you can have "mixed" BC's in some situations we will no doubt see. Restatement (from meta notes): his is a subject that Stak sidesteps. Back in Chapter 5, we learned that if you took a piece of hypersurface C (boundary), and if you are dealing with a 2nd order PDE operator L, then if you specify u and ∂nu on that piece of boundary, it is very likely that your entire solution is determined. There were potential problems if your surface C was tangent to the characteristic curves of L, but Laplace is an elliptic equation and there are no such curves to worry about. Now in Chapter 6 we have learned that if we have a closed hypersurface or boundary, then specifying just u on the entire surface determines the solution (Dirichlet problem). If you tried to also specify ∂nu, your problem would be overspecified and there would be no solution. So somehow the "general theory" of Chapter 5 does not apply to bounding surfaces which are closed. We do know that if such a surface is not closed, then u on that non-closed bounding surface does NOT give you the complete solution, since you can fill in the missing piece with any function you want and for each there will be a different solution. Comment on Theory so Far versus Real World Problems 1. In electrostatics, we often have metal conductors floating around in space. Such conductors have surfaces which are at a constant potential, the reason being that were there any transverse E fields causing potential differences, currents would flow. So "conductor electrostatics" offers no problems of the Dirichlet type other than f=constant. However, the charge distributions which arise in conductor problems cause normal E fields which means cause normal gradients at surfaces. This offers itself to Neumann type BC problems. We could consider electrostatics on insulating objects which can hold charge stably on their surfaces. In this situation, you could have transverse fields and thus you could have some f(θ,φ) function on a surface. You would have to construct it by micro-loading charge at each point on the surface to create the desired potential f on the surface. In this case there would of course be normal gradients determined by f, so you could not do "Cauchy" assignment of the data, and this seems to go with closed surfaces. In both the above examples, we have charges (sources) active. Still, away from the surfaces, we do have the Laplace equation. 2. Steady state heat flow solutions result in stable Laplacian temperature distributions. So you could micro-control a temperature f(θ,φ) at each point on a surface and wait for steady state to arrive. The temperature distribution in a uniform conducting solid object contained by said surface would then be the solution to the kind of Dirichlet Laplace problems we have been considering. 3. In terms of diffusion, you could micro-stabilize a density of some molecules on a surface (that does not attract or repel molecules) and then I suppose in the vacuum space contained by the surface, you would obtain a density n determined by Laplace. Again I think we are talking steady state here. 6.4 Surface Layers (110) Now we are done for the moment with our "boundary value" concerns in the absence of sources, and here we start looking at sources and the corresponding fundy solutions which I just call Green's functions incorrectly. Just as a reminder, here is how you show that G = 1/4πr for n=3 Laplace -δ(r) = 2E : -1 = ∫dV (E) = ∫dS E = ∫(∂E/∂r) dS = ∫(∂E/∂r) rn-1dΩn = Sn(1) (∂E/∂r) rn-1 For n = 3, S3(1) = 4π and E = k/r so ∂E/∂r = -k/r2 and then Sn(1) (∂E/∂r) rn-1 = -4π k and k = 1/4π. This is a long section and I will try to divide it into its various topics. In the text, pencil lines indicate the start and end of topics. A1. Finding u from a source density function q (n=3) (111, 6.30 resulting u) The Laplace is now –2u = q(x) [ called Poisson ] where q(x) is a 3D source density. We can think of u = electrostatic potential and q = charge, if we like. In any event, the solution we know is this u(x) = ∫dξ G(x-ξ) q(ξ) G(x-ξ) = 1/4π|x-ξ| If we look back at Chapter 5, we see that we have a convolution integral here and u = q*G. Writing this out we have u(x) = ∫dξ q(ξ) G(x-ξ) u = q*G dξ means d3ξ in our application 2u = ∫dξ q(ξ) 2G(x-ξ) 2u = q*2G But then q*2G = -q*δ = -q . So using the convolution theory for distributions we start with 2G = -δ and we then prove that 2u = - q, so his point is that everything is rigorous in the sense of distribution theory. That is, it won't matter if q(x) has some delta function aspect to it and is only a symbolic function. B1. Suppose the source density function represents an isolated unit dipole. (112, 6.29) The dipole is two charges of size 1/h separated by distance h. Recall that p = qh and so here p = 1, a unit dipole. Our equation is then -2u = ∂l δ(x-ξ) where l is associated with the unit vector, which is the direction of the dipole, so p = . You can see at once that in this example, q(x) is a distribution (as it is of course in the Green's definition). So what is the solution of -2u = ∂l δ(x-ξ) ? First, in classical physics this would be "the potential of a dipole" which problem has this solution: which we would write here with k = 1/4π and then u(r,θ) = cosθ/4πr2. Stak now solves the equation -2u = ∂l δ(x-ξ) and comes up with this same answer. Here are the details in the Stak notation: usol(x|ξ) = limit h→ 0 of the two single delta charge potentials = Dl(x|ξ) // notation! = limh→0 (1/4π) [ (1/h)/|x-(ξ+h)| - (1/h) /|x-ξ| ] = (1/4π) limh→0 { [ f(ξ+h) - f(ξ) ]/h } = (1/4π) ξf(ξ) = (1/4π) ∂l f(ξ) But we can go compute ξf(ξ) = ξ (1/|x-ξ|) = - r(1/r) where r = x - ξ = /r2 // by trivial calculation using spherical gradient formula = (x ξ)/ |x-ξ|2 So therefore we get ξ (1/|x-ξ|) = (x ξ)/ |x-ξ|2 and Dl(x|ξ) = usol(x|ξ) = (1/4π)∂l f(ξ) = (1/4π) ξf(ξ) = (1/4π) (x ξ)/ |x-ξ|2 = (1/4π)cos(x ξ, )/ |x-ξ|2 = cosθ/4πr2 // NB !!! in terms of the above picture. The dipole is located at position ξ and θ is the angle between the dipole axis and the direction out to the observation point x. In the web solution above, we limit r→ ∞, which is the same as d → 0, because d << r is what we need. In the Stak solution our limit causes the nice form you see above. Note added: Above we have shown that Dl(x|ξ) = (1/4π) (x ξ)/ |x-ξ|2 = (1/4π)cos(x ξ, )/ |x-ξ|2 where our dipole is at location ξ and where is the normal at ξ. In terms of what lies ahead, this normal at location ξ is called n, and we then have Dn(x|ξ) = ∂n E(x|ξ) = (1/4π)cos(x ξ, )/ |x-ξ|2 = k(x,ξ) as on page 119 6.45 where ∂n means ξ and so is associated with the ξ coordinate, not the x coordinate. But suppose we put our dipole at x (with normal ν) and let ξ be some other point. Then we can say Dν(ξ|x) = ∂ν E(ξ|x) = (1/4π)cos(ξ x, )/ |x-ξ|2 = k(ξ,x) as on page 119 B where ∂ν means x and so is associated with the x coordinate, not the ξ coordinate. Of course we can always write E(ξ|x) = E(x|ξ) wherever it appears, since E is symmetric. Application of this note: Look at page 172 C. According to the note above, ∂ν E(x|ξ) = k(ξ,x) A2. Suppose our charge q(x) is spread on a surface σ (simple layer) (112, result 6.30 for u) The result is extremely simple. We just superpose the potential of charges q(ξ) dSξ : u(x) = ∫σ dSξ q(ξ) * 1/4π|x-ξ| n = 3 only // ∫dS q 1/4πR and this is known as a simple layer. Here of course q has dimensions of charge per area. A3. Proof that u of A2 is finite for x actually on the surface (112) The method of proof uses a tangent plane at the point s on σ we approach with x. Our only possible concern is that we might get a divergence from the charge immediately surrounding the point s, since we know that u ~ q/r, it is reasonable to expect possible trouble here. Stak is careful to do things in n dimensions here, but in the 3D case consider σε = (curved) disk around the point s on the tangent plane. As this disk with radius ε is shrunk to 0 (disk becomes a disk), he shows that limε→0 ∫σε dSξ q(ξ) * 1/4π|s-ξ| = 0, which then shows that this dangerous local part of the surface is not causing a divergence in our solution u. His proof is interesting. This picture shows how secγ appears, he refers to the tangent plane area dSt as just dS', so dSξ = dS' secγ As ξ moves closer to s, and assuming σ is a smooth surface, γ → 0 so cosγ gets closer to 1. So it is pretty reasonable to imagine we first get close enough so that sec γ < 1.6 or any number you like, call it 2. So this gets us to result p 112 C. We are still worried that the thing on the RHS could diverge! Polar coordinates then gets you to p 113 A. We are then faced with simply ∫dρ = ρ which of course goes to zero in our limit, QED. It is really the fact that the surface σ has one lower dimensionality that makes this work. A close in-charge is q(ξ) dSξ . As you zoom in, dSξ gets smaller and → 0 so the effective charge is going to zero, and even if you add up a ring of these around s, you still get 0. B2. A dipole layer on surface σ (113 result 6.31 for v) I might have called the density p(ξ), but he uses b(ξ). Again, we just superpose our single dipole solutions and we then have an integral over σ. The angle of the dipole on the surface is assumed to be pointing out normal from the surface, that is what is meant by a "double layer". So then replace = and we have our result as shown in 6.31. The angle is then between the normal at point ξ on σ and the ray from ξ to our observation point x. In other words, just the usual polar angle of the dipole. u(r) = ∫σ dS' q(r') 1/4πR R = r - r' r' is on σ simple layer v(r) = ∫σ dS' b(r')cos(R,n')/4πR2 R = r - r' r' is on σ double layer B3. Proof that u of B2 is finite for x actually on the surface. (113) He just outlines the proof here. For a little disk on the tangent surface at point s, containing our point ξ of interest, we have roughly cosθ = 0 since 90 degrees. But σ is not flat, so really cosθ = constant * |s-ξ| , but then this cancels only one of the denominator factors and we have exactly the same situation as in A3 with the same conclusion. C. Detail about x → s where s is on the surface σ which has a δ source. (113-115) This section is "setup work" for section D which follows. Stak is "modeling" a point s on a surface σ as if that point s were the origin, and as if the tangent plane at s were in the x-y plane and the normal were z, which is he will call x3 here. (1) Simple Layer (113, result p 114 D). Stak imagines a surface = the x1-x2 plane and we approach by letting x3 → 0. We assume we have a symbolic function charge δ(r) at the origin, so u = G = E(x1, x2, x3). On page 114 we have a pointwise function claim for E which is obviously true, where E is the usual 1/4πr Green's. He then wraps this with an integral over the planar surface against a 2D test function to get 6.32 where φ is K2. So now we have result D being true in a distributional sense. This is important because when x1 and x2 both go to 0, this limit blows up in the limit x3 = 0 as well. But this "blowing up" is handled by the theory of distributions. The next section applies to both the "dipole layer" and to the normal derivative of u from "simple layer". (2) Double Layer (114, results 6.33 and 6.34). Recall result top of page 112 which says Dl(x|ξ) = ∂l [ 1/4π |ξ-x| ] = ∂lE. This says that the potential of a Green's dipole is given by the derivative (in the dipole direction) of the Green's monopole. With our coordinates as above, this says Dx3(x|ξ) = ∂3(1/4πr). But ∂3r = x3/r and then ∂3 (1/r) = - x3/r3 and so Dx3(x|ξ) = - x3/4πr3 which he defines as F(x1, x2, x3). He then wraps this with a 2D test function integral in F. Now he wants to show that limx3→0 { x3/4πr3 } = (1/2) δ(x1) δ(x2) because x3/2πr3 is a delta sequence in this sense with parameter x3 . He does in fact use Exercise 5.4 which just says this: if full integral gives a 1, and any integral outside a shrinking core gives 0, then everything must be coming from the center of the core. So limx3→0 (x3/2πr3) = δ(x1) δ(x2), an interesting fact all on its own! Result 6.34 then follows at once if you just integrate this symbolic function equality against a test function. If we approach from x3 < 0, our full integral is -1 instead of +1 in G and that causes the sign change in 6.33 and in 6.34 for this case. Physically what is happening here with this dipole situation? The process of interest is that as x3→ 0, we might hit anywhere in a little disk in the x1x2 plane near the origin. If we miss the actual origin, then our dipole potential is 0 because x3/r = cosψ = 0 and nothing is blowing up to offset this 0. However, if we hit the origin, then we know that our potential must be singular, since v = x3/4πr3 and we have r → 0 as well as x3→ 0. This problem only happens exactly when x1 = 0 and x2 = 0, hence we get something proportional to δ(x1) δ(x2), and the proportionality constant happens to be 1/2. Later when we spread a dipole layer b(r) over this origin region and integrate on a tiny disk around the critical point, we are going to pick up 1/2 b(s) where s is the origin point of our approach. The rest of the surface contributes nothing if it is exactly flat, since x3/r = cosψ = 0. We can now summarize these last two little sections: ( E = 1/4πr = Green's ) (1) limx3 → 0± { E(x1, x2, x3) } = E(x1, x2, 0) in full symbolic function sense. (2) limx3 → 0± { ∂3 E(x1, x2, x3)} = ∓δ(x1) δ(x2)/2 in full symbolic function sense We think of σ as represented here by the x1, x2 plane and we are taking a limit as we approach the surface near r = 0. We assume in (1) that we have a unit delta source at the origin, and in (2) that we have a unit dipole source at the origin. Notice that both results diverge at the origin, the second more violently. (3) (mid page 115). Here is shows that all tangential derivatives of our Green's E are well defined in the distributional sense as we approach the surface σ. Above we treated only E itself and its first normal derivative. Nothing singular happens with these tangential multiindex derivatives because test functions are smooth. But why couldn't you apply this same argument to ∂3E ? I think for the transverse cases, x3 is a parameter, but not for the x3 case, so cannot take parameter limit in that case. D. Extend the Part C analysis from δ sources to layers on the x-y plane (115) Integration on the plane σ always uses variables ξ1 and ξ2. We just use the layer superposition versions of our δ results above. You can think of ξ12 + ξ22 as ρ2 on σ surface. We are going to approach the surface σ along the normal direction which here is x3 (1) Simple Layer u and ∂3u. ( 116 results 6.36 and 6.37) The limit shown in 6.36 is one we really considered much earlier, our point goes onto the surface, and we claim that the integral shown here converges even though our point x gets onto the surface. He uses the symbol u for a function driven by a simple layer. The limit shown in 6.37 is for the normal gradient of the solution u at σ. Here we cannot just take the limit x3→0 through the integral! The reason is that the integrand is very singular in this limit, even as an integral, because the denominator goes to 0 in the limit for part of the integration. BUT, we learned back in 6.33 that we can handle just this situation if we think of the limit of -x3/4πr3 as - δ(x1) δ(x2)/2, rather than being the naive limit 0 that you would at first think. Thus, in 6.37 we pick off the value of the surface layer at the origin. This reminds me of E = 4πσ just above a plane with surface charge σ. So the interesting results here are that the limit of the Laplace potential u is the non-local integral 6.36, whereas the normal derivative of this same thing is the local value shown, depending only on the source right under the approach point.( This local result gets a non-local piece later when σ ≠ plane. ) (2) Double Layer v and ∂3v( 116 results 6.38 and 6.40). The solution v -- as I already noted earlier -- has the exact same form as the object ∂3u above if we just replace simple layer density a with double layer density b and have a + sign out front. Thus we just crib the limit in 6.37, change the sign and do a→b and we have 6.38. So this is nothing new. Now we want ∂3v and this is the first time we have computed such an object. I agree with all the math down through line D where we have the 2D Laplacian appearing inside the integral. If we apply the 2D Green's #2 we can move this to the b function and pick up a boundary integral , but now be careful. We are talking about an integral over our flat σ in the x,y plane but now we are going to imagine that this planar surface is bounded by a closed contour C. The "volume" is the area of the σ surface, and the "area" is the curve C. So now we need to know the outward normal derivative along the curve C of our double layer source function b, so our result is quite complicated now. Our next task is to take the x3 → 0 limit of this mess. As he points out, both terms are convergent in this limit, so we just do it under the integral sign and the result is 6.40, the thing we are seeking. To recap: we quickly arrived at ∂3v = 6.39, but we can again see that if we try to take the x3→ 0 limit, we cannot just take it under the integrand because the integration has a highly divergence region around the origin in that limit, the same kind of problem we had with 6.37. But we found the desired limit and it is the mess shown in 6.40 where we have to retain the nasty C contour integral. E. Now consider the actual curved surface σ. ( 117) Has tangent plane at point s, normal there ν. A nearby point on σ will be called ξ with normal n. We are now going to study our four cases of section F which we can call u, ∂νu, v, ∂νv. [ p 118 top here ] As usual, we break our integrals into a central core ε disk σε around s, and the rest outside that disk σ - σε. (1) Simple Layer u and ∂νu (118 results 6.41 and 6.43) For u, things are not very exciting. We know that as x → s (we approach the surface), the little disk integral goes to 0 as ε→0, so nothing is divergent, so we can just take our limit x → under the derivative. The result is 6.41 (compare to 6.36 but here we are approaching s, not 0) This just says that for any σ you can take your observation point right onto the surface and nothing bad happens. We saw this all before. For ∂νu we get result 6.43. Compare this to 6.37. The ±part is the same, but now we have a new integral over σ which was not present in 6.37. In the earlier result, σ was a flat infinite plane. If we had σ be this plane in 6.43, then cos(90) = 0 everywhere and that term vanishes. Again, it is the "core" which creates the first term in 6.43, and a general σ will still have the second term. This term is what you would have gotten had you naively taken x → s through the integral in 6.42. (2) Double Layer v and ∂νv ( 119 results 6.44 and 6.49). Now for the first time we get a difference between ∂νu and v, and this difference only appears when σ is not flat. In 6.42, ν is the normal at point s and is fixed in the integral. But doing v on p 119A, now n is the local normal to σ and varies over the integration. In the core part, the surface really is flat and we get exactly the same core term, ν and n are the same there. But the residual term still has n in it. Compare this 6.44 to 6.38. Again for flat σ the cos = 0 everywhere on σ. For some reason, at this point Stak pauses to express 6.44 in terms of a defined kernel k to get 6.46, nothing new here at all. But some words are needed to understand how k is really defined, and the confusion is all just in the cosine part. In the numerator of k(s,ξ) we have the dot product of the unit vector of s-ξ and the normal vector at the START of that vector (which is n). For k(ξ,s) the same rule applies, but now the angle is completely different. Thus, there is no connection between k(s,ξ) and k(ξ,s), this is not a symmetric kernel! He then defines a different k (but still calls it k) and rewrites 6.43 using that k. Since k is a symmetric kernel, I suspect we shall be calling upon some Volume I work in the future. For ∂νv we just replay earlier work and end up with 6.49. But a big difference now. We cannot do the limε→0 and we have to keep it in each term! The first has σε removed from the integral, the second has Cε as the boundary of σε . This certainly is crashingly ugly. (3) Two of our four limits above can be expressed using the kernel k notation. Here they are limx→s ∂νu = ∓a(s)/2 + ∫σ dSξ a(ξ) k(ξ,s) // 6.43 becomes 6.47 limx→s v = ± b(s)/2 + ∫σ dSξ b(ξ) k(s,ξ) // 6.44 becomes 6.46 We are soon going to use these equations where the LHS is a prescribed boundary value of something on our surface σ. The integrals remind us of what appears in "integral equations" in Volume 1. Notice how the two objects shown here are NOT the same, though they were the same for the flat σ where of course k was symmetric! Comment: Why ? have we done all this work on evaluating the four objects u, ∂νu, v, ∂νv as we approach the surface σ at point s, the first two cases being for a layer of simple source, the second of dipole source. The source functions recall were a and b on σ. I think the answer is this: when we get around to throwing in our boundary values of solutions again, those are the above four quantities ON the surface! So we need to know how to properly evaluate these four babies. Remember that with a p=2 equation, we might have interest in the solution w and its normal derivative! As he always says, we have our easy ways to get the tangential derivatives. I was not aware that the dipole surface was so important as to get all this treatment, but I'll bet I will find out why. One example would be the heart wall I think. Maybe some bio problems will arise using this. In the above two expressions for limits, here is what we are really saying: limx→s [∫σ dSξ a(ξ) k(ξ,x) ] = ∫σ dSξ a(ξ) limx→s [ k(ξ,x) ] – a(s)/2 Stakgold's effort here is to explain why you get that "extra term" when you interchange the order of limit and integration. He does so rigorously within distribution theory. Exercises (120) Exercise 6.18 In 2D space, find u and v due to layers on σ = C, a curve. (120) For the single layer, we just replace the Green's with no ado and get 121 A. Compare this to 6.41, but here we have not yet taken the limit x→s. For a double layer, the dipole Green's must be found, then that goes into the integral. We repeat the above, usol(x|ξ) = limit h→ 0 of the two single delta charge potentials = Dl(x|ξ) // notation! = limh→0 (1/2π) [ (1/h)ln(1/|x-(ξ+h)|) - (1/h)ln(1 /|x-ξ|) ] = (1/2π) limh→0 { [ f(ξ+h) - f(ξ) ]/h } = (1/2π) ξf(ξ) = (1/4π) ∂l f(ξ) But we can go compute ξf(ξ) = ξ (ln(1/|x-ξ|)) = - r(ln(1/r)) where r = x - ξ = - (1/[1/r]) ∂r(1/r) = - r ∂r(1/r) = r /r2 = /r = (x ξ)/ |x-ξ| So therefore we get ξ (1/|x-ξ|) = (x ξ)/ |x-ξ| and Dl(x|ξ) = usol(x|ξ) = (1/2π)∂l f(ξ) = (1/2π) ξf(ξ) = (1/2π) (x ξ)/ |x-ξ| = (1/2π)cos(x ξ, )/ |x-ξ| = cosθ/2πr so everything is exactly the same except the denominator has only one power of r and 2π instead of 4π. Thus we arrive at p 121 B. We can show that B is fine as x moves onto the contour C as follows. Consider a little interval of C around point x on C, centered on same. Length is 2ε say. In this range do a linear fit to the cosine. Without loss of generality, let x = 0. ∫-εε dξ b(ξ) [k (ξ) ] (1/2πξ) = b(0) (k/2π) 2ε → 0 The integrand is ξ/ξ = 1 and so is non-singular at x = 0, just as he claims. 2D Summary: changes 3D to 2D: (1) ln(R) nor R; (2) 1/R1 not 1/R2 ; (3) 2π not 4π both u(r) = ∫σ dS' q(r') 1/2πln(R) R = r - r' r' is on σ simple layer v(r) = ∫σ dS' b(r')cos(R,n')/2πR R = r - r' r' is on σ double layer Exercise 6.19 In 2D space, compute the two distributional results for E (121) In other words, what are the 2D versions of our 3D results (1) limx3 → 0± { E(x1, x2, x3) } = E(x1, x2, 0) in full symbolic function sense. (2) limx3 → 0± { ∂3 E(x1, x2, x3)} = ∓δ(x1) δ(x2)/2 in full symbolic function sense I will use the 2D coordinates x1 and x3 ! (1) : Item (1) is trivial, repeat page 114 logic. (2) : But the second will be different. Start page 114 F = ∂E∂x3 = ∂[(1/2π) ln (1/r) ] ∂x3 = -x3/2πr2 Again we just lose a power when we do the derivative, see above. Now we need new versions of page 114 G and page 115 H. ∫ (2/2π) x3/r2 dx1 =(2/2π)∫ x3/(x12 + x32) dx1 = (2/2π) x3 ∫dx1/ (x12 + x32) = (2/2π) x3 2 (π/2x3) Schaum p 95 = 1 for x3 > 0 // but = -1 for x3 < 0 But now look at the "outside the core piece" limx3→0 ∫ε∞ (2/2π) x3/r2 dx1 = 0 just because x3 = 0 and r2 non-singular. Thus we argue that we have a delta sequence limx3→0 (1/π) x3/r2 → δ(x1) Then we can argue that F = ∂E∂x3 = ∂[(1/2π) ln (1/r) ] ∂x3 = -x3/2πr2 = - (1/2) δ(x1) in the x3 limit. Then limx3→0 ∂3u = " ∂3u" = limx3→0∫ dx1 a(x1) [-x3/2πr2] = ∫ dx1 a(x1) (-1/2) δ(x1) = - (1/2) a(0) x3 > 0 and the other sign for x3 < 0. So one of our conclusions above is then: limx3→0± ∂E∂x3 = limx3→0± { -x3/2πr2 } = ∓ (1/2) δ(x1) 2D which is the result he quotes in this problem. Exercise 6.20 More 2D results (121) In 2D we got the exact same pick-off function ∓ (1/2) δ(x1) that was ∓ (1/2) δ(x1) δ(x2) in 3D. The constant is the same. So when you go through the details we did in 3D, the 2D results are exactly the same except there is one less power in the dipole or derivative kernel (and 2π not 4π). This is good, so the results are pretty much "form invariant" whether for n = 2 or n = 3. Maybe things are similar for other n. 6.5 Integral Equations of Potential Theory (122) For the first time, we get a clear definition of the Dirichlet Problem for Laplace L: Find the Lu=0 solution u inside a closed boundary on which u is specified. Of course there would then be an interior and an exterior solution to ponder, but I think the exterior is still "interior" in the sense that r = ∞ is the other surface. The Neumann Problem is exactly the same, except it is ∂nu which is specified instead of u. A. The double-layer integral equation method for solving Dirichlet (122 3D) The logic which leads to an integral equation here is interesting. Suppose there were some double-layer function b that made u be the same as our single layer f. Then 6.55 would be true. To meet our BC f, then 6.56 must be true. Then p 122 A must be true. This is like 3.6 Vol 1 page 195 where f is the driving function of Fred inhomo 2. However, all of Chapter 3 was with one variable, and here we are in n=3 space, so this is a subject not yet broached by Stak: integral equations with n variables. Would you call this a Partial Integral Equation? So, suppose you could solve integral equation 6.57 for function b, given BC function f. Then the answer to your problem is given by 6.55 where you insert your solution b. You have thus found a specific double layer b to go on σ which solves your Dirichlet problem with BV f on σ. As we know from Chap 3, at least in 1 dimension, a lot is known about solving integral equations, both exactly and approximately. In our case, the homo equation is 0 = -1/2 b + Kb or Kb = 1/2 b. This homo equation has solution(s) only if 1/2 is an EV. That is why he says that if μ = 1/2 is NOT an eigenvalue, then the solution will be unique -- just the particular solution. But he has not proven to me that a solution of this integral equation exists. The reader is a little rusty on Chapter 3 by this time in his Stak voyage! If the nullspace is empty, Ax = f ought to have a solution for any f. Is our kernel here H-S ? Does not seem so offhand, but maybe, Stak is silent. Stak claims we have already shown that 6.54 has at most one solution. Stak regards it as "known" that our multidimensional 2nd kind Fred equation has a solution, but I cannot remember why I should know that. Perhaps he will have more on this soon. Who is credited with making this connection to the double layer? Maybe Dirichlet? I see that he married the sister Rebecca of Felix and Fanny Mendelssohn ! Amazing! Hold on this question. Question about the reality of the dipole layer Question: Is there in fact a real actual dipole layer on the boundary?? We do know that, were we to put such a dipole layer on the boundary, we would create f(φ) on the boundary and we would create the right solution inside. I think the answer is that such a dipole layer is merely one way to create f(s) but that there are other ways as well. Here are two arguments: (1) For example, I think you could find a third "situation" where the same f(s) and internal solution u is created by a single layer charge density a(s). For example, it seems to me that we could just repeat the above integral equation discussion of this section and say this: u(x) = ∫σ dξ a(ξ)1/4πR (**) limx→s u(x) = ∫σ dSξ a(ξ)1/4πR = f(s) We regard the second equation as a Fred inhomo 1 variety. We solve if for a(ξ) and insert that into (**) to get u. So I would extend the list above to say (1) situation where we impose f(s) on boundary σ (2) situation where we construct a dipole layer b(s) on σ (3) situation where we construct a simple layer a(s) on σ Obviously, we don't expect a(s) = b(s). (2) My second argument is this. We know that for the disk or sphere we can say u = ∫ PK f where PK is the appropriate Poisson Kernel and f is the boundary value. We could rewrite this two ways u = ∫ PK f = ∫ (f PK/G) G = ∫aG a = f PK/G G = monopole Green's u = ∫ PK f = ∫ (f PK/D) D = ∫bG b = f PK/D D = dipole Green's So here are explicit formulas for the equivalent simple and double layer densities a and b that would product the same solution u and which have the same f(s) on the boundary. Of course this only applies to a spherical boundary, not a general boundary σ. So imagine a real sphere and we panel it with a million little plates and hook each to a battery in a network like this So in this method we are really creating a surface charge density a with our million little capacitors and we produce some f(θ,φ) in this manner on the surface. There will of course be some u inside the sphere given by u = ∫ PK f . Alternatively, we could actually place a million little electric dipoles on the spherical surface, each made of two charges and a little stick. Each dipole can have a different p value by varying the charge or the stick length (let's say charge). We arrange this pattern to be b = f PK/D and we obtain the same solution inside as the previous paragraph. So we wonder why Stakgold did not talk about this a(ξ) integral equation method? Probably there is a good reason. [ It comes up in the very next section! ] B. A way to approach the combined interior/exterior Dirichlet problem. (123 nD) (a) first, the exterior problem. In "definition" he defines the exterior problem including u=0 at r = ∞ and f on σ. Solution is called ue for "exterior". Lemma 1 is obvious fact that ue drops off as 1/r or faster, since those are radial solutions that drop off. And of course ∂rue would then be 1/r2 or faster. Lemma 2 is another easily computed fact, just as he presents it: as r→∞, an integral vanishes. Imagine our external region Re bounded on the inside by some surface σ, and on the outside by this Lemma 2 surface at σr at ∞. By Lemma 2, the surface contribution from the infinity surface is 0, so we can make this statement of Green's 2nd Identity for Re: ∫Re dV [ ue2E – E 2ue] = ∫σ dS [ue∂nE – E∂nue] where = - on σ (by which I really mean it is pointing in, not really ) . However, since we are going to combine things below with the "interior" solution, we are going to redefine n right now so that = +, so to speak in the same way, so we replace n → -n in the above equation to get ( several steps here) ∫σ dS [– E∂nue + ue∂nE] = ∫Re dV [ – ue2E + E 2ue] But now on the RHS we can set 2ue = 0 from 6.58, and -2E =δ(x-ξ). So we get ∫σ dS [– E∂nue + ue∂nE] = + ∫Re dV ue δ(x-ξ) = [ue(ξ) if ξ is in Re , 0 otherwise ] This seems to be what 6.64 is saying! Of course for fixed ξ, it you move the surface σ out (just scale it up, say, regardless of shape), eventually ξ will NOT be in Re and then we get 0 and that agrees with Lemma 2, whew! The lemma 2 is only for a sphere, but I think we get the idea. He tells us more directly at the bottom of page 124 how to derive 6.64, but I just did it above. (b) next, the interior problem. Now in 6.61 we state the interior problem. By elementary fiddling, I agree that 6.63 is true. It has the opposite sign sense compared to 6.64 due to the normal n direction for Ri. If we add 6.63 and 6.64, the two ∂nE terms cancel because we claim solution u is continuous at the surface, so ue= ui there. Then we get 6.65 which says the following: both ui and ue are provided by the same equation which, in fact, is just a "simple layer" with layer function being I(ξ) as in p 125 A. I completely agree with this development. This I is just the difference of the normal gradients on the boundary surface between in and out. So if we only knew this gradient difference I think, we would have our complete in and out solution! But all we know is the surface f, we don't know the u or their gradients anywhere, so how does this help? Well, I agree that for a point on the surface, either 6.65 or 6.66 just equals f, so we get the "integral equation" 6.66 which has the form KI = f which is an Ax=y type Fred inhomo 1. Not clear it has solutions for any f. Not even clear that I is meaningful if f is not continuous on the surface σ. But at least this is another path one can perhaps walk down, and he will now do an example. Once you have I, jam into 6.65 and there is your solution for both in and out! Comment: I(ξ) is exactly the thing I call a(s) in my meta notes discussion "Question about the reality of the dipole layer" item (1), where I come up a simple-layer integral equation instead of a double-layer one. So the upshot is that the correct a(ξ) is in fact this I(ξ) thing! Stak goes on to say there are problems with this method when f(s) is not continuous on σ, but claims that the method works well for "the exterior of a thin disk" problem and in fact you must use this method. Summary of 6.65. We imagine a Dirichlet problem with some closed boundary σ and some prescribed f on σ. There is an interior problem on Ri and there is an exterior problem on Re, solutions are ui and ue. He has shown the BOTH these solutions are given by the same integral u = ∫EI dS where I has the form shown in p 125 A. In this equation, points out from the surface σ, so the little arrow would have its tip in Re. It turns out that this 6.65 will get used quite a lot later in this chapter. Example (125) : Simultaneous computation of u in and out for the unit sphere (3D) We are going to try the scheme offered on page 125 where we solve integral equation 6.66 for I, then we put that solution I into 6.65 and out comes our u both inside and outside our boundary. Our geometry is of course the simplest, the unit sphere. We are off and running with 6.67 where primed runs over the spherical surface. We then install the usual 1/R Y expansion from the appendix which is p 126 A and which we are going to jam into 6.67 like so: 4π f(Ω) = ∫dΩ' I(Ω')1/R R = |r-r'| = ∫dΩ' I(Ω') ΣnΣm [ Non/Nmn] Ynm(Ω) Y*nm(Ω') = ΣnΣm Non Ynm(Ω) (1/Nmn) ∫dΩ' I(Ω') Y*nm(Ω') = ΣnΣm Ynm(Ω) (Non Imn) But now expand the LHS and compare 4π f(Ω)= ΣnΣm Ynm(Ω) 4π fmn => 4πfmn = Non Imn Imn = (4π/N0n) fmn Now we know that I(Ω) = ΣnΣm Imn Ynm(Ω) so stick this into 6.65 to get (this is a messy set of steps) u(Ω) = ∫dΩ' (1/4πR) { Σn'Σm' Im'n' Yn'm'(Ω') } R = |r-r'| We now use the horrendous monster A.8 which says this (1/4πR) = ΣnΣm r<n r>-n-1 Ynm(Ω) Ynm(Ω')/[ (2n+1)Nmn] So in we go to get u(Ω) = ∫dΩ' (1/4πR) { Σn'Σm' Im'n' Yn'm'(Ω') } = ∫dΩ' (ΣnΣm r<n r>-n-1 Ynm(Ω) Ynm(Ω')/[ (2n+1)Nmn]) { Σn'Σm' Im'n' Yn'm'(Ω') } = ΣnΣm r<n r>-n-1 1/[ (2n+1)Nmn] Ynm(Ω) Σn'Σm' Im'n' ∫dΩ' Ynm(Ω') Yn'm'(Ω') = ΣnΣm r<n r>-n-1 1/[ (2n+1)Nmn] Ynm(Ω) Σn'Σm' Im'n' δnn'δmm' Nmn = ΣnΣm r<n r>-n-1 1/[ (2n+1)] Ynm(Ω) Imn = ΣnΣm r<n r>-n-1 1/[ (2n+1)] Ynm(Ω) (4π/N0n) fmn Now we pause to note that 4π / { (2n+1) N0n} = 1 // see pencil bottom of page 395 so we resume with = ΣnΣm r<n r>-n-1 Ynm(Ω) fmn and this is in agreement with page 126 C and D ! Remember that "the other r" is r = 1 since it is on the surface of the unit sphere. The result is pretty simple I would say, especially when written out as in C and D: the only in/out difference is which power you put on r. He claims that the final result is valid even when f is discontinuous, such as being 1 on one hemisphere and 0 on the other. So this was a fairly impressive way to solve this problem. Question: Why did our integral equation diagonalize to such a simple form? Notice how the Stak presentation in the above example seems very complicated and mysterious, though it is of course correct. The reader must wonder why we get such a simple result Imn = (4π/N0n) fmn when we process this integral equation. I explain the general idea pretty well on p 45 of my addition theorem thesis where I look at an integral equation that is a "group measure convolution equation". My equation here is 4π f(Ω) = ∫dΩ1 I(Ω1) [1/R]( Ω, Ω1) and the diagonalized equation according to my write-up should be 4πf lmm' = Σm" I lmm" [1/R] lm"m' , where we project things onto our group representation functions Dlmm'(Ω). But my integral here is not the full group measure integration because the second azimuth (call it ψ) is missing. We could consider f(Ω) = f(Ω,ψ) = f(g) and same for other functions, but there is no actual ψ dependence, so these functions are not general functions on the whole group. Then we have the following, where dΩ1dψ1 = dg1 4π f(Ω,ψ) = ∫dΩ1dψ1 I(Ω1,ψ1) { [1/R]( Ω, ψ, Ω1,ψ1) δ(ψ - ψ1) } 4π f lmm' = Σm" I lmm" [1/R] lm"m' ( I am arm-waving here a bit, but am sure this is generally right). We for sure can now say that f lmm' = f lm δm',0 I lmm" = I lm0 δm",0 so we might as well set m' = 0 on the LHS and we then have 4π f lm0 = I lm0 [1/R] l00 => 4π f nm0 = Inm0 [1/R] n00 where in the last step we use Stak's "n" for what we usually call l, the group rep Casimir label. It must then be a fact that [1/R] n00 = N0n . So the reason things diagonalize like this is that we have something that is a special case of a group convolution equation. No doubt I have this written up somewhere, but let's not go look now, we understand the main point. The Neumann Problem in 3D and for general n ≠ 2 (126) Same as Dirichlet but we specify ∂nu = f on our surface σ. The divergence theorem forces an extra condition which is that f00= 0 as shown in 6.69 -- this is because dS ∂nui = dSui is what appears in the surface integral side of the divergence theorem for ui = 2ui = 0. We can carry over lots of our Dirichlet machinery including an integral equation with something like I in it, but Stak does not do that. He claims that the f00= 0 condition is not required for the exterior Neumann problem. If we look at the result p 127 B which would in fact apply in the exterior problem, the surface would be the unit sphere + the sphere at r=∞. If u ~ 1/r and ∂nu ~ 1/r2 and area ~ r2, perhaps this r=∞ sphere makes a contribution to the surface integral, which relieves us of having to require that B be true on the inner unit sphere, thus allowing us to have f00 ≠ 0. Whatever f00 we have, the solution will offset it on the r=∞ sphere. Check this in example below! Example 1: Interior Neumann on unit sphere. Use our usual Ynm expansion and meet the BC and you obtain result F which confirms that f00 = 0 as a precondition for there to be a solution. Plug that amn back into expansion and result is p 128 A (notice the 1/n), valid for ANY f defined on the sphere (with f00= 0). Constant A is the n = 0 term in effect. General expansion p 127 D seems obvious to me, but he will do it as an exercise later. Noted added 10.22.09. In the 3D Neumann interior problem, if ∂nu = 0 on the entire boundary, then all the fnm = 0 and so u = A (a constant) everywhere inside. In the corresponding 3D or 2D Dirichlet interior problem, u = 0 on the boundary implies u = 0 everywhere inside. In the 2D Neumann interior problem ( see p 128 G), we could use 128 G as our general solution form, then ∂ru = B/r + nrn-1 terms + sum nr-n-1 terms. Of this vanishes on σ, then match all eiφ coefficients to 0, but could then stll have B ≠ 0 if σ does not contain the origin. But consider annulus radii a and b. Need B/a = 0 and B/b = 0, so B = 0. So I guess for the 2D Neumann, we still have our theorem: ∂ru = 0 on σ implies u = 0 everywhere inside. I thought I saw this contradicted somewhere, but I guess not. He said there was something like 6.11 for Neumann, but does not appear in this book. I notice today that the fancy books like M&F that I just downloaded don't really have the detail that Stak has on boundary condition problems. They are more general books that cover lots of topics. Example 2: Exterior Neumann on unit sphere. Same starting expansion for u, but r-n-1 powers. Else solve in exactly the same manner and result is p 128 F (notice the –1/(n+1) and no constant A). There is no constant A because of the extra imposed condition u(∞) = 0 which exactly knocks out that A. Also of course we don't allow divergence at ∞. Now check our conjecture above. For large r, our solution is dominated by the n = 0 term which is -f00/(1) r-1, so ∂nu ~ +f00/r2. The r=∞ surface integral is then going to be –4π f00, the - sign because its normal is inward for this surface term. On the inner surface of course ∂nu = f, so its contribution is 4πf00 for the n = 0 term and yes, they cancel. The Dirichlet and Neumann Problems in 2D (128) I agree with expression G for the general form of the 2D problem interior or exterior. Dirichlet 2D For interior Dirichlet, the log and neg power terms go away, and we set A and bn from our BC on the unit circle. For exterior Dirichlet, the log and A and positive power terms go away, and we set the an from our BC. If you require u(∞) = 0, then you must have A = 0 which in turn means condition f0 = 0. Alternatively, you could allow u(∞) = A ("bounded") and then any f0 is allowed. On page 129 he comments that some physical problems that look like 2D Dirichlet don't fit this mold, and example being the cross section of a charged cylinder which is a sort of hybrid 2D/3D problem. I'm sure we shall see this eventually. Neumann 2D Now what about the 2D exterior Neumann problem? On page 129 he puts a circle Cr outside our given closed curve C and shows that the surface integral over Cr is fixed by f and in general won't be 0 and it is independent of r, so must be true for r = ∞, which then gives us our exterior problem. But then he shows that if you don't include the log(r) term in your solution form (ruling it out, say, because at r = ∞ it diverges), you always have the Cr integral be 0 for r=∞, which contradicts the fact that it is determined by f and is in general not 0. So you need the log term! So we have to accept u = const * ln(r) + A as an allowed form for r → ∞, which we state as u/ln(r) → constant (bounded). I guess the constant could depend on θ, so that is why he says "bounded". All very good. [reminds me of Coulomb scattering with an inverse square potential, those extra lot things appear. ] Exercises (130) ______________________________________________________________________________ Exercise 6.21. Interior Dirichlet unit sphere using double-layer integral equation method. This turns out to be a long and very painful problem with many parts. First of all, go back to the very lengthy page 108 Exercise 6.15 which treated exactly this problem. In that exercise, using gymnastics from Appendix A, we obtained the expression p 109 E for the solution u, which involved cosγ in the integrand. In 6.26, Stak just writes this out in full detail. I was happy with that solution, but the pathway was pretty complex. Here we attempt to solve this same problem by a different method, that of the double-layer integral equation 6.57. Part A. The first step is a geometry problem and for that we need a picture: where α = π-ψ. The thing we want is cosψ. But we have an isosceles triangle which we bisect as shown, and in either of the two triangles so formed we can say cos(α) = [1/2 |s-ξ| ] / 1 = |s-ξ|/2 . But we know that cos(α) = cos(π-ψ) = -cosψ so we get our desired result cosψ = - |s-ξ|/2. The kernel k given on page 122B is then this: k(s,ξ) = (- |s-ξ|/2)/(4π |s-ξ|2) = - 1/8π |s-ξ| Our integral equation for double layer b is then (from 6.57) f(s) = -b(s)/2 - 1/8π∫b(ξ)dSξ/ |s-ξ| agrees with p 130 A Here is a non-geometry way to show this: (here θ is the angle between s and ξ ) cosψ = (s-ξ) ξ/R = (sξ - ξ2)/r = (cosθ - 1)/R = - 2 sin2(θ/2)/R = -2 (R/2)2/R = -R/2 but I still needed the above picture! Part B. Now we are on our unit r' = 1 sphere, so this is really f(r) = -b(r)/2 - 1/8π∫b(r')dΩ' / |r-r'| R = |r-r'| where both r and r' are on the spherical surface, as x and ξ in the picture above. From p 397 J or A.12 we know that 1/R = Σn' Pn'(cosγ) so that leaves us with: f(r) = -b(r)/2 - 1/8π ∫b(r')dΩ' Σn'Pn'(cosγ) My inclination at this point is to choose the z' axis in the r direction, so cosγ = cosθ' = z'. But this is wrong because b(r') is not spherically symmetric. Instead, project both sides onto harmonics by applying (1/Nmn)∫dΩ Y*nm(Ω) to both sides, giving fmn= -bmn/2 - 1/8π ∫b(r')dΩ' (1/Nmn)∫dΩ Y*nm(Ω) Σn'Pn'(cosγ) Now expand Pn'(cosγ) using the famous addition theorem, which he states best on p 126 A, Σn'Pn'(cosγ) = Σn" N0n" Σm" (1/Nm"n") Yn"m"(Ω) Y*n"m"(Ω') leaving us with this huge mess fmn= -bmn/2 - T2 T2 = 1/8π ∫b(r')dΩ' (1/Nmn)∫dΩ Y*nm(Ω) { Σn" N0n" Σm" (1/Nm"n") Yn"m"(Ω) } = 1/8π Σn" N0n" Σm" ∫b(r') dΩ' Y*n"m"(Ω') (1/Nmn) (1/Nm"n") ∫dΩ Y*nm(Ω) Yn"m"(Ω) = 1/8π Σn" N0n" Σm"∫b(r') dΩ' Y*n"m"(Ω') (1/Nmn) (1/Nm"n") δnn"δmm"Nnm = 1/8π N0n∫b(r') dΩ' Y*nm(Ω') (1/Nmn) (1/Nmn) Nnm = 1/8π N0n∫b(r') dΩ' Y*nm(Ω') (1/Nmn) = 1/8π N0n bmn So our projected integral equation is now this, fmn = -bmn/2 – (1/8π) N0n bmn = –bmn (N0n/8π + 1/2) but N0n = 4π/(2n+1) so get N0n/8π = (2n+1)-1/2 fmn = – bmn ((2n+1)-1/2 + 1/2) = – bmn ((2n+1)-1 + 1)/2 = – bmn (2n+1)-1 (1 + 2n + 1)/2 = – bmn (2n+1)-1 (2n+2)/2 = – bmn (2n+1)-1(n+1) = – bmn (n+1)/(2n+1) Therefore we find that bmn = – (2n+1)/(n+1) * fmn I Is this correct or not ? I have found something on the web that seems to confirm that it is correct. http://books.google.com/books?id=-bV9Qn8NpCYC&pg=PA101&lpg=PA101&dq=dirichlet+unit+sphere+interior&source=bl&ots=ioiH6jT6eK&sig=jpZXQEnM1vvhPAaBylq-vIEbVMI&hl=en&ei=JOFtSvGnDYSotgOliqDKDg&sa=X&oi=book_result&ct=result&resnum=4 This author's version of 6.54 and 6.55 is this: which tells me that his τ is Stak's b. The author then has an error in his 6.3.16 which he quotes as If we divide this by 4π, we get f/4π = - 1/2 τ() + ∫τ kernel dS so I think his LHS should in fact be 4πf and not f ! In any event, the author goes through lots of gyrations and ends up with this result: which is in exact agreement with my result above bmn = – (2n+1)/(n+1) * fmn so amazingly I am able to confirm this part of my solution with a web probe, I have to offer my thanks to Google for scanning all those books!! Part C. Now we have to take our implied function b(r') and install this into 6.55 which says u(r) = ∫dΩ' cos(r-r',n at r') b(r') / 4π R2 R = |r - r'| but now r is inside the sphere so we cannot use our previous formula for the cosine. So we have to compute this cosine differently. Here is one way to do it: cosψ = (r-r')/R r' /1 since = ' and r' = 1 being on the sphere = (rr' - r'2)/R = (rr' -1)/R = ( rr'cosγ -1)/R = (r cos γ - 1)/R So we conclude that cos(r-r',n at r') = (r cos γ - 1)/R I verified that, when r = 1, this gives our result obtained earlier (set γ = 2β, etc). So we now have u(r) = (1/4π) ∫dΩ' b(r') { (r cos γ - 1) / R3 } where γ is angle between r and r' Warning! Do not confuse the notion of rcosθ = z, a Cartesian coordinate of r, and cosα = z, a shorthand for cosα used with Legendre poly work. From now on, I reserve the meaning of z to be z = cosγ. Thus, the above may be written: u(r) = (1/4π) ∫dΩ' b(r') { (rz - 1) / R3 } The big question: What do we do next? I think we need some fancy formula like A.14. I found lots of these formulas in my Appendix A notes, here are a few of them: rz/R3 = Σn=0 rn z Pn'(z) 1/R3 = Σn=0 rn [(n+1) Pn+1(z) + Pn'(z) ]/z If I subtract the second from the first, I get (rz - 1) / R3 = Σn=0 rn [z Pn'(z) – (n+1) Pn+1(z)/z – Pn'(z)/z ] = Σn=0 rn [z2 Pn'(z) – (n+1) Pn+1(z) – Pn'(z) ] /z = Σn=0 rn [(z2-1) Pn'(z) – (n+1) Pn+1(z) ] /z At this point, I can use Schaum p 147 5th to get = Σn=0 rn [n z Pn(z) - n Pn-1(z) – (n+1) Pn+1(z) ] /z so at least the quadratic z2 has been removed. Now rephrase this last line as = Σn=0 rn [n z Pn(z) - {n Pn-1(z) + (n+1) Pn+1(z)} ] /z Schaum p 147 1st says {} = (2n+1)zPn(z). Then we have this: = Σn=0 rn [n z Pn(z) - (2n+1)zPn(z) ] /z = Σn=0 rn [n Pn(z) - (2n+1)Pn(z) ] = – Σn=0 rn (n+1) Pn(z) and there finally is the n+1 I have been looking for! To summarize: (rz - 1) / R3 = – Σn=0 rn (n+1) Pn(z) z = cos γ γ = angle (r,r') r' = 1 Now that I know the result, I see that I could have gotten it immediately from two of the results I already knew, which were these (just subtract them) [...]-3/2 (rz- r2) = Σn=0 n rn Pn(z) [...]-3/2 (1-r2) = Σn=0 (2n+1) rn Pn(z) I wonder how many of Stak's exercise students came up with this. So, we then have u(r) = (1/4π) ∫dΩ' b(r') { (rz - 1) / R3 } = (1/4π) ∫dΩ' b(r') {– Σn=0 rn (n+1) Pn(cosγ) } = – (1/4π) Σn=0 rn (n+1) ∫dΩ' b(r') Pn(cosγ) Now as usual in we go with Pn(cosγ) = N0n Σm (1/Nmn) Ynm(Ω) Y*nm(Ω') to get = – (1/4π) Σn rn (n+1) ∫dΩ' b(r') { N0n Σm (1/Nmn) Ynm(Ω) Y*nm(Ω')} = – (1/4π) Σnm rn (n+1) Ynm(Ω)N0n (1/Nmn) ∫dΩ' b(r') Y*nm(Ω') = – (1/4π) Σnm rn (n+1) Ynm(Ω) N0n bnm = – (1/4π) Σnm rn (n+1) Ynm(Ω) 4π/(2n+1) {– (2n+1)/(n+1) * fmn } which agrees with p 109 6.25 along with p 109 A. When doing that Exercise 6.15, I regarded the pair 6.25 and A as "obvious". We know we have to have a linear combination of powers rn. And we have also to have a lincomb of Ynm so 6.25 is unavoidable. You then try out p 109 A and you will find that the result is u(1,θ,φ) = f(θ,φ). So that is why things are obvious. Here we obtained this result using the double-layer integral equation approach and we had to thread a very complicated path indeed to get to the obvious result (which of course results in 6.26). Alternative last steps ? go back to this point, u(r) = (1/4π) ∫dΩ' b(r') { (r cos γ - 1) / R3 } where γ is angle between r and r' bmn = – (2n+1)/(n+1) * fmn We could have made use of this fact at this point: (r cos γ - 1) = (r2-1-R2)/2 which just comes from writing out R2 and solving for the quantity shown. We then get u(r) = (1/4π) ∫dΩ' b(r') { (r2-1-R2)/2 / R3 } where γ is angle between r and r' = (1/8π)(r2-1) ∫dΩ' b(r') /R3 – (1/8π) ∫dΩ' b(r')/R This would then be an alternative to using this rule (rcosγ - 1) / R3 = – Σn=0 rn (n+1) Pn(cosγ) which gave us u(r) = (1/4π) ∫dΩ' b(r') { – Σn=0 rn (n+1) Pn(cosγ) } Well OK. In the 2D exercise which follows this one, the "could have" replacement get's us to our answer but that does not seem useful in the 3D case. Let's review the threaded steps: Part A. I verified that cosψ = - |s-ξ|/2 by geometry (s and ξ both on unit sphere) and therefore 6.57 says this [ " double layer b causes f(s) at point s on sphere "] f(s) = -b(s)/2 - 1/8π∫b(ξ)dSξ/ |s-ξ| which I prefer to write this way, f(r) = -b(r)/2 - 1/8π∫b(r')dΩ' /R R = |r-r'| Part B. By projecting onto Ynm harmonics and using 1/R = ΣnPn(cosγ) = addition theorem, this equation became fmn = -bmn/2 – (1/8π) N0n bmn = – bmn (n+1)/(2n+1) which we then solved trivially to get bmn = – fmn (2n+1)/(n+1) which result we were able to verify on the web. Part C. We then wrote our solution u(r) according to 6.55 where r lies inside the unit sphere. This is the potential due to our double-layer b function, u(r) = (1/4π) ∫dΩ' b(r') { (rz - 1) / R3 } z = cosγ γ = angle (r,r') where we installed our computed value that cosψ = (rz - 1)/R. At this point, I had to derive for myself this fancy sum rule (rz - 1) / R3 = – Σn=0 rn (n+1) Pn(z) z = cosγ I installed this into the above for u(r) and then expanded Pn(cosγ) = N0n Σm (1/Nmn) Ynm(Ω) Y*nm(Ω') and I ended up with this final result u(r) = Σnm rn fmn Ynm(Ω) which agreed with p 109 6.25 and A, and hence agrees implicitly with 6.26. ______________________________________________________________________________ Exercise 6.22. Interior Dirichlet unit disk using double-layer integral equation method. (130) The equation shown here is correct for the following reason: Let's first assume that the solution to our Dirichlet problem with single layer f is given by the expression shown for w. We can easily show that it does indeed satisfy Laplace, so the only question is this: does w have the single layer boundary value f? If we take the limit x→ s of our w(x) equation, using our already-proven result p121 6.52, and if we force this limit to be f, we get an integral equation for b. If we can solve this for b, then the w expression with this double layer b was a correct form for the solution of our Dirichlet problem with single layer f. The only difference from the text is that here we have 2D versions of things. The kernel k(s,ξ) is shown in 6.53 page 121. But from part A of Exercise 6.21, we showed that the cosine in the numerator of dipole things is given by cosψ = - |s-ξ|/2 and this is the same for a 2D circle as for the 3D sphere where we look at the plane (great circle disk) containing these two points. So we then get 2D: k(s,ξ) = cosψ/2πR = -R/4πR = -1/4π = constant // constant as claimed! Now let's examine our "integral equation" for this 2D problem: f(s) = ±b(s)/2 + ∫dlξ b(ξ) k(s,ξ) = ±b(s)/2 -1/4π∫dlξ b(ξ) 6.52 where the – sign is for the "interior" problem, the + sign for the exterior. Put this into polar coordinates so we have: (r' = 1) f(1,φ) = ±b(1,φ)/2 -1/4π!Syntax Error, Idφ' b(1,φ') but we only think about b on the unit circle, so just write this as f(φ) = ±b(φ)/2 -(1/2)(1/2π) !Syntax Error, Idφ' b(φ') Because k(s,ξ) does not depend on s, the integral does not depend on φ and is just some number. According to p 94 A+B (the Fourier series in complex form), b0 = (1/2π) !Syntax Error, Idφ' b(φ') so the above really says: f(φ) = ±b(φ)/2 -(1/2)b0 Now apply (1/2π)∫dφe-imφ to both sides to get fm = ± bm/2 - (1/2π)∫dφe-imφ (1/2)b0 = ± bm/2 - (1/2)b0(1/2π)∫dφe-imφ = ± bm/2 - (1/2)b0δm,0 = (± bm - b0 δm,0)/2 2fm = ± bm - b0 δm,0 2fm = ± bm m ≠ 0 2f0 = 0 m = 0 and exterior 2f0 = -2b0 m = 0 and interior Let's consider just the interior problem for the moment. Then we have found bm = -2fm m ≠ 0 b0 = -f0 m = 0 We then have our function b(φ) as follows: b(φ) = Σm bm eimφ = -f0 - 2Σm≠0 fmeimφ Write this as b(φ) = +f0 - 2 Σm fmeimφ = f0 - 2f(φ) // interior We want now to plug this into p 130 C, but first we need cosψ when x is not on the disk, just as we had to face this same issue in Exercise 6.21. The answer there was this: cosψ = (r cos γ - 1)/R γ = angle between r and r'. So we then have w(r,φ) = ∫dφ' b(φ')(r cos γ - 1)/R * 1/2πR R = r - r' r' = 1 = (1/2π) ∫dφ' b(φ') (r cos γ - 1)/R2 A quick picture with r' on the circle at φ', and r inside the circle at φ, shows that γ = φ'-φ . Also, we know that R2 = 1 + r2 - 2r cos(φ'-φ) so our result becomes w(r,φ) = (1/2π) ∫dφ' b(φ') (r cos(φ'-φ) - 1)/[ 1 + r2 - 2r cos(φ'-φ)] = (1/2π) ∫dφ'{ f0 - 2f(φ') } (r cos(φ'-φ) - 1)/[ 1 + r2 - 2r cos(φ'-φ)] (*) This is looking a lot like 6.11, but more messaging is needed. Within w above, let's look at the integral that multiplies (1/2π) f0 : ∫dφ' (r cos(φ'-φ) - 1)/[ 1 + r2 - 2r cos(φ'-φ)] so we could shift variables and write this as !Syntax Error, I dφ (r cos(φ) - 1)/[ 1 + r2 - 2r cos(φ)] Integrand is even, so x2 and 0,π, then split into two terms = 2r !Syntax Error, I dφ cos(φ)/ [ 1 + r2 - 2r cos(φ)] - 2!Syntax Error, I dφ 1 / [ 1 + r2 - 2r cos(φ)] I think you do these by using z = eiφ and changing to a contour, but I see the answer on GR p 366 3.613.2 which handles both our cases. We get 2r !Syntax Error, I dφ cos(φ)/ [ 1 + r2 - 2r cos(φ)] = 2r * πr/(1-r2) - 2!Syntax Error, I dφ 1 / [ 1 + r2 - 2r cos(φ)] = - *2 π/(1-r2) So add these two terms up to get 2r * πr/(1- r2) - 2π/(1-r2) = 2π (r2 - 1)/(1-r2) = -2π So the f0 term in (*) above is then (1/2π) f0 (-2π) = - f0 So at this point we then have: w(r,φ) = (1/2π) ∫dφ'{ f0 - 2f(φ') } (r cos(φ'-φ) - 1)/[ 1 + r2 - 2r cos(φ'-φ)] = - f0 + (1/2π) ∫dφ' f(φ') (-2r cos(φ'-φ) + 2)/[ 1 + r2 - 2r cos(φ'-φ)] But we know that the correct result should be 6.11 which says = + (1/2π) (1-r2) ∫dφ' f(φ')/ [ 1 + r2 - 2r cos(φ'-φ)] To see if we even have a chance, suppose we try f(φ) = δ(φ) so that fm = 1/2π for all m. Then my way: - (1/2π) + (1/2π) (-2r cos(φ') + 2)/[ 1 + r2 - 2r cos(φ')] 6.11 way: (1/2π) (1-r2)/ [ 1 + r2 - 2r cos(φ')] If these are equal, cancel 1/2π and examine: -1 + (2 - 2rcosφ)/R2 = (1-r2)/R2 ? -R2 + (2 - 2rcosφ) = 1 - r2 ? R2 = r2-1 + 2 - 2rcosφ = r2+ 1 - 2rcosφ YES! So I think I am on the right track still. I could argue right now that if they are the same for f = δ, then they must be the same for any f. So now we are motivated. Go back to here: w(r,φ) = - f0 + (1/2π) ∫dφ' f(φ') (-2r cos(φ'-φ) + 2)/[ 1 + r2 - 2r cos(φ'-φ)] Then write (2 - 2rcosγ) = 1 - r2 + R2 So we then have w(r,φ) = - f0 + (1/2π) ∫dφ' f(φ')( 1 - r2 + R2 )/[ 1 + r2 - 2r cos(φ'-φ)] = (1/2π) ∫dφ' f(φ') ( 1 - r2) /[ 1 + r2 - 2r cos(φ'-φ)] -f0 + (1/2π) ∫dφ' f(φ') The last two terms cancel, and we then have our official result 6.11! Summary of the above development: We show first that the "kernel" is a constant, as we were requested to do: 2D: k(s,ξ) = cosψ/2πR = -R/4πR = -1/4π The "integral equation for double layer density b in 2D" is then this: f(φ) = ±b(φ)/2 -(1/2)(1/2π) !Syntax Error, Idφ' b(φ') = ±b(φ)/2 -(1/2)b0 The "solution" for b(φ) is thus (for the lower sign which is for interior) b(φ) = f0 - 2f(φ) // interior and b0 = -f0 So we have successfully "solved our integral equation for b". The next step is to install this solution into our w(r) formula which, for the 2D case, is this w(r,φ) = ∫dφ' b(φ')(r cos γ - 1)/R * 1/2πR where we use the same 3D result that cosψ = (r cos γ - 1)/R. We can rewrite this as w(r,φ) = - (1/4π)∫dφ' b(φ')(- 2r cos γ + 2)/R2 then we can make this trivial replacement (- 2r cos γ + 2) = 1 - r2 + R2 to get w(r,φ) = - (1/4π)∫dφ' b(φ')( 1 - r2 + R2)/R2 = - (1/4π) ( 1 - r2) ∫dφ' b(φ')/R2 - (1/4π) ∫dφ' b(φ') = - (1/4π) ( 1 - r2) ∫dφ' b(φ')/R2 - (1/2) b0 We then insert b(φ) = f0 - 2f(φ) and b0 = -f0 to get, = (1/2π) ( 1 - r2) ∫dφ' f(φ)1/R2 - (1/4π) ( 1 - r2)f0 ∫dφ' 1/R2 + (1/2)f0 The first term is our desired result 6.11, so the last two terms must cancel. We must then have (1/2π) ( 1 - r2)∫dφ' 1/R2 = 1 or !Syntax Error, Idφ' 1/R2 = 2π / (1-r2) and this last result is what you get from GR p 366 3.613.2 with n = 0. QED. What about the exterior solution? Retracing the above steps we have f(φ) = +b(φ)/2 -(1/2)(1/2π) !Syntax Error, Idφ' b(φ') = +b(φ)/2 -(1/2)b0 b(φ) = 2f(φ) - b0 and f0 = 0 // recall Stak saying this is required! and b0 is an undetermined constant which I will just set to 0. Then we have b(φ) = 2f(φ). Now going to the next step above, we had w(r,φ) = - (1/4π) ( 1 - r2) ∫dφ' b(φ')/R2 - (1/2) b0 = - (1/4π) ( 1 - r2) ∫dφ' b(φ')/R2 and if we install b(φ) = 2f(φ), we get w(r,φ) = –(1/2π) ( 1 - r2) ∫dφ' f(φ')/R2 which is what I know is the right answer: same as interior but with a minus sign! If we take the large r limit of this thing, R → r so we get w(r,φ) → – ( 1 - r2)/r2 (1/2π) ∫dφ' f(φ') = f0 which agrees with p 129 A which just says the solution is "bounded" at r = ∞. We find that for the unit circle at least, the value at infinity is the same as at the center, which is the MV of the interior solution. Example (PhL) of 2D Dirichlet: Suppose we set f(φ) = cosφ. Then our interior solution is this w(r,φ) = +(1/2π) ( 1 - r2) ∫dφ' cosφ' /[ 1 + r2 - 2r cos(φ'-φ)] This integral is ∫dψ cos(ψ+φ) /[ 1 + r2 - 2r cosψ] = ∫dψ (cosψ cosφ - sinψsinφ) /[ 1 + r2 - 2r cosψ] = cosφ∫dψ cosψ/[ 1 + r2 - 2r cosψ] - sinφ ∫dψ sin ψ/[ 1 + r2 - 2r cosψ] The second integral is 0 due to parity, leaving us with = cosφ ( ∫dψ cosψ/[ 1 + r2 - 2r cosψ] ) = cosφ 2π r/(1-r2) GR p 366 interior r < 1 so our interior solution is then w(r,φ) = +(1/2π) ( 1 - r2) cosφ 2π r/(1-r2) = r cosφ which pretty clearly matches our boundary value at r = 1. The exterior solution will be this w(r,φ) = – (1/2π) ( 1 - r2) cosφ ( ∫dψ cosψ/[ 1 + r2 - 2r cosψ] ) = - (1/2π) ( 1 - r2)cosφ 2 π r-1/(r2-1) same GR integral location = r-1 cosφ which also matches the boundary value at r = 1. If we pick some azimuth φ, we have w(1,φ) = cosφ and we can draw our solution: As you swing around azimuth with this slice, you can see what you would get. The r = 0 and r = ∞ values are both the MV which is 0. Here is what Maple has to say: restart; > readlib(addcoords)(z_cylindrical,[z,r,theta],[r*cos(theta),r*sin(theta),z]); > f := piecewise(r<.99 or r>1.01, .0); // marker not used! > g := piecewise(r<1, r*cos(theta), r>=1 ,cos(theta)/r); > plot3d(g-5*f, r=0..4, theta=0..2*Pi, coords=z_cylindrical, numpoints=1000); So this is a harmonic function in 2D that meets f(φ) = cos(φ) on the unit circle. Along any ray from the origin it is continuous, but has a discontinuous derivative ∂r at r = 1. Recall that = 3d(g) / |3d(g)| is the normal to surface g(x,y,z) = 0. Our interior surface is z = r cosφ = x, so we have g(x,y,z) = z - x. So we compute n = (z - x ) = - n2 = 2 so = ( - )/ which is in a constant direction, which explains why the interior solution lies on a plane as shown. Exercise 6.23 Neumann on 2D unit circle (130) Just above I did the complete 2D Dirichlet on the unit circle, interior and exterior, for any f(θ). Now Stak asks us to repeat this for Neumann instead of Dirichlet, so ∂ru = f(φ) at r = 1. Fine. But this time we get to use the simple expansion method (not the integral equation method), and here we go: A. Interior. General form of solution is this u(r,φ) = Σn=-∞∞ anr|n| einφ .. including the n = 0 term which includes the A of A+Bln(r) and we throw out the n=0 ln(r) term. Thus ∂ru(r,φ) = Σn=-∞∞ an |n| r|n|-1 einφ ∂ru(1,φ) = Σn=-∞∞ an |n| einφ f(φ) = Σn=-∞∞ fneinφ Compare to get fn = an |n| which requires f0 = 0 as he points out in the text. Our solution is then u(r,φ) = Σn≠0 (fn/|n|)r|n| einφ all done Notice that u(0,φ) = 0. What is the value of u on the circle? u(1,φ) = Σn≠0 (fn/|n|)einφ <u> = (1/2π) ∫dφ Σn≠0 (fn/|n|)einφ = Σn≠0 (fn/|n|) (1/2π) ∫dφ einφ = Σn≠0 (fn/|n|) δn,0 = 0 which agrees with f0 = 0. B. Exterior. General form of solution is this u(r,φ) = Σn≠0 anr-|n|-1 einφ + A + B ln(r) ∂ru(r,φ) = Σn≠0 an (-|n|-1) r-|n|-2 einφ + B/r ∂ru(1,φ) = Σn≠0 an (-|n|-1) einφ + B f(φ) = Σn=-∞∞ fneinφ Compare to get (-|n|-1) an = fn n ≠ 0 B = f0 Our solution is then u(r,φ) = – Σn≠0 fn/(|n|+1) r-|n|-1 einφ + A + f0 ln(r) and as he notes, we have ln(r) divergence of the result which is OK. Review of Chapter 6 up to this point. First, a little review up to this point. In the first 3 sections of this chapter, we solved the Laplace equation with no sources at all, but with the solution prescribed on certain bounding surfaces. This subject is known as "the Dirichlet problem" (without sources). For simple geometries, we could exactly solve the problem. We started with the 2D unit disk, did separation of polar variables, found general forms of solutions, then determining the constants in those forms from the prescribed boundary values. We learned the properties of harmonic functions. In later sections, we learned that these solutions and boundary values could be equivalently generated by certain charge or dipole distributions, but not in these 3 sections. Then in Sections 6.4 and 6.5 we started looking at sources, but only sources distributed on a closed boundary σ, not δ function sources within the volume surrounded by the boundary or on it. These distributed sources are called layers; the monopole and dipole layers are called simple and double layers. By simple superposition of known point-source expressions, we obtained Laplace equation solutions as simple integrals of the form u(x) = ∫σ dSξ c(ξ) K(x,ξ) where c is the layer density and K some kernel function. For the simple layer, the kernel K(x,ξ) is really E(x|ξ) such as 1/4πR which is called a "fundamental solution". In the double-later case, K(x,ξ) is really D(x|ξ) which is the fundy solution for a point dipole source. These integral solution expressions apply in the absence of boundary value impositions on other surfaces σ'. For example, there were no nearby conductors σ' which might develop charge distributions and thus affect the problem solution. We learned in detail how to take the limit of our integral expressions as a point x in the interior approaches a point s on the source-bearing surface σ. In some cases, such as the double layer, you cannot simply take this limit under the integral sign due to the singular nature of K(x,ξ), and in fact a certain "extra term" appears when you properly take this limit. Once we knew how to take the limit, we could then start thinking of u(s) as f(s), a boundary value on σ itself, and we then obtained integral equations relating boundary values f(s) to the surface density c(ξ). In simple geometries, we could solve this integral equation and obtain the same results we got earlier using simpler methods, but the integral equation method applies to all geometries. We understood the notion of "equivalent problems" that are not physically the same. (1) some f is imposed on σ by fiat; (2) a simple or a double layer can be found which actually causes this same f to occur on σ. We obtained a different integral equation for simple versus double layer, though both were driven by f(s) as the inhomogeneous term. The double layer equation has the "extra term" which is that eigenvalue type term appearing in a Fred inhomo 2 integral equation as shown page 195 of Stak I. The simple layer integral equation did not have this extra term and so was a Fred inhomo 1 equation. At this point, for the first time, we dealt with "the Neumann problem" (p 126) which is identical to "the Dirichlet problem" but the boundary value f(s) is for ∂u/∂n instead of u on the surface σ. We found, however, that one must have f00= 0 for the 3D case, otherwise the integral equation has no solution. For Neumann, you take the integral solution expression u(x) = ∫σ dSξ c(ξ) K(x,ξ) and you have to differentiate it once ∂n before taking the limit x→s for the boundary, and that is why we did a lot of work on not just u and v (simple and double layer integrals), but also on ∂nu and ∂nv. That is to say, these latter two are needed for Neumann problems. In both these cases, "extra terms" appear. Our initial work was done in 3D and we considered both interior and exterior problems (to the surface σ). We imposed a requirement that u be bounded at r = ∞ and of course finite at r=0. We found on page 125 that, in the integral expression for u(x) in which c(ξ) is a single layer density a(ξ), a(ξ) = I(ξ) which is the (nameless) difference of ∂nu across the boundary for the internal and external solutions. We used this "method" to find the solution of the unit sphere and we computed I(ξ) as part of the solution. We were happy to see a single integral expression provide both the interior and exterior solution. We then looked at 2D Dirichlet and Neumann problems and a few new features appeared. One was that we needed f0 = 0 for the Neumann, in analogy to the f00= 0 in 3D. But more interesting is that we must maintain the ln(r) solution in 2D which means u is no longer bounded as r→∞. The correct fact is that u/ln(r) should be bounded. NOW for the first time in this chapter, starting in Section 6.6, we are going to consider the Laplace solution in the presence of both a boundary σ and a δ-function point source inside the boundary. As we did with the string in 1D, we require that the solution vanish on the boundary σ. Such a solution is called the Green's Function and is denoted g(x|ξ) as in Volume I. So that is what we are now about to do. But first some comments. Remember that a "fundamental solution" is any solution of LE=δ(x-ξ) in the absence of any boundary condition, while a "Green's Function" must also vanish on the boundary, and we usually write LG=δ(x-ξ) or Lg=δ(x-ξ) in this case. Comment on Integral Equations. In the 1D world of ODE's, in volume I of Stak, we often dealt at the same time with an ODE and with its corresponding integral equation (IE) which always involved a certain function k(x,ξ) called "the kernel" and implied a linear operator called K. Usually the ODE was easier to solve for simple cases, but the IE had nicer "properties" since the integral operator was often bounded or Hilbert-Schmidt. Most of our theorems were proved based on the IE since things are more "controlled" in that framework. For example, we found many theorems relating to the location of the eigenvalues μ or λ of an integral operator. For both the ODE and the IE, there is an eigenvalue problem and the eigenvalues are determined by the imposed boundary values of the solution. Also, the IE world allowed us to prove that the eigenfunctions of the eigenvalue problem formed a complete set in certain common situations, something that was hard to show in the corresponding ODE world. In volume II we are in the world of multiple variables and PDE's. We might have wondered right at the start, had we been curious, if there were an integral equation associated with each PDE. In the above sections, we discovered that there are in fact several "associated integral equations" related to a given PDE. Of course we have specialized in Chapter 6 to only a single PDE, the Laplace equation. In the associated integral equations, the variable in the equations is not the solution u (as it is in the ODE world), but is an ancillary function which we interpret as either a single or double layer density. So at least we have something in the way of an integral equation for a PDE. Perhaps we will learn more about this subject in upcoming sections. Maybe there is an integral equation directly for the solution u whose kernel then would be the inverse of the PDE operator, that is to say, K = L-1. This subject did not occur in Chapter 5, our only other Stakgold chapter involving multiple variables (volume I was all single variable). We certainly know that scattering theory uses integral equations which are associated with a certain potential V(r). The PDE involved here is the classical wave equation or the Schrodinger equation of some sort. In either case, we use Hamiltonian mechanics as our calculational tool and that is how V(r) gets into the picture: K + V(r) = H = E. Perhaps this is only a PDE in QM, not in classical mechanics. Just wandering here. [ In fact, integral equations are going to come up right away in the next section: the variable will be not u, but g, the special Green's Function solution. ] 6.6 Green's Function for L = -2 (130) Definition and Properties: Theorems 1 thru 6. (130) The "Green's Function" g(x|ξ) on an open and bounded region R is the solution not of Laplace, but of Laplace with a δ(x-ξ) source term, AND with g = 0 on a boundary σ fully enclosing the region R. This is like the string thing, where the two ends were tied down to y = 0, but now we are in Rn. If we solve our δ sourced Laplace without regard to boundaries, we get E(x|ξ) = 1/4πR, for example in n=3, and we now see why he is using "E" for the fundamental solution. By itself, of course, this is in general not going to be 0 on any boundary! Whatever the solution is for g, we could of course subtract out E and call the residual v, as he does in 6.75. Then v is a pure Laplace solution which meets the BC v = -E! So think about this. Function v is the solution to a Laplace equation with a prescribed BC v = -E on some boundary. That is to say, we have a specific requirement on the v boundary value (which of course is -E so the sum for g will be 0 on σ). BUT, we know all about solving this problem for v(x,ξ). We have been studying that in the previous section in great detail. It is "the Dirichlet BV problem". Now let's digress a moment to look at the theorems quoted: Theorem 1: g exists and is unique. Well, we know that the solution v exists and is unique, and of course E is what it is, so g must exists and be unique. Theorem 2: The function g(x|ξ) is symmetric. This is not as obvious, so let's do the work he outlines. -2g(x|ξ) = δ(x-ξ) g = 0 on σ // definitions -2g(x|η) = δ(x-η) g = 0 on σ - g(x|η)2g(x|ξ) = g(x|η)δ(x-ξ) // process - g(x|ξ)2g(x|η) = g(x|ξ)δ(x-η) - g(x|η)2g(x|ξ) + g(x|ξ)2g(x|η) = g(x|η)δ(x-ξ) - g(x|ξ)δ(x-η) // subtract ∫R dx { g(x|ξ)2g(x|η) – g(x|η)2g(x|ξ)} = g(ξ|η) - g(η|ξ) // integrate over R The purpose is to get the RHS as shown, which we would like to claim vanishes. Use Green #2 on the LHS and it becomes LHS = ∫σ dS { g(x|ξ) ∂n g(x|η) – g(x|η) ∂n g(x|ξ) } = 0 since both g's vanish on σ, and since we assume no singular behavior of ∂n g(x|-) on σ. QED. But if g(x|ξ) = g(ξ|x), then this must be true as well for v(x,ξ) in 6.75, since it is obviously true for E. So if we put unit charge at B, we get a certain v(A,B). But if we put the unit charge at A, then at point B we see v(B,A). These numbers are the same, v is symmetric. This fact is not really obvious when surface σ is present (on which g = 0). The fact is called the reciprocity principle of electrostatics. Theorem 3: g is positive throughout R. Proof: put a tiny ε sphere around ξ and we know that distant σ has no effect, and we know that g ~ E in this sphere and so g > 0 throughout said little sphere. And of course we know that g = 0 on all of σ. For region R - Rε we know that max must be on the boundary, so the max then must be on that Rε boundary, so we know that max of g over R-Rε = Mε > 0, some number. Similarly, we know min of g over R-Rε = mξ = 0 since this min must happen on the σ boundary. Thus, everywhere in R-Rεwe have that 0 < g < Mε so we conclude that g is positive in R. We don't have g = 0 on R because R is an open region and thus does not contain its boundary σ. Theorem 4: For n ≥ 3, g < E throughout R. Proof: g = E + v. Function v cancels E everywhere on σ, so v must be negative everywhere on σ, since we know E is everywhere positive for n ≥ 3. Since v takes its maximum value on σ, that max value must be a negative number. Thus, v < the negative number throughout R. Thus, v is negative everywhere in R, so g = E+v implies g < E in R. In 2D, E might have either sign in R since ln(1/R), so this theorem then fails. You can think of the potentials of the induced charges on the boundary (of opposite sign) as neutralizing some of the potential of the positive point charge, so this theorem is very intuitive. Theorem 5: The integral operator G implied by g is completely continuous. A "completely continuous" operatior is I think nowadays called a "compact operator". A the end of Chap 2 we got various useful theorems about CC operators. Then in Chapter 3 we showed that H-S operators are CC. I think we know something about the eigenvalues of a CC operator, no doubt Stak will remind us soon. A CC operator maps bounded sets into compact sets. In Rn, compact means closed and bounded. So I guess the key thing is that CC operator maps a bounded set -- whether open or closed -- into a closed and bounded set. This may relate to the nature of the inverse operator, etc etc. I am going to skip this Stak proof because I am not interested in CC right now! This is a typical Stak thing, to put in a proof like this because he did it once and wants to get it into the record. However, I do take note of Stak's definition of the "diameter" d of σ, and the notion of sphere Rd of radius d around ξ as a sphere that contains all of σ, for use in inequality work. These things d and Rd are used below. Theorem 6: For n = 2 and n=3, G is a H-S operator. The double integral of g2 is bounded, he shows this directly for the n=3 case, a short proof using d and Rd. Of course we know that H-S operators are completely continuous. So I guess we have CC for all dimensions, but only H-S for n=2 and 3. Solution of the Dirichlet Problem. (135) Now instead of a δ source, we have a source function q(x) driving our Laplace equation, and we specify a function f on the boundary σ. This is "the Dirichlet problem with a source". We now do the same process we did above in a different context: -2g(x|ξ) = δ(x-ξ) g = 0 on σ // (6.74) -2u(x) = q(x) u = f on σ // our Dirichlet with source problem - u(x)2g(x|ξ) = u(x)δ(x-ξ) // process - g(x|ξ)2u(x) = g(x|ξ) q(x) -u(x)2g(x|ξ) + g(x|ξ)2u(x) = u(x)δ(x-ξ) - g(x|ξ) q(x) // subtract ∫R dx { g(x|ξ)2u(x) - u(x)2g(x|ξ)} = u(ξ) - ∫Rdx g(x|ξ) q(x) // integrate over R Use Green 2 on the LHS, LHS = ∫σ dSx { g(x|ξ) ∂n u(x) – u(x) ∂n g(x|ξ)} = – ∫σ dSx u(x) ∂n g(x|ξ) // g = 0 on σ Thus we have shown that – ∫σ dS u(x) ∂n g(x|ξ) = u(ξ) - ∫Rdx g(x|ξ) q(x) or u(ξ) = ∫R dx g(x|ξ) q(x) – ∫σ dS u(x) ∂n g(x|ξ) where ∂n is wrt x. We can now interchange x ↔ ξ and say u(x) = ∫R dξ g(ξ|x) q(ξ) – ∫σ dSξ u(ξ) ∂ξn g(ξ|x) and finally use the symmetry of g to write u(x) = ∫R dξ g(x|ξ) q(ξ) – ∫σ dSξ f(ξ) ∂ξn g(x|ξ) // agrees with 6.81 and we have also replaced u by f, its prescribed value on the boundary. Repeat the above in 1D: I went off and reviewed all my "integral theorems" and saw how they are all valid in 1D with the correct interpretation of the surface integral as a two point evaluation difference. So I think the above comes out being this: u(x) = ∫R dξ g(x|ξ) q(ξ) – ∫σ dSξ f(ξ) ∂ξn g(x|ξ) = !Syntax Error, Idξ g(x|ξ) q(ξ) – [ f(ξ) ∂ξg(x|ξ) ] |ξ=ba In our famous string problem, we happen to have q(ξ) = δ(ξ-ξ1) and f(b) = f(a) = 0, so this thing ends up saying simply that u(x) = g(x|ξ1) . We would have to redo that problem for random endpoints for the string to see if it really works out right, but we know it will! Now back to u(x) = ∫R dξ g(x|ξ) q(ξ) – ∫σ dSξ f(ξ) ∂ξn g(x|ξ) // agrees with 6.81 So this is our big result! When we solve our operator equation, we get u = Gq + extra surface term as shown. If f = 0 on all of σ, then sure, get u(x) = ∫R dξ g(x|ξ) q(ξ), but this is not the general result. We know g(x|ξ), and we know σ, so we can compute ∂ξn g(x|ξ) , and so here is our complete solution as a function of q(x) and f(σ) !!! We really have a fairly closed form for the solution of any Dirichlet problem at all! Our answer contains two integrals as shown above. The first is the solution of the Green's problem with u = 0 on σ. The second integral is the solution of Laplace with u = f on σ. We spent lots of time solving this kind of problem, but we NEVER wrote the solution in the form shown above with ∂ng sitting in it. Of course we were not talking about Green's functions in those earlier sections. Think of the first integral u = Gq as the particular solution of Lu=δ. Then think of the second integral as a solution of the homo equation Lu=0. We know we are allowed to add a homo solution to a particular solution of Lu=δ and still have a valid solution of Lu=δ. So the second integral is such a homo solution which meets the BC u = f on σ. The particular solution satisfies g = 0 on σ by the definition of g. He claims that ∂ξn g(x|ξ) = - I(x|ξ) is the induced charge on the boundary. What he means is that if you placed a conducting metal shell at σ and a point charge at x, a surface charge density Σ would be induced on the side of that surface facing the inside of the region R. Here is how you would compute this surface charge density Σ. First, put your Gauss pillbox across the boundary with one side embedded inside the metal where E = 0. The Maxwell equation in Stakgold units says E = +ρ (see 6.80 2V = -ρ), the 3D charge density. Applying the divergence theorem to this equation gives ∫dV ρ = ∫dSE or ρdrA = -E(r)A so ρdr = -E and we would call ρdr = Σ, the total surface charge per unit area. So boom, we have -E = Σ. But E = -dV/dr so we have that Σ = dV/dr. Now in 6.81, ξ is a point on the surface σ, so we want dV/dnξ so to speak. We want V to be the potential at ξ. In the Green's world, V is the potential at ξ caused by a unit source at location x. Since g is symmetric, we can write either V = g(x|ξ) = g(ξ|x), where the second form really makes more sense here (potential at ξ due to a unit charge at x). So then Σ = d [g(ξ|x)]/dnξ. So this confirms the claim that ∂nξg(x|ξ) is the metallic surface inside charge density, and then I(x|ξ) is the negative of the induced charge density induced on the metallic surface at ξ due to unit charge at x. And we expect the induced charge to integrate to -1. The idea of "metallic surface" goes with the idea that g = 0 on the surface. It seems to me from 6.83 that I(x|ξ) is also the Laplace potential you get from putting a delta BV spike (with rest 0) at location ξ on the surface σ, like the Poisson kernel thing of 6.11. You then integrate over σ to get the total potential from all such "parts" of the problem. So I now learn that the Poisson kernel I(x|ξ) is not only the Poisson kernel in 6.11 for 2D or 6.26 for 3D, it is also the negative of the charge that would be induced were you to put a metal surface (think tiny but finite thickness) at σ. Meanwhile, we are now off on another subject, eigenvalue problems. Redo the surface charge density argument. I am going to avoid any mention of "Maxwell's Equations" or "electric field E". Just stick to generic potential theorem. (1) Start with region R and boundary σ. Pour liquid "metal" around σ so R is then a cavity in this metal. Ground the metal so it has zero potential. Put a unit point charge inside the cavity at point x. The solution of Laplace inside the cavity is then precisely g(x,ξ). This is just a physical way to picture the problem. (2) Now, consider a standard "pillbox" which straddles the boundary at some location ξ on σ. It is a tiny cylinder which has cross sectional area A. Assume there is some surface charge Σ(ξ) on the inner surface of σ. Then the total charge (source) inside this pillbox is just Σ(ξ)A. As we make A small, we can think of this as a point charge q(ξ) δ(x-ξ) = Σ(ξ)A δ(x-ξ). (3) Inside the cavity we have g(x|ξ) as noted in (1) above, due to point charge at x inside the cavity. But if we move to a point x very close to our pillbox, g(x|ξ) is completely dominated by the effect of the surface charge Σ(ξ). For such a point x, is as if g(x|ξ) were the solution of the equation –2g(x|ξ) = Σ(ξ)A δ(x-ξ). [This is the weakest step in my derivation, but somehow it must be correct. ] (4) Apply the divergence theorem to the function g using the volume of this pillbox, ( = (x)) ∫dx g = ∫dSg The LHS can be written ∫dV 2g over the pillbox. But – 2g(x|ξ) = Σ(ξ)A δ(x-ξ) as noted above in the region of this pillbox, so we get that LHS = ∫dx2g(x|ξ) = –Σ(ξ)A∫dx δ(x-ξ) = –Σ(ξ)A. The RHS is equal to ∫dS ∂n g(x|ξ) where n points out of the pillbox (toward the interior of our region which is embedded in solid metal). Only the inside disk contributes to this integral and it has the value A ∂n g(x|ξ). [ everything is zero inside the metal, and pillbox sides are taken to limit of no side area ] . The normal n is out of the pillbox, but we really want n to refer to normal out of σ so RHS is then really LHS = – A ∂n g(x|ξ) , and specifically it is – A ∂(x)n g(x|ξ). So we then obtain the result from our divergence theorem that – Σ(ξ)A = – A ∂(x)n g(x|ξ). (5) Now we call upon a little lemma: if f(x,y) = f(y,x), then ∂f/∂x = ∂f/∂y. Prove this by writing out the meaning of the deriviative and using the symmetry of the function. Therefore, cancelling the –A, we have shown that Σ(ξ) = ∂(ξ)n g(x|ξ). (6) Our charge density Σ(ξ) is of course a function of the location of the point charge x inside σ, so we should really write it as Σx(ξ), so we have that Σx(ξ) = ∂(ξ)n g(x|ξ). Stak defines this negative derivative (see top of page 136) to be I(x|ξ), so we have Σx(ξ) = – I(x|ξ) = ∂(ξ)n g(x|ξ) . To summarize: if we implement our Green's Function problem on region R by (1) embedding R inside a material in which we can impose g = 0, and (2) placing a unit charge at x inside R, then there will exist a simple layer source density on σ whose value is Σx(ξ) = ∂(ξ)n g(x|ξ) ≡ – I(x|ξ). Thus, I(x|ξ) is the negative of the induced surface charge. (7) There is more to be said here! Suppose we just had our point charge at x in empty space. The Laplace solution would be ux(x') = E(x'|x), the fundy solution. We later in this chapter define v(x'|x) ≡ g(x'|x) - E(x'|x). If we think of the solution g(x'|x) as being due to the sum of the point charge at x and the distributed charge I(x|ξ), then it would seem that solution v(x'|x) would be what you would get if you could somehow freeze the charge I(x|ξ) on the surface and remove the surrounding metal. We know that on the surface we would in this case have v(s|x) = - E(s|x). Going back to g as the sum of the two contributions, we can write g(x'|x) = E(x'|x) + ∫dSξ Σx(ξ) E(x'|ξ) = E(x'|x) – ∫dSξ I(x|ξ) E(x'|ξ) Eigenvalue Problem for L = – 2 (136) * This is an easy two pages, we revisit results from volume 1 for this particular L = -2. Notice that the EV problem includes u = 0 on σ, just as we had the ends of a string tied down at u(a) = u(b) = 0. EV is called λ. Can be simple or degenerate. He then obtains these classic results: (L is symmetric = Hermitian) eigenvalues λ are all real and positive (hence we can talk about μ = 1/λ and μ are all positive) eigenfunctions of different eigenvalues are orthogonal can cast the ODE as an integral equation with μ as EV. We know that when the kernel is symmetric and CC, the eigenfunctions form a complete set, hence that is true here. Two Examples: Rectangle and Disk (139) Example 1: The classic rectangle situation with u = 0 on the perimeter, edges a in x1 and b in x2. He separates variables, and solves the problem in a perfectly clear manner and obtains the EV's and EF's shown in bracket p 139A. These of course are certain modes of a waveguide, probably TE. [ Yes. In the solution, u = V varies across a cross section of the wave guide, so there are transverse E fields. See Jackson p 244 for how k2 and ω2 enter the full problem. ] At top of page 140 he asks "how do we know we have found all the solutions?". He calls upon a Lemma from Vol 1 which says that if your EF's are a product of 1D complete sets, then the 2D set is complete on its domain, and therefore our double sine EF's form a complete set, QED. Notice how this worked: in our original unseparated PDE, we have EV symbol λ. When we separate, a new symbol μ appears. In each "dimension" we have a little EV problem. You pick one of the dimensions (x1) and you find that μ = km, a certain set of values. Then in the other dimension (x2), the corresponding EV is λ - μ = jn , a certain set of values. Then you have λ = km + jn = λm,n as your final EV. In this problem, km= (mπ/a)2 and jn = (nπ/b)2 . Notice that orthonormality of the solutions is this: < um,n| um',n'> = δmm'δnn'. That is to say, both sets of indices are involved, it is a direct product deal of the two dimensions. Example 2: Repeat the above but on a disk with u = 0 on the boundary, the drumhead perhaps. Again we start with λ, and separation causes μ to appear. We pick the azimuthal EV problem first, and we find that μ = n2 . For any integer n, we have sin(nφ) and cos(nφ) as independent EF's which meet our imposed BC at the wrap line that function and derivative be continuous. So each n is doubly degenerate, and you can also take the EF's of course as einφ . At this point, things are a little different from the Cartesian problem. Our radial equation is p 140 C which does not contain just λ-μ but a more complicated combination of λ and μ and r. This is just the way the algebra works out in polar coordinates. So now with μn = n2, we want to impose our u = 0 at r = 1 and find out what λm are going to be (for given n). Our radial equation is Bessel's in z =r with n2 appearing as in p 141 B, and the solutions finite at r = 0 (he did not mention this BC) are the Jn (r) where λ is our desired EV we want to find. The requirement that u(1) = 0 means Jn () = 0 which means must be a "zero" of Jn and these are enumerated as β(n)m. We end up then with λn,m = [β(n)m]2 . We can then normalize the Φ function as in p 141 A, and the R function using known integral p 141 D and the overall result is p 142 A. Orthonormality is again < um,n| um',n'> = δmm'δnn'. The set is proven complete by appealing to the separate 1 dimensional problems and their completeness. And so ends this fairly easy 5 1/2 page section, onward! Green's Function for Unbounded Regions (142) In this section, Stak describes -- as an example of an unbounded region -- the case where we think of Re as the region "external" to the Ri we considered up to this point, so Re lies outside σ and is thus bordered by σ on the inside, and an r=∞ sphere on the outside. He repeats the earlier development, but now we add a new BV condition which is that in the Green's definition, g = 0 at ∞ as well as on σ. This is reasonable since the surface enclosing the Re region has both these pieces. We end up with exactly the same general result which is 6.92, but now we of course are integrating over Re and the normal vector on σ points inward instead of outward. In 6.92, the vector n points outward, so there really is a sign change relative to the second term in 6.81 where the normal also points outward. So the point here is that for such an Re problem, we have a complete closed form solution to our Poisson equation in terms of q(x) charge and f(σ) prescribed on σ (and being 0 at ∞). So this is no harder than the interior problem. What about "the EV problem" on Re ? In Chapter 4 we saw how, when we took (a,b) to (a,∞), the discrete eigenvalues λi moved closer together and formed a continuum ( branch cut in g(x|ξ; λ) ) so that the possibility of a continuous spectrum appeared, related to the operator B-1 being unbounded (see meta notes there). So that is what happens here since we have r = ∞. He says that operator G is no longer completely continuous (compact, since acts on an infinite domain eg). He points out that you can have Re be external to a set of Ri regions and things are the same. But you can also have an unbounded region that cannot be thought of as Re (external to a bounded region). An example would be the infinite wedge. { spend 2.25 hours in AM session today 8.5.09, * to here} Exercises (143) -- eleven problems here Exercise 6.24. State Theorem 4 for n = 2 dimensions. For n=2 we know that E = (1/2π)ln(1/R) which is positive when R<1, but negative when R<1, so this is the difference with other dimensions, fact that E can be negative. g is still positive in R. We are supposed to show that in the n=2 case, g < E + (1/2π)ln(h) where h is the max of R on R . We also know that g = E + v as in 6.75 where v is harmonic in R. So consider v = g - E. We know that v(σ) = g(σ)-E(σ) = 0 - E(σ) = -E(σ). What is the max that -E(σ) can be? Well, – E(σ) = (1/2π)ln(R(σ)) and this is max at some point σ1 on σ where R(σ) is maxed out. Notice that this point also maxes out -E(σ) and v(σ). But this point cannot be further than h from our point ξ, so we know that R(σ1) ≤ h and then we have shown that v(σ) = -E(σ) ≤ -E(σ1) ≤ (1/2π)ln(h). Now σ1 is the point on σ which maximizes v on σ . According to the min/max harmonic theorem, v inside R must everywhere be less than v(σ1) as well. Thus we have v(x) ≤ v(σ1) ≤ (1/2π)ln(h) . But v = g-E so we have now shown that g - E ≤ (1/2π)ln(h) and thus g ≤ E + (1/2π)ln(h), QED. Exercise 6.25. Show that G is H-S in n = 2 dimensions. Keep in mind that the region R is "part of the problem" when you talk about a Green's function. Change R, and you change to a new problem. So, this result was shown for n = 3 top page 135. In n = 2 we get analogously, ∫∫dxdξ |g|2 = want to show this is less than some finite value. dx = dx1dx2 eg We know from the previous problem that g < (1/2π)ln(h/|ξ-x|) so we have ∫∫dxdξ |g|2 ≤ (1/2π)2∫∫dxdξ ln2(h/|ξ-x|) The integrand is positive everywhere (and so is the ln not squared), so if (for the dx integral) we extend the integration volume from R to disk radius h centered at ξ, we only increase the integrand. That is to say ∫dx ln2(h/|ξ-x|) ≤ ∫disk(ξ) dx ln2(h/|ξ-x|) = ∫rdrdφ ln2(h/r) = 2π∫0h rdr ln2(h/r) But this integral is some finite positive number since rln(r2) = 0 at the origin, suppressing a potential ln(r) divergence there. For example, ∫dx x ln2(x) = which is finite at x = 0. So call this finite number K, then we have ∫dx ln2(h/|ξ-x|) ≤ 2πK(h) // K is a function of h Then we have ∫∫dxdξ |g|2 ≤ (1/2π)2∫∫dxdξ ln2(h/|ξ-x|) ≤ (1/2π)2K(h) ∫R dξ = (1/2π)2K(h)A since the area of bounded R is some finite number A. Therefore, since the kernel g is bounded under this double integral over R, we know that the corresponding integral operator G is Hilbert-Schmidt. Exercise 6.26. Show that the kernel k = a(x,ξ)/Rm (m < n) generates a CC operator G. Well, I skipped his related proof on pages 133-135 because I don't have an infinite amount of time, so I will skip this problem as well. Exercise 6.27. Show that g(x|ξ) in 6.90 is symmetric (ie, for the Re external problem). The new thing here is that we are now talking Re, some external region to σ. I will now copy down and modify my proof above for the Ri = R internal region: -2g(x|ξ) = δ(x-ξ) g = 0 on σ and at ∞ // definitions -2g(x|η) = δ(x-η) g = 0 on σ and at ∞ - g(x|η)2g(x|ξ) = g(x|η)δ(x-ξ) // process - g(x|ξ)2g(x|η) = g(x|ξ)δ(x-η) - g(x|η)2g(x|ξ) + g(x|ξ)2g(x|η) = g(x|η)δ(x-ξ) - g(x|ξ)δ(x-η) // subtract ∫Re dx { g(x|ξ)2g(x|η) – g(x|η)2g(x|ξ)} = g(ξ|η) - g(η|ξ) // integrate over Re The purpose is to get the RHS as shown, which we would like to claim vanishes. Use Green #2 on the LHS and it becomes LHS = ∫-σ + Great Sphere dS { g(x|ξ) ∂n g(x|η) – g(x|η) ∂n g(x|ξ) } = 0 since both g's vanish on σ and on the Great Sphere, and since we assume no singular behavior of ∂n g(x|-) on σ. QED. Exercise 6.28. Find the Laplace eigenfunctions for a 2D disk wedge r = 1 This is similar to some exercise done earlier in the book (wedge on p 107). Separate in polar as per page 140 B and C with μ as sep constant. But this time we have Φn(φ) = sin(nπφ/α) so we get 0 when φ = α starting with n = ±1,± 2... We do have the potential n=0 solution which is A + B ln(r), but B = 0 at r = 0 and deal with A later. The eigenvalues are then μ = (nπ/α)2 because this is what gets pulled down when you compute Φ" But why does the limit α = 2π not work right? That is a slit disk, not a full disk, so OK. So in the disk problem we had μ = n2 , here we have μ = (nπ/α)2. We still have our radial equation p 140 C with λ sitting there, and our new μ on board. Let's define (nπ/α) = n' ≠ integer so μ = n'2. We then get p 141 B as our Bessel radial with n ' sitting in there. The solutions are then of the form Jn'(r) and I think we can dispense with the second kind functions due to r = 0. The Bessel zeros are as before, so our eigenfunctions seem to have this form Φn'(φ) Jn'(β(n')m r) = Kn'm sin(n'φ) Jn'(β(n')m r) n' = (nπ/α) with n = ±1, ±2, ... Question: Why on page 107 are the radial functions powers of r, whereas here they are Bessel functions? It is because here we are doing -2u = λu, the EV problem, whereas there we were doing 2u = 0, the Laplace equation. Here we have λn',m = [β(n')m]2 as our eigenvalues. How do we know we have a complete set? Note you cannot just add constant A to the solution since λ ≠ 0. Normalization might be extra work since n' is not an integer. Rescaling some more? Maybe there is a simpler way. Go back to p 141 B which says this: ∂z(z∂zR(z)) - n'2R(z)/z = 0 z = r n' = (nπ/α) We started with the raw equation which has some solution R1(r; λ) call it, equation 140 C. We then said: R2(z;λ) = R1(r;λ) and this got us to p 141 for R2. Now try scaling to variable y = az with some constant a. Then we will have R3(y;λ) =R2(z;λ) = R1(r;λ) The consider: ∂zR2 = ∂zR3 = ∂yR3 * dy/dz = a ∂yR3 in general ∂z = a ∂y z∂zR2(z) = z a ∂yR3 = y∂yR3 ∂z(z∂zR2(z)) = a ∂y (y∂yR3) Our equation above then becomes: ∂z(z∂zR(z)) - n'2R(z)/z = 0 z = r n' = (nπ/α) a ∂y (y∂yR3)- n'2 R3 a/y = 0 Stupid, that got you nowhere! I think we are stuck with those weird Jn' Bessel functions. Can I find confirmation of my solution on the web? Not quickly. Maple does compute Bessel Zeros for any n, not just integers. And it can compute and plot Bessel functions for any n, not just integer, so I guess my solution is not so bad after all. Bateman p 4 says this: Jν(z) = (z/2)ν0F1(ν+1; -1/4 z2) In our case, ν = n' ≠integer, so I have a = n'+1 and 1-a = -n'. For n = 1,2,3,4 we can see that we have no problems with our Jν(z), but for negative n values, it is singular at z = 0 so we would have to rule it out. But then I wonder about this other solution mentioned above. I guess we could just take that to be Nν(z), called Y by Bateman, and yes, that blows up for any= ν. So I think our eigenfunction set is this: Φn'(φ) Jn'(β(n')m r) = Kn'm sin(n'φ) Jn'(β(n')m r) n' = (nπ/α) with n = 1,2,3... As for orthogonality, for two different n' values, the sines will be orthogonal on the range (0,α). For the same n' (non integer) value, the Bessel's will be orthogonal as well, to wit wiki, where here "α" need not be an integer or half integer, etc. So the Bessel functions of any reasonable and non-integral order are orthogonal when scaled by their respective different zeros! As they point out, this would be true for any Hermitian operator, etc etc. Exercise 6.29. Find the Laplace eigenfunctions for a 2D annular ring between r = a and r = b. We are then back to μ= n2 for the φ side of things, and our Bessel solution is now more general, something like this: fmn(r) = An Jn((β(n)m r) + Bn Nn((β(n)m r) and we then have these two BC's fmn(b) = An Jn((β(n)m b) + Bn Nn((β(n)m b) = 0 fmn(a) = An Jn((β(n)m a) + Bn Nn((β(n)m a) = 0 which we can trivially solve for An and Bn. Our solutions are then like p 141 C, with the Bessel J replaced by the above linear combination with An and Bn figured out and installed. Then of course you have to normalize to get the orthonormals. This brings in the N functions I guess. The world seems to be using Y instead of N on these Bessel second kind functions, by the way. Watson called them Y, Jackson calls them N. I can't find any other N users: A&S, Watson, Bateman. But GR use N, quoting from source KU which is Russian. This would be a good question for Jackson, and I wonder if he has changed it in his later editions. Maybe he wanted to avoid confusion with Y as spherical harmonics. // He has not changed it, here is Amazon from his current 3rd edition same as on p 71 of my first edition. The contents shows that most stuff is the same, some chapters have more sections, and he moved his Multipole Fields Chapter 16 into the new Chapter 9, so one less chapter in the book. So how is normalization going to work here? I don't see a rule like that quoted above, but somehow orthogonality must "obtain". Remember that every boundary "shape" provides its own set of eigenfunctions which are complete and must have orthogonality. Here shape is annular ring. Comments on Neumann's. It took me a while to figure out who Neumann was. The name is overshadowed today by John Von Neumann of latter days. The person is Carl Gottfried Neumann 1832-1925. He did the idea of a matrix series even with infinite matrices (Neumann series). He did the Neumann boundary condition idea. http://www.gap-system.org/~history/Biographies/Neumann_Carl.html "Carl Neumann was the son of Franz Neumann who has a biography in this archive. His mother was Bessel's sister-in-law. Carl was born and received his school education at Königsberg where his father was the Professor of Physics. " The father has his own credits. Franz Ernst Neumann 1798-1895, colleague of Bessel! Specific heat thing called Neumann's law. Since he worked with Bessel. One web article says the second kind function is called N "in the German literature". But here at last we are referred to Karl for this thing. And here from "Mathematics of the 19th Century" we have more good detail: So this is some great history. We have Neumann, Weber and Watson all appearing as players. So I think we can establish 1867 as the date of the "Neumann function" and Carl as its founder. Exercise 6.30. Eigenfunctions of Laplace for the unit sphere in 3D. Lots of text here. The starting section is pretty clear, separate into radial and spherical harmonics, get μ = n(n+1) where n = l , obtain radial equation C which of course includes λ and μ. He then changes the function name to Z(z) = z1/2 R(r(z) where z = r and arrives at equation E which shows the Bessel "order" being n+1/2. [ I am adding the corrections in red. ] I recall seeing this somewhere earlier. Yes, on page 53 we were doing Helmholtz which is our EV problem with a δ source added. The radial equation there we made a function change u = r1-n/2 w and we ended up with a Bessel equation with order n/2-1, so not the same as our current situation. So we get equation p 144 E and conclude that R = Z/z1/2 where Z = Jn+1/2(r) so R = (1/r)1/2 Jn+1/2(r) , which in fact does seem to have an r = 0 problem. From AS at r=0 the leading term in the J is (r/2)n+1/2 so for n=0 we do have a divergence at r = 0 in R(r), but maybe that is OK. [ BUT: After the correction, we get R ≈ (1/r)1/2(r/2)n+1/2 = (r)n 2-n-1/2 which is just fine at r=0 even for n = 0, so this problem I found is fixed by the correction. ] We then look for the zeros to get 0 on our circle at r = 1. This determines the λ values as shown in G, no rocket science here. So reinstalling the 1/r factor, he gets eigenfunctions H. The radial functions are orthogonal which is what p 145 I says in terms of labels k ≠ j. The normalization constant is left to the reader. So this "exercise" was just some optional text material Stak wanted to get in. It certainly sounds like an important case to study and get the results down in print. Certainly the radial functions are not pretty. You could express the results in terms of the spherical Bessel functions jn(βr) but you would still have a factor 1/out front, see AS page 437. So yes, the spherical Bessel functions ARE the solutions to this problem. I think Stak just wanted to avoid dealing with another class of functions. This does give at least a meaning to the phrase "spherical Bessel functions". Again, the eigenfunctions form a complete set over a region R that is in n=3 space. Remember that we want u = 0 on the boundary σ. For a rectangular box, we would just get a product of sines. Compute the normalization constant Nmnk ( this section added 5.4.10) Well, let's do a slight variation of this calculation. Let's assume that, for a sphere of radius a, our normalized eigenfunction takes this form, where Y is the Jackson Y, not the Stagold Y: [ warning: this stuff is all wrong, but it is corrected below ] umnk(r,θ,φ) = Cnmk (a/r)1/2 Jn+1/2(βkn+1/2r/a) Ynm(θ,φ) // Stak p 145 H But I think something is wrong here. The reason is this. Page 145 I tells us to compute !Syntax Error, Ir2dr ∫dΩ umnk(r,θ,φ) um'n'k'(r,θ,φ)* =CnmkCn'm'k'!Syntax Error, Ir2dr∫dΩ(a/r)1/2Jn+1/2(βkn+1/2r/a)Ynm(θ,φ)(a/r)1/2Jn'+1/2(βk'n'+1/2r/a) Yn'm'(θ,φ)* = CnmkCn'm'k' !Syntax Error, Ir2dr (a/r)1 Jn+1/2(βkn+1/2r/a) Jn'+1/2(βk'n'+1/2r/a) ∫dΩ Ynm(θ,φ) Yn'm'(θ,φ)* = δn,n' δm,m' CnmkCnmk' !Syntax Error, Ir2dr (a/r)1 Jn+1/2(βkn+1/2r/a) Jn+1/2(βk'n+1/2r/a) // was (a/r)2 where we have used Jackson p 65 for the Y integral. Keep going = δn,n' δm,m' CnmkCnmk' a2!Syntax Error, Idr (ar) Jn+1/2(βkn+1/2r/a) Jn+1/2(βk'n+1/2r/a) // was no factor Now we know from transforms.doc that this is true, [b Jν+1(xν,k)/]-2!Syntax Error, Idx x Jν(xν,kx/b) Jν(xν,k'x/b) = δk,k' // orthogonality and we can set b = a and ν = n+1/2 and x = r and xν,k = βkn+1/2 to rewrite this as [a Jν+1(xn+1/2,k)/]-2!Syntax Error, I!Syntax Error, Idr r Jn+1/2 (βkn+1/2r/a) Jn+1/2 (βk'n+1/2r/a) = δk,k' This is the fact we want to use, but the weight function is missing from the Stakgold integral!!! My gut feeling is that umnk(r,θ,φ) should have 1/ instead of 1/r in front of it in p 145 H. [ correct! ] So I now return to p 144C which I have verified. (pencil check mark). Let's try this z = r and Z = R R = z-1/2Z instead of what he says dz = dr ∂r = ∂z We will then need first: ∂rR = ∂z (z-1/2Z) = ( z-1/2∂zZ – (1/2) z-3/2 Z ) r2∂rR = (z2/λ) ( z-1/2∂zZ – (1/2) z-3/2 Z ) = (1/) ( z3/2∂zZ - (1/2) z1/2 Z ) ∂r(r2∂rR) = ∂z [(1/) ( z3/2∂zZ - (1/2) z1/2 Z ) ] = ∂z( z3/2∂zZ - (1/2) z1/2 Z ) = z3/2∂z2Z + (3/2) z1/2∂zZ - (1/2) z1/2∂zZ - (1/2)2z-1/2Z = z3/2∂z2Z + z1/2∂zZ - (1/2)2z-1/2Z Equation C now reads z3/2∂z2Z + z1/2∂zZ - (1/2)2z-1/2Z - μ z-1/2Z + λ(z2/λ) z-1/2Z = 0 z3/2∂z2Z + z1/2∂zZ - (μ + 1/4) z-1/2Z + λ(z2/λ) z-1/2Z = 0 // next mult by z1/2 z2∂z2Z + z∂zZ - (μ + 1/4) Z + z2Z = 0 Now set μ = n(n+1) so that n(n+1) + 1/4 = n2+n +1/4 = (n+1/2)2 => μ + 1/4 = (n+1/2)2 so we then have z2∂z2Z + z∂zZ - (n+1/2)2 Z + z2Z = 0 Compare this to Stak p 136 and it is Bessel's equation with ν = n+1/2. We can also write ∂z(z∂zZ) = z ∂z2Z + ∂zZ so our equation above is z∂z(z∂zZ) - (n+1/2)2 Z + z2Z = 0 ∂z(z∂zZ) - (n+1/2)2 Z/z + zZ = 0 so Stak's equation E is correct, but I now have to fix things up!!! So my suspicion was correct, and I have fixed the two errors in red! This is a fairly series typo since we are talking about a standard problem! So let's now REDO our norm constant calculation above, doing things right this time: umnk(r,θ,φ) = Cnmk () Jn+1/2(βkn+1/2r/a) Ynm(θ,φ) // Stak p 145 H Page 145 I tells us to compute !Syntax Error, Ir2dr ∫dΩ umnk(r,θ,φ) um'n'k'(r,θ,φ)* = CnmkCn'm'k'!Syntax Error, Ir2dr∫dΩ()Jn+1/2(βkn+1/2r/a)Ynm(θ,φ) () Jn'+1/2(βk'n'+1/2r/a) Yn'm'(θ,φ)* = CnmkCn'm'k' !Syntax Error, Ir2dr (a/r) Jn+1/2(βkn+1/2r/a) Jn'+1/2(βk'n'+1/2r/a) ∫dΩ Ynm(θ,φ) Yn'm'(θ,φ)* = δn,n' δm,m' CnmkCnmk' !Syntax Error, Ir2dr (a/r) Jn+1/2(βkn+1/2r/a) Jn+1/2(βk'n+1/2r/a) where we have used Jackson p 65 for the Y integral. Keep going = δn,n' δm,m' CnmkCnmk' a!Syntax Error, Idr r Jn+1/2(βkn+1/2r/a) Jn+1/2(βk'n+1/2r/a) Now we know from transforms.doc that this is true, [b Jν+1(xν,k)/]-2!Syntax Error, Idx x Jν(xν,kx/b) Jν(xν,k'x/b) = δk,k' // orthogonality and we can set b = a and ν = n+1/2 and x = r and xν,k = βkn+1/2 to rewrite this as [a Jn+3/2(βkn+1/2)/]-2!Syntax Error, I!Syntax Error, Idr r Jn+1/2 (βkn+1/2r/a) Jn+1/2 (βk'n+1/2r/a) = δk,k' So we can now just "keep going" on the above = δn,n' δm,m' CnmkCnmk' a δk,k' [a Jn+3/2(βkn+1/2)/]2 = δn,n' δm,m' δk,k' Cnmk2 a[a Jn+3/2(xn+1/2,k)/]2 => Cnmk = a-1/2 [a Jn+3/2(βkn+1/2)/]-1 = a-3/2 / Jn+3/2(βkn+1/2) ******* Then our normalized eigenfunction is this: (with Jackson Y) umnk(r,θ,φ) = [ a-3/2 / Jn+3/2(βkn+1/2)] Jn+1/2(βkn+1/2r/a) Ynm(θ,φ) The dimension of L-3/2 is consistent with the normalization integral! As a check, let's look at the case n = m = 0 u00k(r,θ,φ) = [ a-3/2 / J3/2(βk1/2)] J1/2(βk1/2r/a) Y00(θ,φ) = [ a-3/2 / J3/2(βk1/2)] J1/2(βk1/2r/a) (4π)-1/2 But Schaum p 138 tells us that J1/2(z) = (2/πz)1/2 sin(z) = (2/π)1/2 z-1/2 sin(z) which clearly has zeros at z = kπ so we know that βk1/2 = kπ. Then we have u00k(r,θ,φ) = [ a-3/2 / J3/2(kπ)] J1/2(kπ r/a) (4π)-1/2 = [ a-3/2 / J3/2(kπ)] (2/π)1/2 (kπ r/a)-1/2 sin(kπ r/a) (4π)-1/2 Meanwhile, Schaum tells us that J3/2(z) = (2/πz)1/2 (sinz/z - cosz) J3/2(kπ) = (2/π2k)1/2 (sin(kπ)/(kπ) - cos(kπ) ) = (2/π2k)1/2 (sin(kπ)/(kπ) - cos(kπ) ) Now J1/2(z) has a zero at z = 0, though it is a square root zero. Are we supposed to include this or not? For the moment, let's NOT include it, and say k = 1,2,3.... In this case, sin(kπ)/(kπ) = 0 and we get J3/2(kπ) = - (2/π2k)1/2 cos(kπ) = - (2/π2k)1/2(-1)k Just for the record, for k = 0 we can use Jackson page 72 which says J3/2(x) = (x/2)3/2 /Γ(ν+1) => J3/2(kπ) = 0 for k = 0. Then we have from above u00k(r,θ,φ) = [ a-3/2 / J3/2(kπ)] (2/π)1/2 (kπ r/a)-1/2 sin(kπ r/a) (4π)-1/2 = - [ a-3/2 (2/π2k)-1/2(-1)k] (2/π)1/2 (kπ r/a)-1/2 sin(kπ r/a) (4π)-1/2 = - [ a-3/2 (π2k/2)1/2(-1)k] (2/π)1/2 (kπ r/a)-1/2 sin(kπ r/a) (4π)-1/2 = - [a-3/2 (π2k)1/2(-1)k] (2/π)1/2 (a/kπr)1/2 sin(kπ r/a) (4π)-1/2 = - [a-3/2 π k1/2(-1)k] (2/π)1/2 (a/kπr)1/2 sin(kπ r/a) (4π)-1/2 = - [a-3/2 π (-1)k] (2/π)1/2 (a/πr)1/2 sin(kπ r/a) (1/4π)1/2 = - [a-3/2 (-1)k] (2)1/2 (a/r)1/2 sin(kπ r/a) (1/2)(1/π)1/2 = - [a-3/2 (-1)k] (a/r) sin(kπ r/a) (1/2) = - [a-3/2 (-1)k] (a/r) sin(kπ r/a) Now let's try our integral : !Syntax Error, Ir2dr ∫dΩ u00k(r,θ,φ)2 = (1/2π) !Syntax Error, Ir2dr ∫dΩ [a-3] (a/r)2 sin2(kπ r/a) = (1/2πa3) !Syntax Error, Ir2dr ∫dΩ (a/r)2 sin2(kπ r/a) = (1/2πa3) !Syntax Error, Ir2dr 4π (a/r)2 sin2(kπ r/a) = 2 !Syntax Error, Id(r/a) sin2(kπ r/a) = (2/π) !Syntax Error, Id(πr/a) sin2(kπ r/a) = (2/π) !Syntax Error, Idθ sin2(kθ) = (2/π) (π/2) = 1 // Schaum p 96 so we think we have our constant factors correct! Exercise 6.31. Eigenfunctions of Laplace for a cone of the unit sphere in 3D. How are we going to get things to vanish on the sides of the cone (centered on the z axis) ? We start by separating variables as in p 144 A. The angular equation is SY = μY. Well, the spherical harmonics do solve this equation even if we restrict to the cone, but they don't vanish on the cone! I tried various approaches and finally I think hit on the right one in Plan F. Plan A. I am thinking of scaling things like this: S(θ,φ) Ylm(θ,φ) = l(l+1) Ylm(θ,φ) S(θ,φ) = -∂θ(sinθ ∂θ) - (1/sinθ)∂φ2 Suppose we define θ' = θ(π/α) so when θ = α, we have θ' = π. Then θ = (α/π)θ' and we have S((α/π)θ',φ) Ylm((α/π)θ',φ) = l(l+1) Ylm((α/π)θ',φ) But this does not seem the right equation because we are really trying to solve with S(θ,φ) there. The sinθ factors are not going to scale properly. Plan B. So let's go to brute force. Separate as on page 394. The azimuthals are the same with integer m as shown on that page. The θ equation is also the same as in A.3. And then A.4 is also correct. We have backed off and we don't know λ at this point, it is no longer n(n+1). We don't care now about divergence at z = -1 because that is not part of our cone. Well equation A.4 is really replaced now with the range cosα < x < 1 and does not go down to -1. So we really have a new ODE here since we have changed the interval. Let me just assume that to avoid blowup at z = +1, we need λ = n(n+1). Then Pn|m|(x) is a solution, but now we have to form linear combinations! That is how you do it. So for a particular n, and λ, we look for solutions that do this: fn(x) = Σm=0∞anm Pn|m|(x) with fn(cosα) = Σm=0∞anm Pn|m|(cosα) = 0 To compute the anm we need the second item below. Then we get Σm=0∞anm Pn|m|(cosα) = 0 But we are going to end up with anm = 0! Plan C. Suppose instead we try this: fn(x) = { Pn|m|(x) – Pn|m|(cosα) } Is this constant a legal solution? Not for n ≠0. So that idea fails. Plan D. Let's go back to A.4 and try scaling the variable. I showed that this linear change of variable at least gets us to -1 < y < 1 x = (1-cosα)/2 y + (1+cosα)/2 But then the factors (1-x2) in A.4 get uglified and we no longer have the Assoc Leg equation. Plan E. Let's try to clearly state what our problem is: -∂x[ (1-x2)∂xΘ ] + m2/(1-x2) * Θ = λ Θ cosα < x < 1 BC's: Θ(1) = finite Θ(cosα) = 0 This does seem to be a well defined EV problem. In Stak volumes 1 and 2 so far, we really were not given any "general methods" for solving an EV problem. All I know is to do a power series around some point like x = 1. Plan F. Suppose we set λ= ν(ν+1) just as a change of variable. Then the ODE without regard to BC's is going to have solutions Pνm(x) and Qνm(x), we know that for sure. We would select the Pνm(x) presuming that the Q functions blow up somewhere, not sure where. Then our problem is to choose ν such that Pνm(cosα) = 0 That is to say, we want to find some ν such that Pνm(cosα) = 0! I could imagine this being possible, as you slowly vary ν the zeros all move, etc etc. Maple can plot these things, and yes the zeros move. I have never seen data on "the zeros in the index" of this function. _EnvLegendreCut := 1..infinity: plot(LegendreP(3,2,x), x= -1..1); // order here is l,m,z Suppose in this example we had cosα = .6. Then perhaps P5.343212(.6) = 0 and then this is an eigenvalue for ν. But perhaps there are several eigenvalues, several values of ν that causes a zero as shown. Maybe an infinite number. 9.20.09 Note Added. I just realized you can simply plot Plm(x) in Maple with l as the variable, and then the zeros are pretty obvious, including the one imagined above which turns out to be around l = 5.16. _EnvLegendreCut := 1..infinity: > plot(LegendreP(x,2,.6), x= 0..100); Well, Bateman does not mention any zeros, nor does AS. But GR has this claim: which is on a different page in my GR. The first item is the one. It suggests that there are zeros only for negative m, but Maple tells me we also have zeros for positive m. Still, we get the idea of an infinite number of zeros. As you increase ν, I think more nodes keep appearing and will sweep across your desired value. OK, so let's make a notation Pνm(cosα) = 0 νm,i = the zeros i = 1,2.3... index λm,i = νm,i(νm,i+1) That is to say, we have equation A.4 for a given m. We find there are eigenvalues as shown above, and no, in general the νm,i are NOT going to be integers! In the limit cosα → -1 then they will become integers somehow. The angular solutions to our problem are then these: eimφ Pν(mi)m(cosθ) m = 0, ±1, ±2 .... i = 1,2,3..... The radial solutions are then going to be as in p 145 H but with n = ν. So this is a pretty horrific result: um,i,k (r,θ,φ) = (1/r) Jν(mi)+1/2 [ β(ν(mi)+1/2k r] eimφ Pmν(mi)(cosθ) where m = 0, ±1, ±2 ...., and i = 1,2,3...... and k = 1,2,3.... The β(ν(mi)+1/2k are the zeros of Jν(mi)+1/2(x), but first you have to find the ν(mi) which are the zeros of Pν(mi)m(cosα) . No doubt I could make Maple find the first hundred zeros here, and for each of those, a hundred zeros of J. What a complete and utter mess! Just an exercise! Exercise 6.32. Eigenfunctions of Laplace for an orange slice of the unit sphere in 3D. (p 144) This solution starts out the same as 6.28 above (wedge in 2D, here we are 3D). Our azimuthal eigenfunctions are going to be Φm(φ) = sin(mπφ/α) μ = (mπ/α)2 where Φ is the function appearing on page 394 of Appendix A. So we get Φ"n(φ) = - μ Φ as shown. In that appendix we had μ = m2 so we are now going to make this replacement: m → m' = mπ/α . We end up then with equation A.4 (p 395) with this new m' installed, same unknown λ left from the radial/angular separation as in A.1 (but now the Y are as yet unknown as well). As in the last problem, we can redefine λ = ν(ν+1). The solution of A.4 is then given by this: Pνm'(cosθ) where ν is as yet unknown, and m' = mπ/α with m = 0, ±1, ±2 etc. I think, however, that now that we have the full range of cosθ, we have to have ν = l to avoid singularities at z = -1. I tried to mess with this in Maple, but find that LegendreP has "bugs" . It cannot evaluate close to z = -1 depending on the m argument, for example here evaluation fails: φ=α This is a long subject on which I will avoid digressing today. For non-integer m, things seem to work, and you do seem to get divergence at z = -1 when l is not an integer. Basically, this issue of l = integer is the same for our "orange slice" problem as for the full sphere problem! I just think Stak makes some misleading comments about z = +1 top of page 395. After all, if we look at a certain websites quote of the AS definition, assocLegendre1 := (mu,nu,z) -> 1/GAMMA(1-mu)*((z+1)/(z-1))^(mu/2)* > hypergeom([-nu, nu+1],[1-mu],(1-z)/2): it seems clear that at z = +1 there is no problem ?? Well OK, there is that outside factor, so I defer this whole topic. So, our "angular" solutions to our current problem are these: ( I set n = l ) Θn(θ) Φm(φ) = Pnm'(cosθ) sin(m'φ) n = 0,1,2,3.... m' = mπ/α with m = ±1, ±2 etc Thus, we have quantized λ to be n(n+1), and so I think our radial solution is the same as for the full sphere. Conclusion: my solution to this problem is given by page 145 H where we replace the spherical harmonic Y there with the above angular function! We still have three indices n, m, k as shown. I don't think there is any limit now on the range of m as there is when m = integer. Note added 3/25/10. First, I think the radial solutions are R = (1/r) Jn+1/2(r) were the are the zeros of these Bessel functions, as shown page 145 H. For the angular solution, I think we have to use this: (where J = 1,2,3....) Θn(θ) Φm(φ) = Pnm'(cosθ) sin(m'φ) n = -m' + J. m' = -mπ/α with m = 1,2,3... The reason is that Pnm'(z) blows up at z = 1 if m' > 0. See "Legendre SL Problem....". Here is an example: Exercise 6.33. Laplace on R where ∂nu = 0 (Neumann) on σ; apply to disk. (145) This is a Neumann instead of a Dirichlet situation. (a) show eigenvalues are real and positive. Page 137 C is still true because σ is not yet used. When we use first Green's on this, the surface term still vanishes because ∂nu = 0 on σ, so result D is still true! This tells us that λ ≥ 0, so EV are real and not negative. (b) What about λ= 0? From D we can see that u = constant is a candidate since grad(u) = 0. In the Dirichlet case, since u = 0 on σ, we got u = 0 inside so λ = 0 was not an EV. But in Neumann, this little argument is not true, so λ= 0 is an EV and u = constant is an EF. (c) solve the unit disk in this Neumann case. Our prototype here is Example 2 on page 140. Our azimuthal functions are Φn(φ) = einφ as top page 141 E. Radial equation is B there and sol is un = Jn(r). But now we do some math: un(r) = Jn(r) un'(r) = Jn'(r) un'(1) = Jn'() = 0 So now we need the zeros not of Jn, but of Jn' . Call these γ(n)k . You can relate J' to a sum of J's so no problem computing these things. So we have = γ(n)k and λ = [γ(n)k]2 k = 1,2,3 enumerator. If λ = 0 we have maybe a special case to worry about, I will ignore that. So here are the solutions: un,k(r,φ) = Jn(γ(n)k r) einφ n = 0,±1,±2 ... k = 1,2,3... λ = [γ(n)k]2 = EV's Now this was for r = 1. Just replace r on the RHS above with r/a to get result for radius a. Done. Exercise 6.34. Repeat last but with "mixed" BC; apply to the disk (145) (a) show EV are real and positive, and show orthogonality. The mixed BC here is hu + ∂nu = 0 on σ where h is some positive function. As before, we start with p 137 C. If we try to arrive at D, we have this surface term do deal with: ∂nu + h(σ)u = 0 ∫σ dS { ∂nu } = – ∫σ dS h |u|2 so that D is now replaced with λ |u|2 = ∫R |u|2 dx + ∫σ dS h |u|2 so now we have an extra term. Since h is positive, this term is positive, and we still conclude that λ ≥ 0. What about λ = 0? I think the second integral is > 0 for ANY non trivial u (given positive h), so even if we try the solution u = constant as our candidate λ= 0 solution, we fail. So λ > 0. Now look into orthogonality. Imagine writing bottom page 137 but with our mixed BC for each line. Note that we would have in the second line, ie + ∂n = 0 on σ. We then get to p 138A as before since BC's not used yet. We now have to show that our new surface term vanishes: surface term = ∫dS { u ∂n – ∂n u } = ∫dS { – u + hu } = ∫dS u {h-} When we say h is "positive" we imply of course that it is real, so {h-}= 0 and we still obtain p 138 B and we obtain orthogonality of EF's which have different EV's. (b) Now set h = constant and solve the unit disk. Angular part is the same and as before we get to the point un(r) = Jn(r) for the radial. But now we need (at r = 1) hu + ∂nu =0 hu = - ∂nu ∂nu / u = -h so un(r) = Jn(r) un'(r) = Jn'(r) Jn'() / Jn() = -h So this is some transcendental equation that yields eigenvalues = ρ(n)k with k = 1,2,3... and then our solution is un,k(r,φ) = Jn(ρ(n)k r) einφ n = 0,±1,±2 ... k = 1,2,3... λ = [ρ(n)k]2 = EV's 6.7 Methods for finding the Green's Function (146) Comment: This question has not really been addressed yet. We did this in Volume 1 for 1D situations: we found a solution on each side and matched them at x = ξ, for example. And earlier in Volume 2 we found some E functions (fundamental solutions) for various L operators, which solutions E were symmetric in r. But now we have the new issue of requiring that g = 0 on some strange boundary σ. Also in the last section we learned how to solve a Dirichlet problem once we know g. (1) The Integral Equation Method First, we do the trick we have done so many times before, to wit, [ think of ξ and t as constants for now ] -2g(x|ξ) = δ(x-ξ) g = 0 on σ // definitions -2E(x|t) = δ(x-t) no BC, "free space Green's" - E(x|t)2g(x|ξ) = E(x|t)δ(x-ξ) // process - g(x|ξ)2E(x|t) = g(x|ξ)δ(x-t) - E(x|t)2g(x|ξ) + g(x|ξ)2E(x|t) = E(x|t)δ(x-ξ) - g(x|ξ)δ(x-t) // subtract ∫R dx { g(x|ξ)2E(x|t) – E(x|t)2g(x|ξ) } = E(ξ|t) - g(t|ξ) // integrate over R, p 146 B LHS = ∫σ dSx { g(x|ξ)∂nE(x|t) – E(x|t)∂ng(x|ξ) } = – ∫σ dSx E(x|t)∂ng(x|ξ) g(t|ξ) – E(ξ|t) = + ∫σ dS E(x|t)∂ng(x|ξ) g(t|ξ) = E(ξ|t) + ∫σ dSx E(x|t)∂nxg(x|ξ) Where ∂n is wrt x, and the dS integral refers to the point x going over σ. Now we do the usual trick of first changing the names of variables: swap the meaning of x and t everywhere: g(x|ξ) = E(ξ|x) + ∫σ dSt E(t|x)∂ntg(t|ξ) Now use symmetry of g and E in three places (not four): g(x|ξ) = E(x|ξ) + ∫σ dSt E(x|t)∂ntg(ξ|t) Now we are supposed to define I(x|ξ) = – ∂nξg(x|ξ) I(ξ|t) = – ∂ntg(ξ|t) g(x|ξ) = E(x|ξ) – ∫σ dSt E(x|t) I(ξ|t) // which is 6.96 x and ξ are inside R and this verifies the claim earlier that if you know I, then you can in theory compute g. Also, comparison to earlier g(x|ξ) = E(x|ξ) + v(x|ξ) => v(x|ξ) = – ∫σ dSt E(x|t) I(ξ|t) Now we take x → s, a point on σ. We can think of I(ξ|t) as the density of a "simple layer" situation as in 6.32 p 114 where we had φ as our layer density. We learned there you can just take x → s right through the integral with nothing strange happening, so you get then g(s|ξ) = E(s|ξ) – ∫σ dSt E(s|t) I(ξ|t) s and t on σ, ξ inside R But g(s|ξ) = 0 since s is on σ, so we now have E(s|ξ) = ∫σ dSt E(s|t) I(ξ|t) s on σ t integrate on σ, ξ inside R // which is p 147 A. Now change variable names s,ξ,t → x,t, ξ : E(x|t) = ∫σ dSξ E(x|ξ) I(t|ξ) x on σ ξ integrate on σ, t inside R This is an aside: Now make the following identifications: (x on σ) E(x|ξ) = k(x,ξ) I(t|ξ) = u(ξ; t) E(x|t) = f(x; t) Then the above says: ∫σ dSξ E(x|ξ) I(t|ξ) = E(x|t) x on σ ξ integrate on σ, t inside R (**) ∫σ dSξ k(x,ξ) u(ξ; t) = f(x; t) x on σ ξ integrate on σ, t inside R where we think of "t" now as a fixed parameter. Look at Vol I page 195 and we see that this is an example of a First Kind Fredholm, but in n-1 dimensions, and where both k and f are derived from E. Now going backwards a bit, start with (**) and change x,t→ s,x ∫σ dSξ E(s|ξ) I(x|ξ) = E(s|x) s on σ ξ integrate on σ, x inside R // which is 6.97 Let's do our identification again E(s|ξ) = k(s,ξ) I(x|ξ) = u(ξ; x) E(s|x) = f(s; x) ∫σ dSξ k(s,ξ) u(ξ; x) = f(s; x) s on σ ξ integrate on σ, x inside R Again, our first kind Fredholm. We know E so we know kernel k and driving function f. So the integral equation is something we should solve for I . We then insert our solution I into 6.86 and we have g !!! Comment: In 1D, the variable x was just somewhere in our interval (a,b) and the integral was over this same interval. But this situation is NOT a limit of our general result above as best I can tell. The 1D limit of the equation above would be f(s; x) = [k(s,ξ) u(ξ; x)] ξ=bξ=a where s = a or b The 2D case might have x in a disk, and then σ would be the bounding circle. So we call this a "Fredholm First Kind" equation, but it is has a different sort of measure than the one on page 195 volume 1. Here we might have x in Rn but the integral is over an n-1 dimensional surface σ. Example: the unit circle. (147) In this example, the integral equation 6.97 is given by p 147 D which we want to solve for I. Notice that arbitrary point x is in fact placed at (r,0). I agree therefore that any function of x and ξ-on-circle is a function of r,ψ as shown, so we can apply that to I. Since I is the radial gradient of g magnitude (minus) I(r,ψ) = I(x|ξ) = – ∂ξn g(x|ξ) = – ∂r g(r,ψ) I agree it is the same at ψ and -ψ so I is even in ψ -- this from the symmetry of the problem. Thus, I agree with the Fourier series expansion on p 148 A. Let's jam that into D: ln(1+r2-2rcosθ) = ∫dψ ln[2-2cos(ψ-θ)] I(r,ψ) = Σm=0 Im(r) ∫dψ ln[2-2cos(ψ-θ)] cos(mψ) Now we rewrite the LHS, using "interesting sum fact 6.22" which I did derive back on page 104, LHS = ln(1+r2-2rcosθ) = -2 Σn=1 (rn/n) cos(nθ) RHS log = ln[2-2cos(ψ-θ)] = -2 Σn=1 (1/n) cos(n(ψ-θ)) This gives -2 Σn=1 (rn/n) cos(nθ) = Σm=0 Im(r) ∫dψ {-2 Σn=1 (1/n) cos(n(ψ-θ))} cos(mψ) = -2Σm=0 Im(r) Σn=1(1/n) !Syntax Error, Idψ cos(n(ψ-θ)) cos(mψ) What is this integral? Write cos(n(ψ-θ)) = cos(nψ)cos(nθ) + sin(nψ)sin(nθ) !Syntax Error, Idψ cos(nψ) cos(mψ) = δn,m π Schaum p 96 !Syntax Error, Idψ sin(nψ) cos(mψ) = 0 because the integrand is odd in ψ so we get -2 Σn=1 (rn/n) cos(nθ) = -2Σm=0 Im(r) Σn=1(1/n) δn,m π cos(nθ) = -2Σn=1(1/n) Σm=0 Im(r) δn,m π cos(nθ) = -2Σn=1(1/n) Σm=1 Im(r) δn,m π cos(nθ) = -2πΣn=1(1/n) In(r)cos(nθ) Comparison then tells us (rn/n) = π(1/n) In(r) => In(r) = rn/π n = 1,2,3... What about I0(r) ? This is not provided by our integral equation above. Instead we go to 6.84 which is true in general for I! In our expansion for I, only the constant term would survive this integral so we get ∫dψ I0 = 1 I0 = 1/2π So here is out complete solution for I I = 1/2π + 1/π Σn=1∞ rn cos(nψ) // which becomes p 149 B Now we get to use our pre-computed identity 6.21 p 103 to write Σn=1∞ rncos(nψ) = (r cosψ - r2) / (1 + r2-2r cosψ) and we then have I (r,ψ) = 1/2π + 1/π *(r cosψ - r2) / (1 + r2-2r cosψ) // which is p 148 G = 1/2π + 2/2π (r cosψ - r2)/ (1 + r2-2r cosψ) = 1/2π * [ 1 + 2(r cosψ - r2) / (1 + r2-2r cosψ) ] = 1/2π * [ 1 + r2-2r cosψ + 2r cosψ - 2r2 ] / (1 + r2-2r cosψ) = 1/2π * [ 1 - r2 ] / (1 + r2-2r cosψ) where ψ is the angle between ξ and x as shown page 147 figure, ξ is on circle, x is inside. If we were to put x instead at φ instead of 0 angle, then angle between ξ and x would be ψ-φ and we get result: I(x|ξ) = I(r,φ ; 1,ψ) = 1/2π *(1 - r2) / (1 + r2-2r cos(ψ-φ) ) // which is p 148 H. Now we jump back to page 135 where we had the general Poisson solution p 135 D. If we want to do a Laplace solution with no source, only the second term survives, and this is where 6.83 comes from on page 136. So if we stick in our I above we get u(r,φ) = ∫dψ f(ψ) 1/2π *(1 - r2) / (1 + r2-2r cos(ψ-φ) ) and we see that I is in fact that Poisson kernel thing in 6.11 (and not some charge density). [ But elsewhere I show it can be interpreted as a surface charge density in embedding metal if a point charge is put at point x = (r,φ). ] So in this example, we have used the "integral equation method" 6.97 together with various kick-ass identities to reproduce our Laplace-in-circle solution 6.11, a typical Stakgold tour de force. Now we want an explicit expression for g itself! I agree with p 148 J, then p 149 A and B. B is just a variant of what I showed above (see reference to p 149 B bold above). So now we have the following mess for 6.98: g = E - ∫dψ { 1/2π * Σn=1 rn cos(nψ-nφ) /n } {1/2π + 1/π* Σm=1 r0m cos(mψ-mφ0) } The isolated 1/2π term contributes nothing because ∫dψ cos(nψ-nφ) = 0 for n = 1,2,3... This leaves g = E - 1/2π2 * Σn=1 rn/n Σm=1 r0m∫dψ cos(nψ-nφ) cos(mψ-mφ0) We can shift the integration variable ψ' = ψ-φ0 to get n(ψ-φ) = n(ψ'+φ0- φ) so we have ∫dψ' cos n(ψ'+φ0- φ) cos(mψ') ∫dψ cos n(ψ-Δφ) cos(mψ) Δφ = φ-φ0 cos(n(ψ-Δφ)) = cos(nψ)cos(nΔφ) + sin(nψ)sin(nΔφ) so again we are left with this for the integral cos(nΔφ) π δm,n and so then g = E - 1/2π2 * Σn=1 rn/n Σm=1 r0m∫dψ cos(nψ-nφ) cos(mψ-mφ0) = E - 1/2π2 * Σn=1 rn/n Σm=1 r0m{ cos(nΔφ) π δm,n } = E - 1/2π * Σn=1 rn/n r0n{ cos(nΔφ) } = E - 1/2π * Σn=1 (rr0)n/n * cos(nΔφ) and once again we invoke 6.22 p 104 to write this as = E - 1/2π * { - 1/2 ln [ 1 + r2 r02 – 2r r0 cos(φ-φ0)] } = E +1/4π * ln [ 1 + r2 r02 – 2r r0 cos(φ-φ0)] // which is p 149 D = 6.99 So this second term is really v(x,ξ) which we add to E to get g. We have our explicit answer for g on the unit circle (which implies g = 0 on the unit circle). Then in steps E,F,G,H we are able to write this Green's solution as the sum of E from our source, E from an image source, and a constant, as shown in result H. The image source is at 1/r0 (outside the unit circle) whereas the real source is at r0. Contradiction? As usual, I am confused. If I take two sources of opposite sign , the positive one at the origin of a coordinate system, I would write (from p 50 5.105) -V = 1/2π * ln r - 1/2π * ln R where r2 = x2 + y2 R2 = (x-d)2 + y2 where the negative source is located on the x axis at position (d,0). So -2πV = ln(r/R) For V to be a constant on some σ, we need r/R = constant on σ, or r2 = KR2. So we need x2 + y2 = K [(x-d)2 + y2] x2 + y2 = Kx2 + Kd2 + Ky2 -2Kxd x2(1-K) + y2(1-K) + 2Kxd = Kd2 For example, if K = 1, this becomes 2xd = d2 or 2x = d or x = d/2 which is the well known vertical line between the two charges. In the limit K = 0 we get x2 + y2 = 0 which is a circle of radius 0. Define a = K/(1-K). Divide by 1-K to get x2 + y2 + 2axd = ad2 Now complete the square x2 + 2axd + (ad)2 - (ad)2 + y2 = ad2 (x + ad)2 + y2 = ad2 + a2d2 = ad2(1+a) = ρ2 And indeed, this is a circle centered at (-ad,0) with radius ρ2 = ad2(1+a), ρ = d . Here is our picture: Now in the limit K = +∞ we get a = -1 and this ρ = 0, so we get a circle of zero radius centered at the origin. For K = -∞ we get a = +1 and circle is then centered at x = -d, not at origin. So my confusion is this: how on earth can the equipotential surface be a circle centered at the positive charge as this Stakgold discussion is claiming? [ but that is NOT what he is claiming! ] Resolution: In my work above, I have the origin O located at the positive charge. My coordinates were shall we say (r,θ) or (x,y) relative to this origin. But Stak is working instead with an origin located at the center of the above circle and the +1 charge is located at point ξ which is distance r0 = ad to the right of this center (and his image charge is located point ξ*, a distance 1/r0 to the right of this center.) So let's call his origin O' with coordinates either (r',φ') or (x',y'), We then know that x' = x + ad y' = y So let's rewrite my locus equation above in these coordinates: (x + ad)2 + y2 = ρ2 ρ2 = ad2(1+a) my origin x'2 + y'2 = ρ2 ρ2 = ad2(1+a) Stak origin Regardless of origin choice, we see the main interesting fact: in this 2D dipole problem, the equipotential curves are all circles! As we adjust parameter K and hence a, the origin of the circle moves relative to the +1 charge, and the radius changes, but it remains a circle. Now what happens if we put the -1 charge at distance ad+d = 1/r0? We can see that d = 1/r0 - r0 . We can adjust our K and hence our a to get ρ = 1 and then we have a unit circle. And on this circle we have -2πV = ln(r/R) r2 = x2 + y2 R2 = (x-d)2 + y2 Since this V is constant on the circle, we can just pick some nice point on the circle for evaluation. Let's pick then the intersection of the circle with the x axis. There we have (my coord system) r = ρ – ad R = d-r = d - ρ + ad -2πV = ln(r/R) = ln[ (ρ-ad) / (d - ρ + ad) ] But we know that ρ-ad = 1 - r0 d - ρ + ad = (1/r0 - r0) - 1 + r0 = 1/r0 - 1 So we get -2πV = ln(r/R) = ln[ (1 - r0) / (1/r0 - 1) ] = ln [r0(1 - r0) / (1-r0)] = ln r0 and hence V = - 1/2π * lnr0 on this unit circle. If we then add the term + 1/2π * lnr0 as shown in 6.100, we cause V = 0 on this unit circle, and then we have arrived at the Stakgold 6.100 including its image charge interpretation. So the contradiction has now gone away. Plots of the 2D equipotential lines. r := ( x^2 + y^2)^(1/2): R := ( (x-1)^2 + y^2)^(1/2): > v := -ln(r/R); with(plots):implicitplot({seq(v=N/2,N=1..6)},x=-0.6..0.3, y=-0.5..0.5, grid = [200,200], scaling=CONSTRAINED); So here you see some examples of the perfectly circular equipotentials. For the record, Maple can also do this with the following call, and the curves are a little smoother as shown on the right. with(plots):contourplot(v, x=-0.6..0.3, y=-0.5..0.5,contours = 10,grid = [200,200], scaling=CONSTRAINED); Comment on the 3D dipole equipotential lines. Let's consider (+1) and (-1) unit charges separated by d, and look only in the z=0 plane which then matches our earlier picture but now we have a slice through a 3D surface. We have drawn a sphere of radius ρ here, but ignore it! Then we have 4πV = 1/r – 1/R where r2 = (x2 + y2) R2 = [(x-d)2 + y2] so that (4πV)2 = 1/r2 + 1/R2 - 2/(rR) // ±V are now both included due to squaring (4πV)2 - 1/r2 - 1/R2 = -2/(rR) r2R2(4πV)2 -R2- r2 = -2rR [r2R2(4πV)2 -R2- r2]2 = 4r2R2 // ±r are now both included due to squaring [r2R2v2 -R2- r2]2 – 4r2R2 = 0 v = 4πV Let's have Maple write this out in terms of x and y. I defined symbols a,b,c,e as shown below, and then the equation following is a restatement of the last equation above, and I define the expression to be f: a ≡ r2 b ≡ R2 c ≡ v2 e ≡abc f ≡ [e -b- a]2 – 4ab = 0 Maple code is then as follows. I set up all the symbols and then asked Maple to expand everything so I could see the equation of the curve in Cartesian coordinates. restart; a := x^2 + y^2; b := (x-d)^2 + y^2; c := v^2; e := a*b*c; f := (e - b - a)^2 - 4*a*b; g := expand(f); h := sort(g,[x,y]); We then find that the LHS or our equation LHS = 0 is this function 0 = h(x,y; v,d) which gives us a very horrible 8th degree polynomial in x and y for our equipotential surface. Our equation h(x,y) = 0 in fact includes the solutions to four problems, not just one problem! We can summarize those four problems as follows: ±4πV = ±1/r – 1/R 4πV = 1/r – 1/R + – // our intended problem inner left curve 4πV = –1/r + 1/R – + // polarities reversed inner right curve 4πV = 1/r + 1/R + + // both charges positive outer left curve 4πV = – 1/r – 1/R – – // both charges negative outer right curve For this reason, plotting the above 8th order polynomial mixes all the four problems together and the plots are not so useful. Here is an example: v = 3.75 I just wanted to show that we really do have 8th order curves here. If we are interested in plotting, it is better to go back and do it this way: (here I set d = 1 again) > restart; > r := ( x^2 + y^2)^(1/2): > R := ( (x-1)^2 + y^2)^(1/2): > v := 1/r - 1/R: > with(plots):implicitplot({seq(v=N/4,N=1..8)},x=-2..2, y=-2..2, grid = [200,200], scaling=CONSTRAINED); Notice the use of the "seq" trick to get 8 plots on a single graph. Here you see first of all the strange non-circular shape of these 8th order curves. The innermost curve has v = 2.00, the outermost curve has v = 0.25. The innermost curves are "far" from the -1 charge and are pretty close to being circles. On the left side, for the outer curves, the -1 charge at x = 1 has little effect, and we have almost circles there as if there were no -1 charge. But on the right, the circles are compressed toward the +1 charge for this reason: Because of the larger negative contribution from the -1 charge, we have to move closer to the +1 charge to maintain a given value of v. Just for fun, here are the same plots if both charges are positive: v := 1/r +1/R: with(plots):implicitplot({seq(v=N,N=1..6)},x=-2..4, y=-2..2, grid = [200,200], scaling=CONSTRAINED); Here you see for large v that each charge has its own curve, but around v = 4 these curves merge into a single curve, and of course for small v these will become large circles. Now back to Stakgold 6.100 on page 149: The summary is that we are trying to find the Green's function inside the unit circle, meaning g = 0 on that circle. We find that the solution is the E function for the positive charge, plus another E function for a -1 charge placed at a certain "image point" location, and then a constant added. If we knew ahead of time that the equipotentials for a pair of 2D charges were circles, it might have occurred to us to look for a solution by adding an image charge and adjusting its position so as to get an equipotential circle of radius 1. The image charge is of course not really there. It's just that the 2-charge problem gives the same solution within the unit circle as our Green's function problem for g. That is to say, within the unit circle if we put our unit + charge at some location ξ, and examine the potential at position x (also inside the unit circle), that potential is given by g(x|ξ) where we require that g = 0 on the circle. Anticipation of the next section: Suppose we scale the negative charge by k = 1/(ad) = 1/r0. This is the same as replacing 1/R by k/R = 1/(adR). Probably k is really ρ/ad where ρ is the radius of your desired sphere. So with ρ = 1, we replace R → adR (2) The Method of Images You basically guess a set of image charges and their positions and sizes, and then show that your guess gives the correct results. Examples Example 1: Green's function for the unit sphere (3D) (150) Earlier we did this in 2D with the result being 6.100 where we found an image point charge -1 located at a certain inverse point. We did that using the horrendous "integral equation method". Now we are going to do here the analogous 3D situation. We have a charge inside our sphere at position ξ, and we guess an image charge of size A sitting outside that ξ* where |ξ*| = 1/|ξ|. Then our v(x,ξ) is as shown in C. We then quickly arrive at F which is our candidate solution to this problem: the potential at x due to a point charge at ξ and vanishing on surface σ. We made it vanish by construction in D for any point s on σ, so we know that our guess was right. Now just to convince us even more, he computes our friend ∂nξg(x|ξ) appears in our general equation u(x) = ∫R dξ g(x|ξ) q(ξ) – ∫σ dSξ f(ξ) ∂ξn g(x|ξ) // q = 0 in source-less Dirichlet so we suspect ∂nξg(x|ξ) [ ≡ -I(x|ξ) ] will be minus the Poisson kernel as shown in our 3D result way back in 6.26. He laboriously computes this derivative and in 6.103 it all falls out that the resulting I(x|ξ) is exactly the Poisson kernel in 6.26. I am skipping the math detail here because I don't think it contains anything new for me. Comments: So the solution here is sort of 4πV = 1/Q - k/R. For this problem, I want to use a different picture where now the sphere is centered at the origin and we are looking at the z=0 slice, and where I am now calling the one of the distances Q instead of r as used formerly, because r has a new meaning. To adapt this to our 3D case, imagine the +1 charge point to be ξ = (r0,0), the point on the sphere to be r, and the negative charge out at ξ* = (d+r0,0) *. In order to get 0 potential on the spherical surface, we need to make two choices (1) set d so that d+r0 = ρ2/r0 (2) scale the negative charge by factor ρ/r0. This means that our image charge has size - ρ/r0 and is located at (ρ2/r0, 0). Since we want V = 0 on the sphere, our equation is simply 1/Q = k/R where k = ρ/r0. We don't need Maple to do the algebra which shows that the V=0 surface is in fact the sphere or radius ρ shown. 1/Q = 1/R' = > Q = R' = (r0/ρ) R => Q2 = (r0/ρ)2 R2 which says [(x-r0)2 + y2] = (r0/ρ)2 [(x - ρ2/r0)2 + y2] so it is not hard to imagine we might get a circle. We don't need Maple for this: ro2[(x- ρ2/r0)2 + y2] = ρ2[(x-r0)2 + y2] ro2[x2 + ρ4/r02 - 2xρ2/r0 + y2] = (ρ2x2 + ρ2r02- ρ22xr0 + ρ2y2) ro2x2 + ρ4 - 2x ρ2r0 + r02y2 = ρ2x2 + ρ2r02- ρ22xr0 + ρ2y2 ro2x2 + ρ4 + r02y2 = ρ2x2 + ρ2r02 + ρ2y2 x2(ρ2 - ro2) + y2(ρ2 - ro2) = ρ4 - ρ2r02 = ρ2(ρ2-r02) x2 + y2 = ρ2 Now, suppose we throw in z2 in all the right places. That means we add +z2 in both opening []. This means that z2 carries through just as does y2 and we will then get this result x2 + y2 + z2 = ρ2 Conclusion: Draw a sphere of radius ρ centered at the origin. A distance r0 < ρ out on the x axis put a charge +1. A distance ρ2/r0 out on the x axis put an image charge – (ρ/r0). Let v be the sum of the potentials of these two charges. One will then find that v = 0 on the sphere! Therefore, the sum of these two potentials is the Green's Function for the sphere. This is the content of equation 6.102 where the ratio 1/|ξ| is our ρ/r0, since he has ρ = 1. Example 2: Green's function for a half space (3D) Here we have our +1 charge sitting somewhere in the right half space at position ξ. The distance of this charge from our grounded plane is ξ1> 0. The "region" R for this problem is the half great sphere on the right. Our BC is that g = 0 on this entire boundary, which includes the plane and the great sphere. The image charge of course goes at ξ* = (-ξ1, ξ2, ξ3) and has value -1. Now, what can we calculate in this situation? Now, since we have guessed what g is, we can compute ∂nξg(x|ξ) on the plane. This is of interest for two reasons. First, it tells us the surface charge density on the metal plane which I have been calling Σ (if we were to put a metal plane down the middle and a point charge at ξ: Σx(ξ) = ∂(ξ)n g(x|ξ) ≡ – I(x|ξ) ). This result is shown in p 152 D where he has now put the charge out at x, and ξ is a point on the plane. He makes this confusing swap x↔ξ because I(x|ξ) is defined as shown. We can see that this distribution peaks up on the plane at a point directly under the charge (but is everywhere negative and integrates to -1) plot3d(1/(1 + x^2 + y^2 )^(3/2), x=-4..4, y=-4..4); The second reason we want to know I(x|ξ) is that we can the jam I into 6.83 u(x) = ∫R dξ g(x|ξ) q(ξ) + ∫σ dSξ f(ξ) I(x|ξ) second part only to find the solution to this Laplace problem when an arbitrary potential f(ξ2, ξ3) is specified on the plane. Recall that I(x|ξ) is always the "Poisson kernel" for a problem on some region R with some boundary σ on which g = 0. Thus we have the double integral p 152 E for our general solution u. Then using a limit exercise we did earlier in the book (page 15 5.2), we can take the limit x1→ 0 and show that u → f on the plane, as shown in G. I am pretty sure this double δ thing is exactly the same limit we did more recently on page 115, (6.33), he could have used that as a more recent reference. So yes, the half-space bounded by a plane is a good "standard toy problem" to work with. I have already displayed the equipotentials for this "dipole problem". You get those oblate things around each charge, and in the limit they becomes the plane boundary so that g = 0 on this plane. We know that the gradient -u of our potential is always perpendicular to the equipotential lines, so in this problem the E field at the plane is perp to the plane and proportional to the Σ we computed above. This picture shows the equipotentials in the z = 0 slice if the charge is placed at the origin and the vertical plane is placed at x = 0.5 (not shown). It just happens that my picture is this way, with the region of interest on the left half of the plane instead of the right. Example 3: Neumann's function for a half space (3D) Same problem but Stak refers to the Green's function as h instead of g, and we want ∂nh = 0 on the plane instead of g = 0 as in the last problem. If we guess that the image charge now is in the same location but has +1 charge, what do we see? The form for h is indeed "even in x1" with the plane at x1= 0. Any function f(x) that is even has f'(0) = 0, which is just what we want here ( = 1). Thus, h' = 0 on the plane. If we were to draw our equipotential lines for this problem, they are still those distorted circle things, but now both sets of rings have the same V, not opposite V. This means that we have the same E field arrows as well, except the arrows due to the image charge are all opposite direction. The upshot is that the total E field at the plane is parallel to the plane (no longer a metal plane!) . That is what you expect if ∂nh = 0 at the plane, since ∂nh is the E field perp to the plane. To match the above, our Neumann problem equipotential lines would be the left half of our earlier picture where zero slope at the dividing plane says ∂nh = 0. Example 4: Line charge parallel to a plane (3D) You would first write this out in 3D cylindrical coordinates with z parallel to the line charge. Symmetry tells us that every slice z = z1 will have the same solution. In the remaining coordinates ρ,φ the Laplacian is the same as a 2D one in 2D polar coordinates. If we want g = 0 on the plane, we "guess" that we put an image line charge over on the left, opposite charge, same distance from the plane. This is what 6.106 says. We could then compute I(x|ξ) following the steps of example 2, and he claims that the resulting "Poisson kernel" is as shown in 6.107 (without proof, yet) and there he is integrating this kernel against an assumed potential prescribed value f along the plane in a direction perp to the line connecting the two line charges. As proof, suppose we started with 6.105 where direction ξ3 would be along the line charge. If we prescribed that f did not depend on ξ3, we would have a 3D problem with two point charges and an f that does not vary in this direction. The solution to this problem in a plane through the two point charges seems to be the same as the solution to the line charge problem in any slicing plane. I would have to ponder this for a while, to accept that the two problems have the same solution, though it does seem reasonable. // OK, I am sure that is correct. Since I have considered the "dipole in 2D" above already, I know what the potential lines will look like on the active side of the plane: they are all circles! (3) The Method of Full Eigenfunction Expansion (153) First, a little volume 1 review. On page 259 we looked at the driven equation with L = -∂x2-λI, and of course we were in a 1D world. We had the EV problem Lφ = 0 and we found some φn and λn. Then on page 261 we arrived at a certain "bilinear series" for the Green's function 4.8 with λn-λ in the denominator. We regarded λ as a parameter of g and obtained various fancy results. We even thought about analyticity in λ of g, for example. Now fast forward to the current context where now L = -2 so in effect we have λ = 0. We expect to get our bilinear series with λ = 0, and that is what 6.108 says, but now we are in a 3D world (or nD). Stak provides us with a very quick derivation of 6.108 in steps A,B,C,D. In A, we have just taken our Green's defining equation, multiplied by an eigenfunction un* and integrated. In B we are able to move the 2 to the other side because both g and un vanish on the bounding surface -- this uses Green's Theorem. But then of course we know 2 on un and this brings in λn. So result C is A + B + our knowledge of 2 on un. Then in D, we think of expanding g on the complete set un and all we need is the expansion Fourier coefficients. But those are provided by C, all done! [ Notice that this is a multidimensional Fourier Expansion on the volume R.] Recall our earlier finding that for the Laplacian, all λn are positive, so we don't have to worry about a divide by 0 situation here. You can see that EF's with the smallest eigenvalues have the largest weight in g! This is really quite an elegant formula I think, this 6.108. On page 154 top he quotes this eigenfunction expansion for two simple problems: the rectangle a,b and the circle r0. These results look extremely reasonable. You of course have to have normalized φn and you have to know the λn and we did all that earlier for these two problems. The reason it is sin sin is that you have g = 0 on all of the boundary, by the way. You don't have sin sinh, for example. (4) The Method of Partial Eigenfunction Expansion (154) This section opens with a very clear description of the generalized 2D "curvilinear quadrilateral" problem of the type I encountered in many earlier exercises. We assume our problem is separable in the two dimensions of this quadrilateral, as would be indicated by a possible writing of the total PDE as in p 154 A, and the trial solution as in B. If we set either factor to a constant like λ or μ, we have a 1D eigenvalue problem going across its direction of the quadrilateral. The idea is to pick one or the other of these two eigenvalue sets and regard those eigenfunctions as a complete set as your starting point. This is why this method is called "partial EF expansion" because we only do EF's in one of the two variables. Comment on Separation of Variables. In the text here, Stak says that this whole method only works if the homogeneous operator L(s,t) in the two variables s and t "separates" as shown in p 154 A as the sum of two operator terms, each involving only its own variable. Only then can you have a "separation constant" λ and only then do you get the eigenfunction systems A and B shown on page 155. What Stak fails to say is what happens when you throw in that delta function when you start talking about a Green's Function solution. He shows this in three examples, but never explains the general case, so let's do that here. At this point we have something like this for our full Green's function system: { (1/a1(s))L1(s)+ (1/a2(t))L2(t) } g(s,t) = f(s0, t0) δ(s-s0) δ(t-t0) (*) where the "delta function" separates into a product of φ functions as shown, which will always be the case no matter what coordinates we pick, sort of a Jacobian rule I guess. This does NOT tell you that g(s,t) is separable! Suppose we pick the #1 system and have (1/a1(s))L1(s)φλ(s) = λ φλ(s) and we go off and solve this EV problem to get the λ and φλ(s) functions. Then let's TRY this form of a solution for g, ie, let's try "expanding on these EF's φλ" : g(s,t) = Σλ gλ(t) φλ(s) If we insert this into (*) above we get Σλ gλ(t) λ φλ(s) + Σλ φλ(s) [(1/a2(t))L2(t) gλ(t)] = f(s0, t0) δ(s-s0) δ(t-t0) (**) Now we know that from the φλ EF problem there will get some "expansion of the delta function" that looks like this: (assuming operator L1 is Hermitian or whatever is needed) (1/2πi) dλ g(x|ξ; λ) = – δ(x-ξ)/s(x) = – Σn φn(x) n(ξ) – ∫dν φν(x) ν(ξ) // (4.95) which in our current case will just write symbolically like this: δ(s-s0)/s(s0) = Σλ φλ(s) λ(s0) If we insert this into the RHS of (**) we get Σλ gλ(t) λ φλ(s) + Σλ φλ(s) [(1/a2(t))L2(t) gλ(t)] = f(s0, t0) δ(t-t0) Σλ φλ(s) λ(s0)/s(s0) But since we know φλ(s) is a complete set, we set its "coefficient" to zero: gλ(t) λ + [(1/a2(t))L2(t) gλ(t)] = f(s0, t0) δ(t-t0) λ(s0) /s(s0) We can rewrite this as: [(1/a2(t))L2(t) + λ ] gλ(t) = { f(s0, t0) λ(s0) /s(s0)} δ(t-t0) = K(s0, t0, λ) δ(t-t0) So now we see the MAIN POINT: you end up with Green's Function system in the other variable t. The numbers s0, t0, λ are just some constants here. So you have reduced things to a 1D Green's Function system in the other variable. This only works if the original total L = { (1/a1(s))L1(s)+ (1/a2(t))L2(t) } was separable as shown in {...}. Long Comment on Example 1 (rectangle in 2D) to come: In my earlier work on such quadrilateral problems looking for u, we had a 0,0 BV situation in one dimension, and a 0,f situation in the other direction. So we used an EV problem in the 0,0 direction with some μ eigenvalue, and of course that same EV appears in the 0,f equation. But, for a simple L = 2 in 2D, this meant that solutions were "sine in one direction and sinh in the other direction". The sinh was good enough to satisfy the 0 of the 0,f. See Exercise 6.13 above in these notes. If we try to apply that method to our 0,0 and 0,0 problem, we get the 0,0 and 0,f result which is this: (from above) u(x1, x2) = Σn=1∞ bn [ sinh(nπx2/a)/ sinh(nπb/a)] sin(nπx1/a) bn = (2/a) !Syntax Error, Idx1 f(x1) sin(nπx1/a) But now if we set f = 0, we get bn = 0 and then u = 0, so this method only produces the trivial solution of the rectangle problem! Fascinating. We know there is a non-trivial solution g, but this method is not finding it! I wonder why not? We know that the sine part is correct in the first variable. Here is my attempt to explain this, again looking at Exercise 6.13 above. There, we assume that the solution had this form: u = Σn pn(x1)qn(x2) pn= sin(σnx) qn = sh(σnx) This "form" is a sum of terms such that each term is the product of two functions such that one function solves the first ODE with σn and the second function solves the second ODE with this same σn. For f(θ,φ) ≠0 this gave a non-trivial solution, and since the solution is unique, we found the solution and solved the problem and we were done. But for f = 0, this FORM of the solution must not be valid. If we look at the correct g solution in 6.109, for example, it basically has this form g(r | r') = Σnm φnm(r) nm(r')/λnm φnm(r) = pn(x) qm(y) This is of course a "complete form" because it expands over the eigenfunction set in each direction. There are lots of terms where the two functions have different eigenvalues n ≠m. Now, in the Example 1 below, we are going to call the variables x and y which I will combine as r. We will start off with fm(x) as above with the same sin σm as above. We will then expand [ let r' = (ξ,η) ] g(r | r') = Σm gm(y; r') fm(x) . // 6.113 and p 155 C We will then show that gm(y; r') must solve a certain Green's equation p 156A (containing both parts of r') which we solve by Volume 1 methods to get p 156 E. It is the y part of r' (η) which plays the usual Green's source point role, and ξ is just tagging along as a parameter. Notice things like y< = min(y,η). So the form of our solution here is a single sum of the product of a regular 1D eigenfunction fm(x) in one of the variables, and then a 1D Green's function for the other variable gm(y; r'). One nice feature is that the result is a single sum and not a double sum. [ I will clarify this method on the next section in my own notation.] On page 157 Stak talks about another method where we do as above to get gm(y; r'), and then we expand that thing on qn(x) as well to get g(r | r') = Σmn gmn( r') pm(x) qn(y) This is just a double Fourier projection of our problem. Comparison with our original double expansion then shows that gmn( r') = pn(x') qm(y') / λnm which is exactly what we see in p 157 C. If I have understood this correctly, I do not understand why Stak says "which is not typical" because it would seem that this method always gives the same double sum as the generic double expansion. Review of the Actual Method We try to find a solution of the form ( r = (x,y) ) g(r | r') = Σm gm(y; r') fm(x) where gm(y; r') = ∫dx g(r | r') fm(x) Here we have just done a 1D Fourier analysis on one of the two variables of r (namely we have replaced the variable x with discrete index m) The next step is to "process" the original PDE which is -2g(r|r') = δ(r-r') . Here is that processing: (a) multiply the PDE by fm(x) and then integrate over x over its natural range for its "direction". For the rectangle that would be from 0 to a. This gives -∂x2g - ∂y2g = δ(x-x') δ(y-y') - fm(x)∂x2g - fm(x)∂y2g = fm(x)δ(x-x') δ(y-y') - ∫dx fm(x)∂x2g - ∫dx fm(x)∂y2g = ∫dx fm(x)δ(x-x') δ(y-y') = fm(x') δ(y-y') (b) Now do parts once on the first term - ∫dx fm(x)∂x2g = + ∫dx ∂x fm(x)∂xg + [fm(x) ∂xg] |endspoints = ∫dx ∂x fm(x)∂xg where the "parts" vanishes because fm(x) vanishes at both endpoints. Now do a second parts to get - ∫dx fm(x)∂x2g = ∫dx ∂x fm(x)∂xg = – ∫dx ∂2x fm(x)g + [∂x fm(x) g]|endpoints = – ∫dx ∂2x fm(x)g where this time it is g that vanishes at the end points. But we know that fm(x) is an EF of -∂2x in its dimension, so we can replace – ∂2x fm(x) = λm fm(x). Then our first term has become – ∫dx ∂2x fm(x)g = – λm ∫dx fm(x) g(r|r') = – λm gm(y; r') (c) We write the second term as - ∂y2∫dx fm(x)g(r|r') = - ∂y2 gm(y; r'). (d) Therefore our processed equation has this form: – λm gm(y; r') - ∂y2 gm(y; r') = fm(x') δ(y-y') // which is p 156 A This is just a 1D Green's equation with an extra multiplier function fm(x') on the right which we know (in theory) how to solve from our Volume I efforts. Example 1. Green's function for the rectangle. (155) We are working in 2D, σ the boundary of a rectangle a,b as usual. If we pick the (x,a) problem as our 1D EV problem, we get sin(mπx/a) as the φm(x). So we expand our solution Green's g on this as shown in p 155 C where gm(y) is an as-yet unknown coefficient which depends on y as a parameter. Since this is your basic sine transform, we can write the inverse as in 6.112. But at this point we know neither g nor gm. We follow instructions to process the full PDE to get result p 155D. Let's examine the two terms on the LHS of D. -(2/a)∫dx sin(mπx/a) ∂x2g(x|y) = He says to do a double parts on this, and I can call upon my rule ∫dV ψ(∂i2φ) = ∫dV (∂i2ψ) φ + ∫dS J where J = {ψ (∂iφ) – (∂iψ) φ } where i = x, ψ = sin(mπx/a) and φ = g(x|y) . The surface term will then have J = { sin(mπx/a) (∂i g(x|y)) – (∂i sin(mπx/a)) g(x|y) } but for any fixed y, both sin(mπx/a) = 0 and g(x|y) = 0 at x = the two endpoints, so surface term = 0. We then have -(2/a)∫dx sin(mπx/a) ∂x2g(x|y) = -(2/a) ∫dx ∂x2sin(mπx/a) g(x|y) = +(2/a)(mπ/a)2∫dx sin(mπx/a) g(x|y) = (mπ/a)2 gm(y) The second term is trivially -∂y2gm and the RHS stays as is, so we arrive at p 156 A. But this equation is an ODE Green's type equation for the function gm(y) where η is just some fixed point in the delta. Looking at 6.112 where y is a parameter, we know g(x|y) = 0 at y=0 and b so gm(0) = gm(b) = 0 and we have here a 1D Green's equation therefore, with a constant thing on the RHS. Both ξ and η are fixed values, the variable in this ODE is y. So we now go off and use our Volume 1 methods to solve for the Green's function gm(y| η; ξ). You do this by going away from y = η so you have the homo equation but its solution is the sinh type thing not the sin thing due to the signs shown. On each "side" of η, we choose a solution which vanishes at that end, so C looks reasonable. Question: why is this continuous at y = η? Consider the point y = η-ε. At that point the solution has the form sh(η-ε)sh(b-η). But at the point y = η+ε, we have sh(η)sh(b-[η+ε]). As you let ε→0, these two are indeed the same, hence continuity. I will now skip some math and assume he has gotten the Green's correctly. The usual g' jump condition nails the constant C. It does seem odd to have ξ always hanging around as a parameter. The variables are y and η of this thing gm(y|η). So this little Green's is given in p156 E and then we have our answer for g in 6.113. Well, the solution of 6.111 really has the form g(x,y|ξ,η) where we are showing the coordinates in 2D of the two "points". So fine, we expect the result in 6.113 to contain ξ. Having done all these gyrations, one wonders why this result is useful. Well, for one thing it is a single sum and not a double sum (note that 6.108 is symbolic for what here is a double sum). Also, in the double some "complete form" we have sine in both variables, whereas here we have sinh in y. Maybe that is helpful. The answer is going to lie in how many terms of the series you have to compute to get a good answer! Suppose |y-η| is relatively large, and to be explicit, suppose y << η. Then have a...y......η...b gm(y) ~ sh(mπy/a)sh(mπ[b-η]/a) ~ emπy/a [emπb/a e-mπη/a - e-mπb/a emπη/a ] ~ emπy/a emπb/a e-mπη/a ~ emπb/a emπ(y-η)/a Now throw in the denominator sh(mπb/a) ~ emπb/a and you get gm(y) ~ emπ(y-η)/a ~ e-mπ(η-y)/a = e-mπ|η-y|/a in my case η> y. Then sequential series terms will have ratio e-π|η-y|/a which is the claim he makes in p 156 F. So if |η-y| > a say (ball park), then series is converging fast, maybe need only the first term (or two). He then states the same solution where we start with the OTHER complete EF set, the sin(nπy/b). Then results 6.114 and p 157 A are obtained by just doing x↔y, a↔b, ξ↔η. Then we can use this form for a good expansion if it happens that |x-ξ| is large. This does seem pretty obscure, but maybe we can get an example. Now p 157 sees a slight subject change. Suppose in 6.113 we expand our friend gm(y) on the y EF set. Then we define those coefficients as gm,n and we have the transform B. Not surprisingly, the double projection completely diagonalizes the PDE giving C which is an explicit expression for the numbers gm,n. This takes us to the double sum 6.115 which not only looks like our old double-EF expansion, but in fact is identical to 6.109. But he says that in general you tend to get a double sum which is not the same as your double EF expansion, so one of the sums must have some functions that are not EF's. As noted above, this mystifies me. I would think there would be a unique double sum that is the solution. We continue now mid page 157. How can we use our newly found Green's function g to solve our old problem "the Dirichlet problem on the rectangle with f on the top edge"? He pulls out the combination of D and E which I agree with, then he computes I as -∂ηg assuming η > y which seems to be good enough. In D we are integrating only along the top edge of the rectangle where f ≠ 0. So we compute our normal derivative of g to get p 158 B (but he left out a factor which confused me for a while!) and then when we install this I into p 157 D, we get final result 6.117. Happily, this agrees with the solution he quoted and which I derived in Exercise 6.13, page 108 and see notes above. I note that you could call the factor multiplying f(ξ) a "Poisson kernel" type of object. That is, if we put f(ξ) = δ(ξ) as a spike on the top rectangle edge, this I think would be the Laplace solution. If we were to put a point charge at location (x,y) inside the rectangle and if we were to put a metal perimeter with u = 0 on the entire perimeter σ, then -I would be the charge distribution induced on that perimeter! Of course you have to evaluate this I thing on the perimeter. All four edges would have charge. Would make another interesting Maple plot, or an exercise. I could repeat my plots above that I made for the disk. Mid page 158: he comments that the limit shown should work out right, and in fact I did do something with that in my Exercise 6.13 solution see notes abobve. Bottom page 158: A Practical Application. We specialize to a square rectangle b = a and to the case where we apply f = 1 across the top (and f = 0 on the other three edges), and this gives u as in p 158 D. You can think of this as a steady state heat problem where you have temperature reservoirs to hold the edge temperature to 1 on the top and 0 on the other three sides. Question: what is the solution temperature in the middle of the plate? We get an excellent answer from just the very first term in the sum p 158 D. He shows the first two terms in p 158 E which again I agree with. This term is .253, and yes we got the exact answer of 1/4 in an earlier exercise. So he is showing how sometimes things converge fast. He has an amusing comment that applied mathematicians use the single-term series approximation as their "favorite weapon". I might add that theoreticians do this a lot too, such as the trace approximation, see elsewhere. Earlier solution gives same result. Now, we solved this rectangle problem earlier in Exercise 6.13 above and got this result: u(x1, x2) = Σn=1∞ bn [ sinh(nπx2/a)/ sinh(nπb/a)] sin(nπx1/a) bn = (2/a) !Syntax Error, Idx1 f(x1) sin(nπx1/a) If we set f = 1, then we find that bn = (2/a) !Syntax Error, Idx1 sin(nπx1/a) = 4/(nπ) n = odd = 0 if n = even so our solution was u(x1, x2) = Σn=1,3,5.. 4/(nπ) [ sinh(nπx2/a)/ sinh(nπb/a)] sin(nπx1/a) If we go to the middle of the rectangle, we have u(a/2, b/2) = Σn=1,3,5.. 4/(nπ) [ sinh(nπb/2a)/ sinh(nπb/a)] sin(nπ/2) and if further we make the rectangle be square, we get u(a/2, a/2) = Σn=1,3,5.. 4/(nπ) [ sinh(nπ/2)/ sinh(nπ)] (-1)(n-1)/2 = (4/π) Σn=1,3,5.. [(-1)(n-1)/2/n ] [ sinh(nπ/2)/ sinh(nπ)] and this replicates the series we have on page 158 D or E. So we are not getting something new here really. Interesting sum identity. It must then be true that Σn=1,3,5.. [(-1)(n-1)/2/n ] [ sinh(nπ/2)/ sinh(nπ)] = π/16 but I don't know how to show this. It is true that [ ] = (1/2) [ 1/cosh(nπ/2) ] = sech(nπ/2)/2 so we then would need to show that Σn=1,3,5.. [(-1)(n-1)/2/n ] sech(nπ/2) = π/8 At this point, we could use n = 2m+1 to get Σm=0,1,2 and also (-1)(n-1)/2 = (-1)m . Then we have to show Σm=0,1,2 ((-1)m /[m+1/2]) sech( π (m+1/2)) = π/4 I know this sum is correct because Maple tells me it is so: evalf(sum (((-1)^m/(m+1/2))* sech(Pi*(m+1/2)), m=0..20)); but I don't know how to prove the sum is correct analytically. I fiddled a bit with this, but I can't go off and derive everything in the world, so I now let this go after several hours. The proof will surely be some contour thing. Example 2. Green's function for the unit circle. (159) We just turn the same cranks and push the same buttons as in Example 1. Our full 2D Green's equation is first stated in 6.118 [ like 6.111 for rectangle]. If we start with the azimuthal EF's of φ, we can expand as in p 159 B. [ like 6.112 for rectangle]. We process the PDE and end up with a Green's equation D just for the expansion coefficient gn(r; ro) [ in the rectangle we had r' = ξ,η. Here we would have r' = ro, φo but we have just set φ0 = 0 to simplify things. the PDE was p 156 A for the rectangle. ] In this problem, we always have to treat the n=0 term separately due to its A + B ln(r) nature, no problem. We solve for the radial Green's in the usual 1D manner: find solutions on each side that meet BC's and that match at the middle, which here means r = r0. Assume unknown constants, and we have p 159 E and F. By applying the usual g' jump condition, the constants become known as in p 160 B and C. Notice the 1/n factor. Next, plug these solutions for gn(r; ro) back into expansion p 159 B and the result is p 160 E -- note of course the azimuthal expo in the second term sum. The two n sides can be folded to the positive side to get an explicitly real version which is p 160 F. We then add back the φ0 in the usual manner, and we get our final answer for "the Green's function for the unit circle" in 6.119 where r< = min(r,r0) etc [ compare this to rectangle results p 156 E and 6.113: we have a single sum, we have the original EF's that we started with, and then multiplied by a Green's for the other variable.] Now if we want to use this to solve Dirichlet on the unit circle, we need that "Poisson kernel" thing (called I) which is -∂ng, so he uses our new g to compute -∂ng and the result is p 160 J. He does not do the next step, but you can replace the sum in there using 6.21 p 103. You end up with our famous result 6.11 which contains the original "Poisson kernel". Of course 6.11 was derived by just doing an expansion on einφ and doing the trivial sum shown bottom page 94, but here we get the same result very indirectly by first computing g. We could alternately start with a base set of EF's for the radial EV equation which I know are some kind of Bessel J functions. Then we could repeat this entire analysis with the variables reversed. In this problem, the two variables r and φ are NOT on an equal footing, so this does require some work, and Stakgold saves that effort for us to do in Exercise 6.39 to come on page 167. I did not verify all the equations in this section in the interest of time. There is nothing "new" here. It is of course an excellent example, and the tie back to known results is very good. Example 3. Green's function for two horizontal plates with a line charge (2D strip). (161) The idea is that we put a plane at y = 0, and another plane at y = a. These are perp to the viewer who is looking down the +z axis into paper. Between these two horizontal planes we have a z-parallel line charge which passes through the xy plane at x=0 and y=η. Since the problem has no z dependence, as usual we treat it as a 2D problem. We look only at the z = 0 plane, where the line charge appears as a point charge and the two plates appear as a horizontal strip. We want the Green's function g for this situation (g = 0 on both those horizontal boundaries), and then we want to use that to solve Dirichlet. It would seem that there would be an "image method" solution to this problem. You could also think of this as the rectangle problem re-centered at the origin, then take limits as a/2 → +∞ and -a/2 → -∞ (ie, a→∞). So having set this problem up, let's turn the crank. The new wrinkle here is that we have an infinite interval in one of our dimensions. In this dimension, we get into Sturm-Liouville analysis due to the infinite endpoints, perhaps some Weyl's theorem and pseudo eigenfunctions and continuous eigenvalue spectra. So Stak's approach is going to be to take this infinite direction as the basic EV set. The EF's are eiαx and the EV spectrum is the entire real axis (the discrete EV's converged!) What this means is that our "opening" expansion and recovery formula [ analogous to p 155 C for the rectangle] involve not a sum but an integral, and it is just the Fourier Transform. What might have been called gα(y) is instead called gα(y) = (1/2π) g^(α,y) which then allows Stak to use our former FT notation. OK, now we have handled the x direction, and we come up with a Green's equation in the y variable for gα(y) which appears in D, which we get by repeating our usual processing steps on the original PDE. Notice that η appears in the δ function. This time instead of showing the various steps involved in getting the g^(α,y) Green's function, he just quotes the result ion p 162 A where y< = min(y, η). The final Green's g is then an integral of this thing times the expo e-iαx (instead of a sum). He then in C writes this out more formally, showing all the right variable names, as in g(x,y | ξ,η). As in the previous examples, we now want to compute I so we can do the Dirichlet problem. In previous examples, I was a sum, but here it is an integral as in 6.125. This thing is evaluated on the top edge of our strip where we are going to impose some f(x). Then our Dirichlet answer is 6.126 with our "Poisson kernel" I. At bottom of p 162 Stak starts over with the EF's in the y direction which are the usual sines. So the starting expansion is then p 162 F. The usual processing then gives us a Green's equation for gn(x) as shown in p 163 A. Ie, the starting EF's were in y, so our Green's is now in x. As before, he just quotes the solution of this Green's equation in p 163 B. I has expo decay behavior in both x directions. Then we jam this into our starting expansion and we get g as in p 163C. This result is then "written out" in proper coordinates as in 6.127. We are invited to compare this to the g we just got doing things in the other order. He comments on a way to show that the two results are equal! He then goes on to compute I doing things this second way and gets 6.128. He claims that this form is clumsy to work with and not very convergent, so he messages this bad I into the form D with E, then F, and finally we get to p 164 6.129. This is, he says, the "best form" of the Poisson kernel for this problem (the strip) and of course you just stuff 6.129 (which has no sums or integrals) into 6.126 and you have your Dirichlet problem solution for the strip with some f(x) along the top edge and u = 0 on the lower edge. Comments: These three examples were rather tedious, but arm the student with a place to look should he or she need to do such a calculation some day. Again, these were examples of using the "partial EF expansion method" to compute a Green's function g. We did this for R = rectangle, disk, and strip. In each case, once we got the g function, we computed I to get a formula for the Dirichlet solution. (5) Complex Variable Method for 2D (164) There are several components to the theory here. First we have the Riemann Mapping Theorem which is mentioned at the start of Ahlfors Chapter 6 (which I have not read yet, having read only Ch 1 thru 3), which says this: in the complex plane, given any simply-connected region R, and given an arbitrary point z0 lying in R, there exists an analytic mapping f(z) (1-1) which will map R into the unit disk, and will map z0 into the origin. It does not say how to find this mapping f. However, we are invited to think of the mapping as being in two parts that you concatenate. Assume your original point is z0 inside R. You first find a mapping f1 that takes R into the unit disk and takes zo in R to some location z0' = f1(z0) inside the unit disk. Once you have done that, there is an obvious second mapping f2 which maps the unit disk into itself and brings z0' to the origin, so that f2[z0'] = 0. That mapping f2 is claimed to be this: f2(z) = eiα (z-z0') /(1-z0') 0' ≡ 1/z0' Recall Ahlfors's discussion of pairs of symmtry points where he had z and z* = R2/ where R is the radius of the circle in question which here is R=1. So in our case, the point symmetric to z0' is 1/ 0' which is of course along the same ray as z0'. Stak refers to z* as an "inverse point". First let's show that the unit circle goes into itself. We can compute (set z0' = a here) |f2|2 = [ z - a - z + a ] / [ 1 - z - a + za ] where note that Stak has made two offsetting sign errors! Now, if z is on the unit circle, z = 1 and |f2|2 = 1, so we map unit circle into unit circle. Now we just have to make sure the interior goes to the interior. Since z0' is in the interior, and since f2(z0') = 0, this is verified. So this thing "shuffles" the inside of the unit disk in some manner. Then we can take f(z) = f2[ f1(z)] as our full Riemann mapping. So f(z) = f2[ f1(z)] = eiα (f1(z)- f1(z0)) /(1- f1(z) 1(z0)) = eiα (s(z)- s(z0)) /(1- s(z) (z0)) where Stak uses f1 = s, and thus we arrive at p 166 B. Now this is all very nice, but what is the connection to Green's Functions? Well let's imagine we have found our analytic mapping function f(z) which we will call f(z,z0) since it also depends on z0. Then here is The Big Claim: g(r|r0) = - (1/2π) ln | f(z,z0) | r = (x,y) z = x+iy etc (a) Now, we know that F ≡ ln(f) = ln(|f|eiφ) = ln|f| + iφ, so ln|f| is the real part of F = ln(f) which is analytic since f is analytic (I guess we stay away from f = 0 which means z = z0 which is the usual Green's singular point). But Ahlfors showed us that both the Re and Im parts of an analytic function like F are harmonic, ie, they satisfy Laplace. Thus, our g(r|r0) satisfies Laplace. (b) Next, suppose z lies on the boundary σ of R. Then f(z) lies on the unit circle so |f| = 1 which then means that ln|f| = ln(1) = 0. Thus, g = 0 on the boundary σ. Things are looking good! (c) Now our Green's equation says -2 g(r|r') = δ2(r-r') where I now think of r = (x,y). It is probably a lot of work to apply 2 to ln|f|, so we instead use the distribution meaning of δ and we want to show that ∫d2r 2 g(r|r0) = -1 for "volume" = ε disk around r' as ε→ 0 But ∫d2r 2 g(r| r0) = ∫d2r g(r| r0) = ∫dS g(r| r0) = ∫dS ∂r g(r| r0) = ε ∫dφ ∂r g(r| r0) where dS = dl = Rdφ = εdφ for our little disk. Now write (z-z0) = εeiφ for point z lying on our little integration circle. Then write [ want to distinguish r = (x,y) from z = x + iy ] g(r|r0) = - (1/2π) ln | f(z,z0) | ≡ G(z,z0) ∂r g(r|r0) = ∂εG(z,z0) = - (1/2π) ∂ε ln | f(z,z0) | = - (1/2π) 1/|f| * ∂ε |f(z,z0)| (0) But now write this Taylor expansion for z on the circle near z0 , and use f(z0,z0) = 0, f(z,z0) = f(z0 + εeiφ,z0) = f(z0,z0) + εeiφ f '(z0,z0) + ... = εeiφ f '(z0,z0) => (1) ∂ε f(z,z0) = f '(z0,z0) eiφ and (2) |f(z,z0)| = ε | f '(z0,z0) | and (1) => (3) ∂ε |f(z,z0)| = | ∂ε f(z,z0) | = | f '(z0,z0)| Therefore, inserting (3) into (0) and then using (3), we get ∂r g(r| r0) = ∂εG(z,z0) = - (1/2π)(1/|f(z0,z0)|)* | f '(z0,z0)| = - (1/2π) |f '| / |f| = -1/(2πε) so our little integral is then ε ∫dφ ∂r g(r| r0) ≈ ε∫dφ { -1/(2πε)} = ε { -1/(2πε)} ∫dφ = ε { -1/(2πε)} 2π = -1 QED This last part is not really obvious until you write all this out, some of which Stak did. Example 1: What is the Green's Function for the unit circle itself? We want g(r|r0) where r0 is some point in the unit circle. Here we have in effect f1 = 1 ( meaning z0' = z0) and we then have f = f2 = s, so g(r|r0) = - (1/2π) ln | f(z,z0) | = - (1/2π) ln | f2(z,z0) | = - (1/2π) ln |(z-z0) /(1-z0)| = - (1/2π) ln |z-z0| /|1-z0| = - (1/2π) ln |z-z0| /(|z-1/ 0||z0|) = - (1/2π) ln |z-z0| /(|z-z0*||z0|) = - (1/2π) ln |z-z0| + (1/2π) ln |z-z0*| + (1/2π) ln |z0| and this agrees with our earlier result p 149 H where |z0| = r0 and z0 = ξ if you like and z0* = ξ*. There Stak did use * to indicate the inverse point, but he did not call it by that name there. Then in the current section, he does NOT use * for the inverse or symmetry point! It seems pretty clear that, for a circle of radius a, you replace z → z/a and z0→ z0/a and this gives us page 166 C. History: George Green invented Green's functions around 1825. Self-taught and isolated, his work was not noticed until 1833. So we must have had the concept of a Green's Function around that time. http://greenfunction.blogspot.com/2007/06/george-green-miller-mathematician-and.html "To quote Sir Edmund Whittaker’s authoritative history of Theories of Aether and Electricity – "It is no exaggeration to describe Green as the real founder of that "Cambridge school" of natural philosophers of whom Kelvin, Stokes, Lord Rayleigh and Clerk Maxwell were the most illustrious members in the latter half of the nineteenth century" " "Green’s mathematics was nearly all devised to solve very general physical problems. His first interest was in electrostatics. The inverse-square law had recently been established experimentally, and he wanted to calculate how this determined the distribution of charge on the surfaces of conductors. He made great use of the electrical potential, and gave it that name, and one of the theorems that he proved in this work became famous and is known as Green’s theorem." I downloaded a copy of his famous paper, but it is long and hard to find the Green's Function in it. Probably there was no δ notation at the time, etc etc. So OK, Green's function's circal 1828. Now we come to the mapping part of all this. Concerning the Riemann Mapping Theorem: "The theorem was stated (under the assumption that the boundary of U is piecewise smooth) by Bernhard Riemann in 1851 in his PhD thesis. Lars Ahlfors wrote once, concerning the original formulation of the theorem, that it was “ultimately formulated in terms which would defy any attempt of proof, even with modern methods”. Riemann's proof depended on the Dirichlet principle (whose name was created by Riemann himself), which was considered sound at the time" "The first proof of the theorem is due to Constantin Carathéodory, who published it in 1912. His proof used Riemann surfaces and it was simplified by Paul Koebe two years later in a way which did not require them." So I guess we can say the connection was made around 1851. Enough on this. Exercises (166) Since this doc is now 126 pages long, I am going to continue it in a separate document because things are getting a little sluggish and the chances of a "problem" I know are larger with large documents. I may decide later to combine these two documents.