stakgold chap 7 raw 1

Stakgold Chapter 7 notes PhL 3.22.11 7.1 Introduction 2 The basic specification of a heat or wave equation BV problem. 2 Comment clarifying the BC's: 3 Analogy with electrostatics. 3 Green's Theorem for heat and wave operators. 4 What if region R is unbounded? 4 7.2 The Causal Green's Function for Heat Conduction 4 Specification of the Green's Function problem and its alternate form. 5 So where exactly do we want g = 0 in our heat conduction version of potential theory? 7 The free space solution for g, connection to Schrodinger Equation. 7 Facts about the adjoints L* and g*. 8 Review of Green's Function role in potential theory. 8 Starting the derivation of 7.14 8 Revisit the free space solution for g. 9 Continuing the Derivation of p 199A and (7.14) 9 General Comments on What has Happened Here. 12 Application to problems with an Infinite Rod (200) 14 Meaning of "a rod". 14 No surface integral for n = 1 in 7.14. 14 In (7.15) we get the full-rod problem over all of R 14 Stak's comments about u depending continuously on the initial data. 15 Infinite rod with dipole impulse at the center. 16 Teapot tempest in the dipole problem. 16 [p 202 rewrite the general infinite rod solution] 16 Exercise 7.1: Get (7.14) for The Neumann case ∂nu = h on σ. 17 Exercise 7.2: Get (7.14) for The Radiative case ∂nu = h on σ. 17 Exercise 7.3. Interpreting the BC's as additional source terms. 19 7.3 Methods for finding Causal Green's Functions (all 1D space) 20 A. The Method of Images and other Trick Methods. 20 [ p 204 The Half Rod Green's Problem, g=0 at end (image method)] 21 [ p 205 Half Rod with u=0 on the end (ie, h(0,t)=0) , but u(x,0) = f(x). ] 22 [ p 206 The Above solution for large t. ] 25 [ p 206 Time Laplace method for half-rod with f = 1 ] 25 [ p 206 The f=1 half rod solution for small t ] 25 [ p 207 The f=0 half rod with left end driver h(t). ] 25 [ p 207 taking the t = 0+ limit of the above half rod solution.] 26 [ p 208 a Subtraction Method for the half rod ] 26 [ p 209 the half rod with insulated end, Neumann ∂ 27 [ p 209 Green's for the half rod with radiative BC ∂ 27 [ p 211 Green's for an insulated wire ring ] 28 [ p 212 Green's Function for a finite rod of length l with g=0 at each end] 28 B. The Method of Expanding in Spatial Eigenfunctions: h=0, initial f(x) (213) 29 [ 214 Causal Green's for the above.] 30 [ 214 Repeat the above for Neumann case ] 30 [ 215 R unbounded.] 30 [ 215 The general u(x,t) problem by a non-EF method] 30 C. Examples of the Eigenfunction Expansion Method (216) 33 Example 1: Green's Function for finite rod by EV method u = 0 at both ends. 33 Example 2: The half-rod with radiation from left end (a continuous λ problem) 33 Radiation explained: 33 What is the relevance of Lλ = δ(x-ξ) p 217B? 33 First, what do we know about this in Potential Theory? 33 Apply this to Volume II page 217 for his "second method". 34 So what was the purpose of the Lλ equation shown as p 217 B? 34 What is the "first method" Stak refers to on page 218 top? 35 D. The Laplace Transform Method (218) 37 7.4 Uniqueness and Continuous Dependence on the Data 39 [224] The Maximum Principle (Theorem 2) 40 [225] The Minimum Principle 40 [226] Theorem 4 is the thing about continuous dependence 40 7.5 Miscellaneous Heat Conduction Equation Topics 41 The Semigroup Connection. 41 [228] The Backward Problem. 42 The Diffusion Interpretation of the Heat Conduction Equation 45 Asymptotic Formula for The EV's of -2 45 A 3D Composite Medium Heat Conduction Problem 49 The Stefan Problem (237-238) 50 Exercise 7.4 Solve the ring problem by doing Fourier in the x variable. 51 Exercise 7.5 Redo the Stefan problem with ice at u = V < 0. 53 Exercise 7.6 Freezing of water in a cylindrical tank. 57 Exercise 7.7 Obtain the Weyl formula for an n=3 region. 59 Exercise 7.8 Obtain the Weyl formula for Neumann and Radiative BC. 61 Exercise 7.9 Another way to show max and min of harmonic u both lie on the boundary. 61 Exercise 7.10 Finite rod with dipole source in center. 62 Exercise 7.11 Average Temperature is a constant in time. 64 Exercise 7.12 Add a cu linear term into the heat equation and resolve in EF method. 66 Exercise 7.13 Maximum principle for fancy diffusion equation. 67 Exercise 7.14. Show Theorem 1 page 223 for Neumann and for Radiative BC's. 68 Exercise 7.15 Show Theorem 1 page 223 for the fancy diffusion equation 7.108 if c ≥ 0. 68 Exercise 7.16. The half-rod driven from the left end by various methods. 68 Exercise 7.17. Deriving the Weber Transform 68 Exercise 7.18. Practice with the Weber Transform 69 7.1 Introduction The basic specification of a heat or wave equation BV problem. We are going to work in n+1 dimensions (t,x) and we are going to consider two equations. One is the parabolic heat-conduction equation, and the other is the hyperbolic wave equation. In each case, we are going to think of a physical region of real space R (we allow this to have n dimensions, but usually we will be concerned with n = 3 or 2 or 1 ). We will look inside this space at some function u(t,x) defined over all of R, and at each point x in R this function "evolves" with time t. This evolution is represented by the n+1 dimensional cylinder shown on page 195 whose sides are parallel to the t axis in n+1 dimensional hyperspace. It is really a hypercylinder. I will use n = 3 in the following discussion. The two boundary value problems are clearly stated, 7.3 and 7.5. For the heat case, things are driven by a time-varying 3D heat source function q(t,x) out in the region R. We ignore t < 0 completely. The function u(t,x) is the function which solves the heat conduction equation, called the response function. We have two boundary conditions for the heat situation. One is that u(0,x) = some prescribed function f(x) at time t=0. The other is that the full function u(t,x) is forced to a prescribed value h(x,t) on the cylinder walls, x on σ. In potential theory, this would be the Dirichlet problem and u would be the potential. It would be possible instead to specify ∂nu(t,x) on the walls of cylinder and this would be the Neumann problem. Comment clarifying the BC's: It is good to reinforce this BC situation: if R were a disk, then the initial value condition f(x) specifies the solution over the entire disk at t = 0 ( which is the entire bottom surface of the cylinder), while the boundary value condition h(x,t) specifies the solution just at the perimeter of the disk, but at all times (which locus is the entire walls of the cylinder). One might wonder if you also need to specify the solution over the entire disk at time t = T ( which is the entire top surface of the cylinder), with the idea that you require specification on a closed surface! When we soon look at Green's theorem on our closed surface (for the heat equation L) we will find it to say !Syntax Error, Idt !Syntax Error, Idnx [ v(t,x) (∂t - 2)u(t,x) - u(t,x) (- ∂t - 2)v(t,x) ] = ∫ dnx [u(T,x)v(T,x)- u(0,x)v(0,x)] + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nv(t,x) - v(t,x)∂nu(t,x) ] (7.6) when applied to two arbitrary functions u and v. You see the full volume integral on the LHS over the volume of the cylinder. On the right we have all the surface terms. The last term is an integral over the sides of the cylinder. The first term has two pieces, and yes, the first piece with u(T,x)v(T,x) is one over the top end of the cylinder, and u(0,x)v(0,x) is over the bottom end. So in this general statement, yes, we need to deal with surface terms on ALL surfaces of the cylinder. Now, our next development step is to say that u and v are not just arbitrary functions, but that u is our desired solution function to 7.3 and v is the Green's Function v(x) = g(x0|x) of 7.8 with arguments swapped. Now since v(x,T) = g(x0|x,T) = g(←), we see that our causality boundary condition on g ( which we have not mentioned yet!) causes v(x,T) = 0 since we assume that the point x0 with its t0 lies in the cylinder, so g = 0 since it is propagating backwards in time from any spacetime point on the top surface. This then makes the top surface integral disappear, and then we have only surface integrals over the bottom and sides of the cylinder!!! Analogy with electrostatics. In the heat world, u is temperature, and ∂nu(t,x) on the boundary is heat (thermal energy) flow into our out of the boundary (whereas q(t,x) is a volume source in the interior). The analogy between electrostatics and static heat flow is presented by me in "potential analogies.doc" located in the Stak folder. The electrostatic analogy of q would be the 3D charge density ρ (sourcing field lines) and ∂nu is like ∂nV being surface charge density σ, also sourcing field lines. But here we are talking the time-varying heat conduction equation, and only its static limit can be identified with potential theory. Time varying electrostatics gets us into Maxwell's Equations which is how that static world moves in time, and those equations are certainly not the same as the heat conduction equation. ( they are more wave like than heat like). As for the wave equation, the boundary value problem as presented in (7.5) is exactly the same, except since second order in time we have to specify u(0,x) = f1(x) and ∂tu(0,x) = f2(x). The BC on the sides is the same as for the heat case. Doubtless there are Dirichlet, Neumann and Mixed possibilities for the BC on the cylinder side here as in the heat case. Stak says he will concentrate on practical stuff, but will talk a bit about existence of solutions etc. Green's Theorem for heat and wave operators. He then writes down Green's Theorem for the two equations. This is really just the divergence theorem in n+1 dimensions applied to the vector Jμ = (J0,J) which is the (non-unique) surface current associated with a differential operator L. In our electrostatics analogy, we call this thing a current in the sense that if you integrate it over a boundary to get the total outflow of charge, you get the time change in the total charge enclosed, based on the continuity equation ∂μJμ = 0 or J = -∂tρ. But in our application for operator L, J is just the "vector thing" that you get in your Green's Theorem on the surface integral side of things. I proved this general Green's Theorem (page 40, (5.72) and (5,73)) in great detail here: "parts and greens.doc" // result stated in Section 5 " attempts to prove the full p 40 Green Theorem.doc" The basic idea is that the current surface term is the "parts" that results after you do lots of parts integrations to swing operator L from one side to the other side so the LHS is the integral of u(Lv) - (Lu)v and the RHS is then the parts terms. There are lots of parts terms because L has multi-index k and n variables etc etc. Then I derive the specific case for n+1 dimension heat and wave L (7.6) and (7.7) here: "N dim - Geometry, Div Thm, Green.doc" and the results are stated in the next section where we need to use them. The second doc is in a bad state, What if region R is unbounded? Stak's final comment on page 197 concerns what happens if your spatial region R = Rn happens to be infinite. This is represented by the elliptical boundary in the page 197 drawing. His plan is to sort of regulate this thing by also drawing a math sphere of radius r about some point, and treat R as the intersection of these surfaces. In my notation your σn in Green's Theorem consists of the two pieces he shows, one piece being on R, the other being on your math sphere portion inside R. Of course the plan is then to go to the limit r→∞ and assume that your response function u decays appropriately far away so things can be controlled in the limit. 7.2 The Causal Green's Function for Heat Conduction Comment: I have redone this section more generally in "BV problems for Laplace and Heat.doc", where we have Dirichlet, Neumann and Radiative cases all treated uniformly, and where I constantly compare the heat situation to the Laplace situation of potential theory. Things are clearer to me now than when I first wrote the notes below. Specification of the Green's Function problem and its alternate form. Back in Chapter 5 on page 59 it was shown that you could replace the Green's function problem 5.135 with the homogeneous alternate problem 5.136 which is a simpler problem, where C = Hu. I am completely happy with this, see Chap 5 raw notes. [ but I show it a foot below anyway! ] Now 7.8 is like 5.135, and 7.8a is like the simpler problem 5.136. Instead of using u and C, we use g everywhere, but note that in 7.8a we only have t > 0. The time origin is shifted here to t0 instead of 0, no big deal. Notice in the statement of our Green's Function Problem 7.8 that we want g = 0 for t<t0 but we also want g = 0 on the sides of our famous cylinder. What happens at t = t0 has to be determined by solving the δδ driven PDE. and we want g = δ(x-x0) inside the bottom (t=0) where the δ source is located at some point x0 (previously called ξ)! I am going to use my own notation as in the last doc above so that bold x = a point in R, possibly on the boundary σ of R x0 = (t0, x0) = location of the driving spacetime pulse δ ( ie, at t=t0 and Green's space point at x0 ) g(x|x0) = g(t,x|t0,x0) So here is (7.8) where in general t,t0 can be anywhere in (-∞,∞) : (∂t- 2) g(x|x0) = δ(t-t0)δ(x-x0) // driving impulse is at t = t0 and x = x0 g(x|x0) = 0 t<t0 // no response at negative time since causal g(x|x') = 0 x on σ, x' in R // vanish on cylinder walls (7.8) Since I use this elsewhere, here is a suggestive picture to set the stage for later. In this picture I have no integration arrows at the "1" location, so it basically say g = 1*g, and it shows that g = 0 on the boundary (assuming the Dirichlet case). Later we will restrict to all times positive, but 7.8 is valid for all time, so just imagine there is no t=0 marker in the above Feynman diagram. The main idea is that we are going to soon think of both endpoints of the arrow as being at positive time, so they are both in our "cylinder" which runs t = 0 to t = T. The propagator g(x|x0) has the sense of g(←) only (non-zero only in this direction), by which I mean t > t0. If you try to run it backwards, you get 0. This is the causality idea being built in at the start. So this heat evolution world, g is NOT going to be symmetric, so we must be more careful: g = 0 when the first argument has spatial x on the boundary σ of R. Now here is the alternate problem (7.8a), where now g is only defined for t > t0 ( remember that THIS g is like the u in Chapter 5, previous g is like C). Stak seems to add an extra condition that both t,t0 must be positive, not sure why he does this. So here is (7.8a). Although g and h are identical for t > t0, I will write 7.8a in terms of the function h : g(x|x0) = θ(t>to) h(x|x0) g = h for t>t0 (∂t- 2) h(x|x0) = 0 // homogeneous! h(x|x') = 0 x on σ, x' in R limt→t0+ h(x|x0) = δ(x-x0) // spatial point source on the cylinder bottom. (7.8a) In the following section, I will show how 7.8a arises from 7.8, it is pretty simple. We don't say here whether h exists for t < t0 because it won't matter, since in the end it is g that we care about. The Connection between 7.8 and 7.8a, general form for g, compare to potential theory. (1) First, we define g as the solution of this equation when t > t0 and we say g = 0 otherwise (∂t - 2)g(x,t ; x0, t0) = δ(x-x0)δ(t-t0) with BC=0 (whichever it is, Dirichlet etc) Thus, we might write g this way, to make the forward propagation idea explicit, g(x,t ; x0, t0) = θ(t-t0) h(x,t ; x0, t0) At this point we have no idea what the functions g (or h) are. Depends on shape of σ. Now consider: (∂t - 2)g(x,t ; x0, t0) = (∂t - 2)[ θ(t-t0) h(x,t ; x0, t0)] = θ(t-t0) (∂t - 2) h(x,t ; x0, t0) + δ(t-t0) h(x,t ; x0, t0) (*) Suppose we could find a function h(x,t ; x0, t0) such that (∂t - 2) h(x,t ; x0, t0) = 0 h(x,t0+ ; x0, t0) = δ(x-x0) with BC=0 (**) Then if we insert that h into (*) above we find that the red term vanishes and so (∂t - 2)g(x,t ; x0, t0) = δ(t-t0) h(x,t ; x0, t0) Now the usual delta rule says that δ(t-t0)F(t) = δ(t-t0)F(t0) so we then get (∂t - 2)g(x,t ; x0, t0) = δ(t-t0) h(x,t0+ ; x0, t0) = δ(t-t0) δ(x-x0) This then gives us a candidate for g(x,t ; x0, t0) that meets all the conditions of our original problem! Those conditions were: (1) the δδ-driven PDE; (2) g vanish for t<t0; (3) g respects the BC. So perhaps 7.8a provides an easier way to find g than 7.8 since you have a homo PDE. Later on page 225 Theorem 3 Stak shows that the general u heat problem has a unique solution, and therefore so does the g heat problem. If we have found a solution by the h method, that must be it. (2) Now let's jump ahead just for a moment where we will learn that we can write h in this manner h(x,t ; x0, t0) = Σiφi(x)φi(x0) e-λt where -2φi = λiφi with φi = 0 on σ It is easy to show that this h solves the h problem stated above. Here are all three parts: (1) h(x on σ,t0 ; x0, t0) = 0 because φi(x on σ) = 0 so BC = 0 (2) h(x,t0 ; x0, t0) = Σiφi(x)φi(x0) = δ(x-x0) (3) (∂t - 2) h(x,t ; x0, t0) = Σiφi(x)φi(x0) (-λi)e-λt + Σiλiφi(x)φi(x0) e-λt = 0 Of course we know without doing any more work that if g(x,t ; x0, t0) = θ(t-t0) h(x,t ; x0, t0) and h satisfies the conditions shown above (**), then (∂t - 2)g(x,t ; x0, t0) = δ(x-x0)δ(t-t0). (3) Now in potential theory, we have a slightly different situation: -2g(x; x0) = δ(x-x0) BC = 0 g(x; x0) = Σiφi(x)φi(x0)/λi -2g(x; x0) = -2 [Σiφi(x)φi(x0)/λi] = Σiλiφi(x)φi(x0)/λi = Σiφi(x)φi(x0) = δ(x-x0) So where exactly do we want g = 0 in our heat conduction version of potential theory? If we look at the form 7.8, here is what we see. First, we have a Green's point source somewhere inside the cylinder at location x0, t0 which is a point at some height t = t0 in the cylinder. Right off the bat we know that g(x|x0) = 0 at every point inside and on the surface of the cylinder below t = t0, based on the g(←) idea. Thus, we know that g = 0 on the bottom surface. Second, we have g = 0 on the boundary σ of R at all times, which means g = 0 on the walls of the cylinder. So the bottom part of the cylinder walls are set to 0 twice, so to speak, by our conditions. Therefore, we can think of the requirement being that g = 0 on the bottom and sides of the cylinder, whereas in potential theory we had g = 0 just on the perimeter of our disk. It is really the walls of the cylinder that is the important part, and generalizes potential theory. Now, you can ask what g might look like on the slice cylinder right at t = t0 (a slice in which lies in the middle of our cylinder running t = 0 to t = T.). This is where 7.8a comes in. If you take the limit of g from above, you find that g → δ(x-x0) on this slice t = t0. But if you were to take the limit from the bottom, I think you could get g → 0 (if the limit exists). The free space solution for g, connection to Schrodinger Equation. Remember that 7.8 and 7.8a are different problems, but they have the same solution in region t≥t0! If the region R is all of space, we know the solution g = C because we computed it back on page 60 in Chapter 5. Note that this is very much like Saxon's free particle propagator in quantum mechanics (page 61 Saxon) and the reason is that the Schrodinger equation is the same as the heat equation but I suppose with t → it or something like that. In Stak, this free space full region thing is called the Causal Green's Function for Heat Conduction. So again, g = 0 on the boundary σ at all times, which is just like potential theory, but now we have two extra conditions as well. The first says g(x|x0)=0 for t < t0, and the second says g = 0 on the cylinder slice t = t0 except at point x0 where we have δ. We have no condition on the top of the cylinder where evolution is ongoing. Really 7.8 is more like potential theory than 7.8a since in 7.8a we have a homo PDE and the δ condition on the bottom. Facts about the adjoints L* and g*. Now, by goofing around on the top of page 199 with a backwards-in-time version of g called g* which is the solution to the "adjoint problem", he shows that g*(xμ, yμ) = g(yμ, xμ). The reason he needs to worry about g* is that the operator (-∂t-2) appears in our huge ugly Green's theorem gizmo 7.6, and we need to know what this does to g. Of course if L = (∂t-2) then L* = (-∂t-2), so the latter really is the formal adjoint operator. (7.11) is obtained by negating all times in (7.8): t → -t and t0→ - t0. The first line's time δ of course stays the same then. So t < t0 becomes -t < -t0 or t > t0 . In the last line, we still have the first argument be on σ. The proof p 199A that g*(x|y) = g(y|x) requires use of 7.6 and I am going to skip this proof for now, but this is a key result, I don't think the proof is too hard. This result says G*xy = Gyx = GTxy so that operator G* = GT in matrix language, where G is an integral operator generated by g(y|x). In potential theory we always had L = L* and we ended up with G* = G = GT (symmetric). Here in heat theory we do not have L = L*, so not surprising that the G property is a little different. Now comes the big moment of truth with 7.6 becoming 7.14, but let's hold the phone a moment and go back to potential theory. Review of Green's Function role in potential theory. In potential theory for the equation we have a Green's Theorem 5.73 p 40 which I write here ∫R dx ( v 2u – u 2v ) = ∫σ dS (v ∂nu - u∂nv) Now if we assume that 2u = q(x), some Poisson source, and we set v = g (which vanishes on σ), this Green's Theorem becomes our familiar result (for the Dirichlet situation g = 0 on boundary) u(x) = ∫R dξ g(x|ξ) q(ξ) – ∫σ dSξ f(ξ) ∂ξng(x|ξ) which is the instruction for using the Green's Function g to get the solution of any Dirichlet problem with prescribed potential f(ξ). The solo u(x) term comes from the fact that -2g = δ, and the v∂nu surface term vanishes since g = 0 on the surface, and the other surface term is what you see! Starting the derivation of 7.14. So, on page 199 Stak does this same thing for "heat conduction theory" instead of potential theory. He sets u = u (solution of 7.3 with source function q) and he sets v = g, and the 7.6 becomes p 199 A which he rewrites as 7.14. Here you see the locations of the heat source function q, the initial temperature distribution in our space R as prescribed by f(x), and finally the cylinder side function h(x,t) which is sort of the new kid in town: we allow the prescribed "potential" on the boundary to vary in time! If you know the three functions q,f,h and if you know your Causal Green's function g, then this equation tells you u(x,t) everywhere and at all time. As the time integrals in 7.14 show, u(x,t) is only affected by "past history" just as you expect. I am pretty sure Feynman would represent this equation with a Feynman Diagram where g is a free space heat propagator, though I don't know quite how that would work. Let's go back to our definition of g stuff (∂t- 2) g(x|x0) = δ(t-t0)δ(x-x0) // driving pulse is at t = t0 and x = x0 g(x|x0) = 0 t<t0 // no response at negative time since causal g(x|x') = 0 x on σ, x' in R (7.8) Revisit the free space solution for g. If it happens that R = En , then entire space, then we know g(x|x0) exactly! Recall in potential theory how we had - 2g(x|x0) = δ(x-x0) and we said if there are no local boundaries, then g(x|x0) = 1/4π|x-x0| and we called this the "fundamental solution" E(x|x0) = 1/4π|x-x0|. We solved this first with x0 = 0 and in spherical coordinates with no boundaries we wanted to look for a solution with no angular dependence and we then found that E(x|0) = 1/4π|x| = 1/4πr. We associated this with V = q/4πr for a point particle in isolation (in Stak units). Now a similar thing happens for our heat equation if there are no boundaries. We start with t0 = 0 and we show that there is a "fundamental solution" for g that we call C which is this C(x,t|x0,0) = H(t) [4πt]-n/2 exp(-|x-x0|2/4t) Then we do our usual time shift and this becomes g(x|x0) = C(x,t|x0,t0) = H(t-t0) [4π(t-t0)]-n/2 exp(-|x-x0|2/4(t-t0)) (7.10) You might think of this as the "heat potential" (so to speak) of a "unit point heat source", just as we think of V = 1/r as the "electrostatic potential" of a "unit point charge". At some point we ought to find an integral which is like our ∫dV σ/R integral in electrostatics, but it will be ∫dV q C. This is what the first term of 7.14 says in the absence of boundaries since then g = C, BUT we still have the second term in 7.14 even if there are no boundaries. My inclination is to call this C thing a "free space heat propagator". If there are boundaries present, we will be adding to this thing another piece which (1) solves the homo equation (∂t- 2)f = 0, and which (2) causes things to be right on the boundaries. But in any problem where you have just the entire space, this thing C is the propagator. In non-rel QM the heat equation is the Schrodinger equation with t→it and the above thing really is called the "free space propagator" (Saxon p 61, already noted above). In the limit that t → t0, the C object becomes just δ(x-x0), and we see the general idea that as time moves forward, g "spreads out" in space. Note what happens as t-t0 gets large as t goes into the future. The exp part is never larger than 1 and the other factor drops off as 1/ (t-t0)n/2 so basically C(x,t|x0,t0) → 0 as t→∞. This is true of derivatives like ∂xC(x,t|x0,t0) as well. We are happy to see "influence" fading away as we move into the future, the propagator being the carrier of influence. For n = 3 space, it is t-3/2, so not inverse square in time. Continuing the Derivation of p 199A and (7.14) We start with 7.6 which I quote from "N dim - Geometry, Div Thm, Green.doc" [ 4.15.11 I cannot find this doc! !Syntax Error, Idt !Syntax Error, Idnx [ v(t,x) (∂t - 2)u(t,x) - u(t,x) (- ∂t - 2)v(t,x) ] = ∫ dnx [u(T,x)v(T,x)- u(0,x)v(0,x)] + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nv(t,x) - v(t,x)∂nu(t,x) ] (7.6) We set v = g, the solution to problem (7.8a), and we set u = u, the solution to heat problem (7.3). We need to be more specific with v = g: we want v(x) = g(x0|x) so it is the 2nd argument that is in v. Then we have (notice we write ∂n = ∂nx in one term to make clear what we mean! ) !Syntax Error, Idt !Syntax Error, Idnx [g(x0|x) (∂t - 2)u(t,x) - u(t,x) (- ∂t - 2) g(x0|x) ] = ∫ dnx [u(T,x) g(x0|x)- u(0,x) g(x0|x)|t=0] + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) - g(x0|x)∂nu(t,x) ] Now let's make these insertions: (∂t - 2)u(t,x) = q(t,x) (- ∂t - 2) g(x0|x) = (- ∂t - 2) g*(x|x0) = δ(x-x0)δ(t-t0) // from 7.11 line 1 g(x0|x) = 0 on σ With just these three items we get !Syntax Error, Idt !Syntax Error, Idnx [g(x0|x) q(t,x) - u(t,x) δ(x-x0)δ(t-t0) ] = ∫ dnx [u(T,x) g(x0|x)- u(0,x) g(x0|x)|t=0] + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) ] or !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) - u(t0,x0) = ∫ dnx [u(T,x) g(x0|x)- u(0,x) g(x0|x)|t=0] + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) ] or changing all signs - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = - ∫ dnx [u(T,x) g(x0|x)- u(0,x) g(x0|x)|t=0] - !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) ] We are assuming that T > t and t0 and T is some large value. But g(x0|x) = 0 when t0 < t so that g = 0 when t > t0. I think it is safe to assume that ∂n g(x0|x) = 0 as well in this range! Thus we can adjust both the upper time endpoints to get this - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = - ∫ dnx [u(T,x) g(x0|x)|t=T- u(0,x) g(x0|x)|t=0] - !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) ] Now g(x0|x)|t=T = 0 because we just said above that g(x0|x) = 0 when t > t0 and T > t. So that knocks out another term and we have - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = + ∫ dnx [ u(0,x) g(x0|x)|t=0] - !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) ] Now we replace u(0,x) = f(x) and on the boundary we have u(t,x) = h(t,x) - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = + ∫ dnx f(x) g(x0|x)|t=0 - !Syntax Error, Idt !Syntax Error, I dSn [ h(t,x)∂nx g(x0|x) ] or u(t0,x0) = !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + ∫ dnx f(x) g(x0|x)|t=0 – !Syntax Error, Idt !Syntax Error, I dSn [ h(t,x)∂nx g(x0|x) ] // agrees with p 199 A Now I will to the x ↔ x0 swap manually to get u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idnx0 g(x|x0) q(t0,x0) + ∫ dnx0 f(x0) g(x|x0)|t0=0 – !Syntax Error, Idt0 !Syntax Error, I dSn0 [ h(t0,x0)∂nx0 g(x|x0) ] I have been maintaining n as the number of spatial dimensions, but I will not drop that notation and also swap arguments for h and q to match Stak, and the above becomes u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(x0,t0) + ∫ dx0 f(x0) g(x|x0)|t0=0 – !Syntax Error, Idt0 !Syntax Error, I dS0 [ h(x0,t0)∂nx0 g(x|x0) ] // agrees with (7.14) Comment: I have redone all of the above more generally in support doc "BV problems for Laplace and Heat.doc", where we have Dirichlet, Neumann and Radiative cases all treated uniformly. General Comments on What has Happened Here. One might (and should) ask: How does this all differ from just stating Green's Theorem in n+1 dimensions and calling one of those dimensions "time" t ? (a) First, what do we mean by "Green's Theorem" ? The web definition is that "Green's Theorem" usually means just Stokes Theorem for n=2. This meaning is completely unrelated to what we are going to talk about here. Here we shall use Stakgold's definition of "Green's Theorem" (he calls it that on page 40) which we can express in both integral and differential forms: ∫dV (vLu - uL*v) = ∫JdS vLu - uL*v = J Green's Theorem It is the divergence theorem which shows these two forms to be the same. Here J is the (non-unique) surface current which is associated with operator L. As an example of Green's Theorem in 1+n dimensions we can consider L = (∂t - 2), and we find that J = [ uv – ( v u – u v) ] from Stak p 41. Again, this is Green's Theorem in a 1+n dimension case. Stakgold states this in equation (7.6). (b) In "N dim - Geometry, Div Thm, Green.doc" I consider the divergence theorem in m dimensions ∫V dV A = ∫S dSA The Divergence Theorem I then think of m = 1 + n and I name the first coordinate t, and I introduce relativity "like" notation (but with a diag(1) metric tensor!) xμ = (t,x). Breaking off time in this way, the divergence theorem can be applied to Stakgold's "cylinder" as a volume in 1+n dimensions and we obtain !Syntax Error, Idt !Syntax Error, Idnx ∂μAμ(t,x) = ∫ dnx [A0(T,x)- A0(0,x)] + !Syntax Error, Idt ∫ dSn (x) A(t,x) This applies to any vector Aμ = (A0,A) we like. It is nothing more and nothing less than the divergence theorem in 1+n dimensions applied to our cylinder as volume. (c) Going back to L = (∂t - 2) and the Green's Theorem in 1+n dimensions and the fact that in that theorem we have J = [ uv – ( v u – u v) ] , we can write J in our 1+n dimensional notation this way Jμ = (uv, uv - vu) If we then take this to be vector Aμ in our n+1 dimensional divergence theorem, we get !Syntax Error, Idt !Syntax Error, Idnx ∂μJμ(t,x) = ∫ dnx [J0(T,x)- J0(0,x)] + !Syntax Error, Idt ∫ dSn (x) J(t,x) So this equation is the 1+n Divergence Theorem applied to a certain 1+n dimensional vector. If within that special vector Jμ we select v(x) = g(x0|x) (notice that we have reversed the "usual order" of g's labels here), we end up with the famous equation 7.14 in the Dirichlet boundary case (see "BV problems for and Heat.doc" for the other cases). Here is the above equation written out in this case: !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) - u(t0,x0) // volume contribution = ∫ dnx [{u(t,x) g(x0|x)}|t=T – {u(t,x) g(x0|x)}|t=0] // top and bottom contributions + !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) - g(x0|x)∂nu(t,x) ] // contribution of the sides (d) So, in answer to the opening question, yes, equation 7.14 is in fact just a statement of the 1+n dimensional Divergence Theorem applied to the vector Jμ which appears in Green's Theorem in 1+n dimensions, with u and v inside Jμ set to certain functions. If we then ask "In the 1+n dimensional divergence theorem, we have to have a closed boundary in the surface integral, so why is the cylinder top excluded from this boundary in 7.14?" The answer is that we do have a fully closed boundary and it just happens that the top gives no contribution due to the unidirectional nature of the particular function g(x0|x) we have selected for v(x) inside Jμ. That is to say: since t0 is a time within the cylinder, if x is a point on the top surface, then g(x0|x)}|t=T = 0 from the causal fiat definition of the function g. The bottom and sides of the cylinder are then the only boundary contributions to the surface integral you see in 7.14. The bottom is f(x), the sides are h(x,t). (e) We put the causality rabbit in the hat when we defined g to be causal by fiat in 7.8. One result is that this allowed us to have the alternate form 7.8a. Another result of this causality fiat definition is that we can find u(x,t) which solves the system 7.3 by just plugging q,f,h into 7.14. Had we defined g in some other way, we would not have obtained solution 7.14. So let's try for a grandiose Theorem. Theorem: IF we define g with the causal fiat added, THEN we discover that we are able to find a solution (namely 7.14) to the PDE system 7.3. Thus, system 7.3 is a "reasonable class of problems" to study since it has a solution (which we later show is unique). The boundary conditions are neither overspecified nor underspecified, so we have well-posed Cauchy Data. Perhaps there are other classes of problems that can be solved by other methods. Since causality was installed by fiat into g, we find that u(x,t) is only affected by u(x',t') at earlier times. But this aligns with our expectation of what happens in a physical heat problem in the real physical world, and this will also arise in the wave equation. In both cases the concept is evolution (from the past into the future). One could imagine one of those "other classes" of problem where we do not build causality into g, and then we would find that u(x,t) was affected by u(x',t') in both the past and the future. This would be a fine class of problems to study, but this class of problems would not include real world evolution problems where we expect the "arrow of time" to have one direction. This class would include problems where t is just another spatial dimension and we have L = Laplace, and of course the solution u(x,t) at some point is influenced by u(x,t) on the entire 1+n dimensional boundary and that would include values of t in the "future and in the past". Stakgold is not asking why the arrow of time goes only in one direction, nor does he even mention this question. He just says that, given that it does, problems of the class 7.3 are associated with real world problems in heat flow. Regarding the time arrow direction, for heat and diffusion I think it is statistical mechanics and the notion of going in a direction that lowers the Gibb's free energy that determines time's direction, but I am very rusty on that subject. Application to problems with an Infinite Rod (200) Meaning of "a rod". A "rod" means that n = 1 for spatial dimensions. So a rod is something that is just a line or line segment along the x axis, say. An infinite rod is the entire x axis. Such a rod has no boundary σ at all. A semi-infinite rod running (0,∞) I call a "half rod" below, and it has a boundary consisting of the single point x = 0. A finite rod has a boundary consisting of the two endpoints. This 1D model can also apply to a "thin rod" of finite cross section if you say it is "insulated" so no heat can flow out the sides of the rod. A requirement here is that the rod is thin enough so there is no variation of u on a cross section of the rod. The validity of this assumption would have to be checked for a given problem. The time change of things has to be slow enough to get equilibration across the cross section at any time instant. No surface integral for n = 1 in 7.14. In n = 1 applications, recall that the Divergence Theorem has no integral on the RHS, it is just parts integration in 1D. ∫V dV A = ∫S dSA // divergence theorem = Gauss's theorem (1) Schaum 22.59 !Syntax Error, I dx ∂xA = ∫S dSA = ∫S dS A(x) = A(b) - A(a) The "surface S" here is just the two points a and b. So the divergence theorem is then "calculus". For that reason, 7.14 has no spatial integral on the RHS dSn , it just has the integrand sitting there. In the case of the infinite rod, not only is there no dSn integral, but there is also no boundary so the entire last term in 7.14 is not present. In (7.15) we get the full-rod problem over all of R1 with some t = 0 temperature distribution f(x). This is the problem whose Green's Function g = C we already solved for general n as shown in (7.10) and we now quote that g solution for n = 1 in p 200A. The full rod's cylinder picture consists of an infinite vertical sheet lying in the x-t plane starting at x=0. One could regard the two vertical edges at x = -∞ and x = ∞ as being the "surface of the cylinder", but these surfaces are ignored I think because the integrand of any contribution they might make would be zero. For example, as explained below, we would have for the vertical bounding lines at x=-∞ and x=∞: [ u(t0,x0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(t0,x0) ] = 0 when x0 = ± ∞ So in effect, the full rod has no boundary of interest, and therefore the last term in 7.14 can be neglected, so we then have this as the applicable (7.14), where u(0,x0)=f(x0), the initial condition of the infinite rod: u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(t0,x0) + !Syntax Error, Idx0 {u(0,x0) g(x|x0)}|t0=0 and if q = 0, we simply get u(t,x) =!Syntax Error, Idx0 f(x0) {g(x|x0)}|t0=0 but g(x|x0) is the C shown in p 200A and we set t0 = 0 in that so that C(x,x0)|t0=0 = 1/ exp(-(x-x0)2/4t) and then we have u(x,t) = 1/ !Syntax Error, I dx0 f(x0) exp(-(x-x0)2/4t) // agrees with (7.16) This then is the complete solution to our problem! [ this is the Hello World of heat conduction problems.] It is a superposition of free space propagators starting at x0 and weighted by f(x0). If you were to set f(x0) = δ(x0-x1), you get g. And if you take the limit t→0+ we get u(x,0+) = f(x) since the C propagator becomes a δ function (I proved this in Chap 5, and it seems pretty natural). Here is a plot of the solution when f(x) = δ(x-2) This is a picture of our famous "n=1 free space propagator" C(x,t: x0=2,t0=0). Along the t=0 line it is a delta function δ(x-2), and this delta then spreads out in x, and weakens as t increases. I think any student should have this picture in mind! Comment added 12.8.11. Suppose you start the entire rod at some temperature u0. Then what happens? u(x,t) = 1/ !Syntax Error, I dx0 u0 exp(-(x-x0)2/4t) = u0/!Syntax Error, I dx0 exp(-(x-x0)2/4t) = u0/!Syntax Error, I dx exp(-x2/4t) = u0/ ( 2 ) Maple = u0 hurray! Stak's comments about u depending continuously on the initial data. The idea there is that we want fk → f to imply solution uk→ u where k indicates some sequence. We usually write this as saying for given ε we can find δ so that |fk-f| < δ => |uk-u| < ε . You can maybe think of f as variable x here, and u = u(f) in a functional sense. Then if u is bounded over f, we know it is also continuous. You would ponder || u || ≡ maxf |u(f)| and he shows using p 200C that in fact || u || ≤ ||f|| so that yes, u is bounded as a function of f and therefore the solution u is a continuous function of f. Infinite rod with dipole impulse at the center. How now wants us to consider the special case not f(x) = δ(x-ξ) [ I just did this above] but f(x) = -δ'(x) which is the little dipole impulse at x = 0. This is little hot/cold instant shot on the rod at x = 0, I could imagine that it might do nothing. He blindly inserts this f(x) into (7.16), does the obvious parts integration, and ends up with (7.17). He claims logically that this -∂xC(x|0) for t>0 solves the homo heat equation. [ Recall in potential theory that a dipole source results in a solution u = ∂l E, so that is basically what is happening here as well, u = - ∂xC. ] So here is a picture of the infinite rod solution with f(x) = -δ'(x): Teapot tempest in the dipole problem. Now where does p 201 A come from? I think in the Stak mind, this is coming by taking the t→0 limit of the left expression in 7.17 for fixed x > 0 and u ~ t-3/2 exp(-x2/4t) → 0. He is "shocked" that we seem to get this limit u → 0 but I am not shocked because this limit operation does not include the x=0 region. He notes that if you go to t→0+ constrained by x = ±2 you get ±∞, so clearly things are not continuous in this limit. So at page bottom, he takes the careful distribution limit of (7.17) and in fact finds in this manner that u → -δ'(x) which is of course the f(x) we started with. This seems a tempest in a teapot to me. Comment: The solution to the infinite rod with the -δ'(x) initial condition is u(x,t) = -∂xC(x,t|0,0) as shown in 7.17. We showed above that C→0 as t→∞, so we get u→0 in this limit. This is what we expect from work shown below because the total heat in our δ' pulse is 0, and this just spreads out and the average temperature of the rod starts at a non-uniform zero and ends up at a uniform u = 0. [p 202 rewrite the general infinite rod solution] Now we are going to define z = (x0-x)/[2] dz = dx0/[2] [2]z = (x0-x) x0 = x + [2]z u(x,t) = 1/ !Syntax Error, I dx0 f(x0) exp(-(x-x0)2/4t) // this is 7.16 for the ∞ rod problem = 1/ !Syntax Error, Idz [2] f(x + 2z) exp(-z2) = 1/!Syntax Error, Idz f(x + 2z) exp(-z2) // agrees with 7.19 Now !Syntax Error, Idz exp(-z2) = , so you can "see" that the limit as t→0 is f(x) taken out of the integral and then the rest is just 1. I agree. Rest of section is fine. This last is just a restatement of our 1D rod problem solution 7.16. Exercise 7.1: Get (7.14) for The Neumann case ∂nu = h on σ. [ I really did 7.1 and 7.2 in my separate doc "BV problems...." , but leave them here as well. ] The slightly larger question here is this: mimicking the derivation of (7.14), how are things different if we insist not that g = 0 on σ, but that ∂ng = 0 on σ, sort of Neumann versus Dirichlet. Our BC this time is not that u = h on σ, but that ∂nu = h on σ. I traced the steps above, and the only difference is that we now keep "the other term" in the surface dS integral part, and this term has the opposite sign. The result we got above was this u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(x0,t0) + ∫ dx0 f(x0) g(x|x0)|t0=0 – !Syntax Error, Idt0 !Syntax Error, I dS0 [ u(x0,t0)∂nx0 g(x|x0) ] // agrees with (7.14) which I then modify as just discussed to get u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(x0,t0) + ∫ dx0 f(x0) g(x|x0)|t0=0 + !Syntax Error, Idt0 !Syntax Error, I dS0 ∂nx0u(x0,t0) g(x|x0) or u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(x0,t0) + ∫ dx0 f(x0) g(x|x0)|t0=0 + !Syntax Error, Idt0 !Syntax Error, I dS0 h(x0,t0) g(x|x0) // which agrees with Ex 7.1 [ See "BV problems for and Heat.doc" for the general case written out. ] Exercise 7.2: Get (7.14) for The Radiative case ∂nu = h on σ. This is the same idea, but we have the locally mixed BC (the one Sneddon calls radiation BC). Way back if we keep both surface terms we would have gotten - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = + ∫ dnx [ u(0,x) g(x0|x)|t=0] - !Syntax Error, Idt !Syntax Error, I dSn [ u(t,x)∂nx g(x0|x) - ∂nx u(t,x) g(x0|x)] We now have ∂nu + θu = h on the boundary σ, and we want then to use a Green's function which has the property ∂ng + θg = 0 on the boundary. So write [..] = u∂nx g(x0|x)- ∂nxu(t,x)g(x0|x) = u∂nx g(x0|x)- ∂nxu(t,x)g(x0|x)-θug+θug = = [ u∂nx g(x0|x) + θug] - [∂nxu(t,x)g(x0|x)+θug] = u [ ∂nx g(x0|x) + θg] - [∂nxu(t,x) +θu] g(x0|x) = u [ 0 ] - h(t,x) g(x0|x) = - h(t,x) g(x0|x) So the above then becomes - !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + u(t0,x0) = + ∫ dnx [ u(0,x) g(x0|x)|t=0] - !Syntax Error, Idt !Syntax Error, I dSn [- h(t,x) g(x0|x)] and we have then u(t0,x0) = !Syntax Error, Idt !Syntax Error, Idnx g(x0|x) q(t,x) + ∫ dnx [ u(0,x) g(x0|x)|t=0] + !Syntax Error, Idt !Syntax Error, I dSn h(t,x) g(x0|x) and we can then swap everything to get u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(x0,t0) + ∫ dx0 f(x0) g(x|x0)|t0=0 + !Syntax Error, Idt0 !Syntax Error, IdS0n h(t0,x0) g(x|x0) // result for Exercise 7.2 and this is the same form as the result to Exercise 7.1, but of course h has a different meaning in this case. Also, if we let θ→0, then we get exactly the result of Ex 7.1 in all respects. [ again, see doc quoted above where I do all the cases ] Exercise 7.3. Interpreting the BC's as additional source terms. The double Heaviside p 203 A is just fine. What is 2(uHR(x))? If we go back to page 14 Example 5, we can set = uHR so this example then says < 2(uHR),φ> = !Syntax Error, Idx φ 2u + ∫σ dS ( u ∂nφ - φ ∂nu) // this is exactly 5.16 = !Syntax Error, Idx φ HR 2u + ∫σ dS ( u ∂nφ - φ ∂nu) // φ(xσ) in the dS integral Now if we were to formally set the test function φ(x) to δ(n)(x-x1), where x1 is some point in R, then we could get a symbolic version of the above equation: -1 L-1 L-n L-n L-1 2(uHR(x1)) = HR(x1)2u(x1) + ∫σ dS [ u ∂n δ(n)(xσ-x1) - δ(n)(xσ-x1) ∂nu ]) and now replace x1 by x, some point in the volume, to get 2(uHR(x)) = HR(x)2u(x) + ∫σ dS [ u(xσ) ∂n δ(n)(xσ-x) - δ(n)(xσ-x) ∂nu(xσ) ]) Now we have to imagine somehow that ( n is number of dimensions, and n is the normal direction) δ(n)(x-xσ) = δ(n-1)(x-xσ) δσ(ξ) ∂nδ(n)(x-xσ) = δ(n-1)(x-xσ) ∂nδσ(ξ) ∫σ dS f(x) δ(n-1)(x-xσ) = f(xσ) where δ(n-1) is along the boundary surface and δσ(ξ) is perp to it, then somehow we get 2(uHR(x)) = { HR(x)2u(x) - [∂nu(xσ)] δσ(ξ) + u(xσ) [∂n δσ(ξ)] } This looks very close to p 203 B but he has a - sign where I have a + sign on the last term. It seems to me that the two surface terms always have a different sign as they do in 5.16. Ignoring this sign problem, it is then easy given p 203 B to get (7.21). He neglects to state what I think is his main point here: That if you "extend" the definition of u to v as shown in p 203 A, so that v is then defined in all space, not just inside the volume R, then you find that in addition to your normal driving term q which you see in p 203 E for the heat equation for u, you have a bunch of extra terms. These extra terms are "effective sources" for this problem viewed in the larger sense of v in all space. There are three extra terms that arise from the boundaries. These is an f(x) "initial condition" term, and there is a "Dirichlet boundary condition term" with h(x,t), and there is a "Neumann boundary condition term" which involves ∂nh(x,t) on the boundary. [204] He continues on p 204 writing (7.21) as p 204 A where p is our RHS with many terms of 7.21. He then uses the Green's function inversion idea to get B which he calls "one solution". This is a double integral of g times the RHS of 7.21. After some fiddling he gets v given by D, and I can then compare this to what I got on my "BV problems doc" u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idnx0 g(x|x0) q(t0,x0) + !Syntax Error, Idnx0 {u(0,x0) g(x|x0)}|t0=0] + !Syntax Error, Idt0 !Syntax Error, I dSn0 [ u(t0,x0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(t0,x0) ] (7.14 gen) and we have in effect derived this 7.14 gen result in an extremely unpleasant manner. His last term in D is really the one where I have a sign question and a general question. With his sign, it is this: + ∫dt0 ∫dVξ g(x|x0) h(ξ,t) [∂n δσ(ξ)] In some manner I don't buy, he is able to move this ∂n over to the g factor (somehow it misses the h) and this causes a sign change as we see in D. This entire Exercise is extremely wobbly I would say, and I think he would agree. But the point is clear: you can interpret (7.14 gen) as having three extra terms which are "sources" roughly on a footing with q(x,t) which arise from the initial condition f(x), from the Dirichlet boundary situation h(x,t) and from the Neumann boundary situation ∂nh(x,t). You can only "see" these things as sources if you enlarge your view to volume which contains R. I could imagine this demonstration being done in a much clearer and less hazy manner, perhaps in a paper somewhere. But the point is made. 7.3 Methods for finding Causal Green's Functions (all 1D space) [204] A. The Method of Images and other Trick Methods. The Half Rod. We now come to the "half rod" which is a little more interesting than the full rod. The "cylinder" in 1D requires a few words. At each time, the "cylinder" is a half axis (0,∞) whose boundary is the two points x = 0 and x = ∞, so you could say that the "surface of the cylinder" in this case is two vertical lines, one at x = 0 and one at x = ∞, while the cylinder itself is the flat vertical semi-infinite sheet connecting these two boundary lines. The boundary line at x=∞ is ignored on the grounds that any heat solution will be 0 there, as discussed above for the full rod, so the effective boundary is just the vertical line at x = 0. It is on this line that the boundary condition would be specified as h(0,t). Then on the horizontal line at t = 0 we would have the initial condition f(x) specified. How does our u(x,t) equation appear in this special 1D case? We start with our general u(x,t) u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idnx0 g(x|x0) q(t0,x0) + !Syntax Error, Idnx0 {u(0,x0) g(x|x0)}|t0=0] + !Syntax Error, Idt0 !Syntax Error, I dSn0 [ u(t0,x0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(t0,x0) ] (7.14 gen) As in the divergence theorem in 1D, the surface integral here is just a point at x0 = 0, so the above reads u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(t0,x0) + !Syntax Error, Idx0 {u(0,x0) g(x|x0)}|t0=0] + !Syntax Error, Idt0 [ u(t0,0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(t0,x0) ] (7.14 gen, half rod) and I now change from u(t,x) notation to u(x,t) to agree with Stak, and insert the f and h functions: u(t,x) = !Syntax Error, Idt0 !Syntax Error, Idx0 g(x|x0) q(t0,x0) + !Syntax Error, Idx0 f(x0) { g(x|x0)}|t0=0] + !Syntax Error, Idt0 [h(0,t0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(tx0,t0) ] (7.14 gen, half rod) [ p 204 The Half Rod Green's Problem, g=0 at end (image method)] For this Green's Function definition, see the system 7.22. We have just the first term of the last line above since g = 0 on σ, but we also set h = 0 as usual for a Dirichlet Green's system so the entire last line above is gone. Then f = 0 as well, so the second line is gone, and q = δδ so we just get g = g, at least things are consistent. In passing, note that for the half rod, the normal direction at that boundary is = -. So our statement of the Green's system is 7.22 where we have a δ(x-x0)δ(t) driving pulse in the first line. But now we suddenly consider a different problem: a full rod problem with the pulse above, and with the negative of that pulse at x = -x0 (an image or mirror pulse). We know how to solve this problem because we know the full-rod propagator which is C, so we superpose the two pulse solutions as shown p 205A. The claim is this solution, called w in B, satisfies all three lines of (7.22) and therefore, appealing to uniqueness, this must be the solution to the half-rod problem, and we then have (7.23). This is really an "image method" solution. Let's compare to the electrostatics problem of a point charge near a grounded plane where we only care about x ≥ 0, the side where our charge is located. We consider an alternate problem which is this point charge plus a negative point charge at the mirror point over in x < 0. For x > 0 this second problem's solution solves . For x = 0 the second problem solves the requirement V = 0 on the plane. But these are in fact the requirements of the first problem, so by uniqueness, this second problem solution part with x ≥ 0 (which is trivial to write) is also the solution of the first problem. We never think about what happens at x < 0. So good work, (7.23) is the Causal Green's for the half rod problem, our first non-trivial Causal Green's function computation! And here is a Maple plot for x0= 2 setting the 4π=1: Things to notice: (1) get decay at any x0 as t increases; (2) as the plotting routine gets close to t=0, we get closer to seeing the δ(x-x0); (3) the edge at x=0 is glued to the floor at all time as the BC requires. Of course if we plot for x<0 as well, we see our solution plus the image activity going on over there which is an artifact in terms of "our problem", but is the correct solution for the double impulse on the full rod. [ p 205 Half Rod with u=0 on the end (ie, h(0,t)=0) , but u(x,0) = f(x). ] Our formula from above then has just the first two lines, but we assume no volume sources so q = 0 and then we just have u(t,x) = !Syntax Error, Idx0 f(x0) { g(x|x0)}|t0=0] // this is 7.25 with 7.23 as g I now repeat the above obvious letting old notes remain: Our half rod has some initial temperature distribution f(x) along the rod, it could be anything reasonable. But, the left endpoint of the rod is going to be held at u=0 for all time, this is our σ BC. This would seem to require that f(0) = 0 to have things be consistent. Since we already know the half rod Green's, we jam stuff into (7.14) to get u(x,t). In passing, note that the dS0 integral in the third term there is missing for the 1D rod problem! That third term is just the time integral of ∂n0g times u = h(x0= 0,t) but in our example we are holding this at 0. So the point is that the third line is 0. Since no sources q, the first line is also 0, and we have only the second line left which gives (7.25) which is rewritten as (7.26). We can see that at x = 0, the two integrals cancel and we meet our boundary condition that u=0 at the left end for all time. And at t = 0+ we see the 2nd integral vanish and the first just gives f(x). Stak likes always to "check things" after he finds a solution. [ Stak shows this last fact after changing variables.] In the special case f(x) = 1, the solution is (7.27) with the error function erf(z) which Maple tells us looks like this: [ z = x/(2) ] so again we see that u = 0 when x = 0 for all time. Physical interpretation? The entire half rod starts off at constant temperature u = 1, but the left end is held against a temperature reservoir u = 0. As t increases, for any fixed x we have z decreasing, so u eventually drops below 1, below 1/2, and tends to 0. In other words, the cold source being held on the left end pushes its coldness down the rod. Eventually the entire rod will be at u = 0, but this till take an infinite amount of time. This is our first actual heat flow problem in this chapter! Here are some plots of (7.27). If we could really get to t=0, the sheet would be glued to u=1 almost all the way across the rod, as shown in the lower picture where t = (0,1) only. This also shows how we are glued down to u=0 along the line x=0. The f = 1 glue line is non-uniform of course. As usual, for large t everything decays to temperature u=0 since heat all leaks out the left end of the rod. Lower right picture repeats lower left but out to t = 10 instead of t = 1 with 10x more surface points, and is rotated so you can see that back of the surface as it falls away from f(x) = 1. Comments on the Error Function qua heat solution. The erf has this definition (omitting the overall constant (2/)) erf(z) = !Syntax Error, Idu e-u and I will now show that erf[(1/2)x t-1/2] is a solution of the homo 1D heat conduction equation. This takes a surprising amount of work. We start off with, ∂z erf(z) = e-z ∂z2 erf(z) = ∂z [∂z erf(z)] = ∂z[e-z] = e-z (-2z) Here ∂z can be regarded as either a partial or a full derivative. Now suppose z = z(s,w,u...) . Then we know that (here ∂s is a partial and this is a partial derivative chain rule) ∂s F[z(s,w,u...)] = (dF/dz) ∂sz(s,w,u...) which I write as (now ∂s is partial and ∂z is either partial or full) ∂s F(z) = ∂sz ∂zF(z) So I justify my sloppy notation only inasmuch as we can regard all ∂ as partials. Then: ∂s e-z = (∂sz) ∂z e-z = (∂sz) e-z (-2z) ∂s erf(z) = (∂sz) ∂z erf(z) = (∂sz) e-z ∂s2 erf(z) = ∂s[∂s erf(z)] = ∂s[(∂sz) e-z] = (∂s2z) e-z + (∂sz) ∂s e-z = (∂s2z) e-z + (∂sz) [(∂sz) e-z (-2z)] = e-z[ (∂s2z) - 2z (∂sz)2] If we set z = (1/2)x t-1/2 then ∂tz = -(1/4)x t-3/2 ∂xz = (1/2) t-1/2 ∂x2z = 0 ∂x2 erf(z) = e-z[ (∂x2z) - 2z (∂xz)2] = e-z[ 0 - 2{(1/2)x t-1/2} (1/2)2 t-1] = e-z[ - (1/4) x t-3/2] ∂t erf(z) = e-z (∂tz) = e-z[ -(1/4)x t-3/2] Since these terms are the same, erf(z) is a solution of the 1D homo heat equation. QED. Now, we might wonder: what is the general solution of (∂t-∂x2)u(x,t) = 0 ? I only know the theory of ODE's, not PDE's, so I don't know how to answer that question ! Stak has never taught us the theory of partial differential equations. I can see that a solution is u(x,t) = A + B erf[(1/2)x t-1/2], but that is probably not the most general solution. Stak's book is on "boundary value problems" and not " partial differential equations", so I would have to consult something like C&H on this question, which I won't do right now. It is on the list. [ p 206 The Above solution for large t. ] Really we are just doing a large t expansion of the error function (meaning a small z expansion) and our solution for the f = 1 half rod above is now (7.28). This tells us that a point on the rod reaches u = 1/2 when x/(2) = 1/2, or t = x2. If a point of interest is twice as far down the rod from some other point, it takes 4 times more time to act on that point. [ p 206 Time method for half-rod with f = 1 ] If we take our heat equation (7.24) and replace t with s using the usual Laplace transform, the PDE becomes (7.29) and I verified on scratch (and on page bottom marked by α) that the solution of (7.29) is (7.30). Now if we wanted to verify that we get the (7.28) erf for large t by this Time Laplace Method, we would have to examine the INVERSE Laplace Transform to that shown in p 206 A, evaluated for large t. Stak expends 4 pages in his Appendix B figuring out what the large t expansion is for the inverse of a Laplace Transform. He does this by looking at the Mellin-Barnes vertical contour, moving it to the left, and studying the pole and cut residues it picks up. I will write (have written) notes on this Appendix B elsewhere, but the end result for inverse is (B.10) given some f(s) that we are inverting. So presumably, if we use this (B.10) with f(s) = (7.30), we will recover the large t expansion (7.28). [see Schaum page 169 32.110 for the inverse being erf! ] He is showing us this because there will no doubt come a time when our only method of solving a problem is to do this Time Laplace transform, obtain f(s) = his (x,s), and it will be so messy that we cannot look up the inverse Laplace, so we will use (B.10) to get an asymptotic expansion for our solution u(x,t). Notice that the Time Laplace method converts our PDE in x,t into an ODE only in x, and that ODE is in fact the 1D Helmholtz. With more space dimensions, this is still true: our heat conduction equation becomes the inhomo (driven) Helmholtz equation BUT we have k2 < 0. In the case of (7.29), we have k2 = -1 and driving function 1. I am not sure this thing is called Helmholtz when k2 < 0, but maybe. [ p 206 The f=1 half rod solution for small t ] This is a much easier expansion, and the result is the expansion page 207A. So in this example, he has shown the reader how to get both large t and small t expansions for the result, and if the result is too complicated to evaluate, this would be very useful. [ p 207 The f=0 half rod with left end driver h(t). ] We apply our general half-rod form above to get: u(t,x) = !Syntax Error, Idt0 [h(t0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(tx0,t0) ] (7.14 gen, half rod) Physically this is a "frozen half rod" (u=0) and we apply a time varying heat source at the left end. We expect this heat to wander to the right in some reasonable manner. Looking at our generic heat solution (7.14), we see that both the first and the second term now vanish (since q = 0 and f = 0), so the entire u(x,t) will come from the third term which contains ∂ng. As I noted above the integral dS0 is absent in this 1D spatial problem and σ = just the point x0 = 0. We know g for our half rod problem, so we have to compute ∂ng need for this third term of (7.14). So he computes ∂ng in p 207 B (I verified it fully). We stick this into (7.14)'s third term and set x0 = 0 since that is σ, and our solution is (7.35), expressed as an integral over h(t0) dt0. He may do it soon, but I could imagine a time pulse h(t0) = δ(t0 - 0) at x=0 and then we seem to get u(x,t) = x exp(-x2/4t)/[ t3/2] which would describe something moving to the right down the rod. This is not C because there is an extra power of t in the denominator. Here is a plot Our end-driven pulse spreads to the right on the rod, but also decays fast in time. The right picture shows the spread out to the right at early times before the decay sets in. The entire left front edge was at 0 just after the spike hit, bit we are too late to see that. Stakgold would have loved these plots and I am sure that current day books on the subject are loaded with such plots, but I don't know the titles of those books. Note added 5.20.11. What is the solution of our frozen half-rod problem if h(t,0) = 1? Well, it is the integral shown in 7.35 where set h(t0) = 1. We need then to do this integral !Syntax Error, I dt0 (t-t0)-3/2 exp(-a/[t-t0]) a = x2/4 Maple does this integral as follows so our solution 7.35 becomes u(x,t) = - (4π)-1/2 x / (x/2) * [ -1 + erf(x/2)] = [ 1 - erf(x/2) ] At t=0 we have erf(∞) = 1 so u = 0 everywhere, the initial frozen rod. At t = ∞ we have erf(0) = 0 and the entire rod is then heated up to u = 1 and that is the steady state solution. [ p 207 taking the t = 0+ limit of the above half rod h(t0) = δ(t0 - 0) solution.] We want to show u(0,t) = h(t), but we naively seem to get u = 0 from 7.35. But he shows that the integrand is really a delta sequence so in fact the correctly done limit gives h(t) as shown in p 208 B. Similar to teapot thing above. [ p 208 a Subtraction Method for the half rod ] We saw that the x→0 limit of our driven cold rod problem above is a bit tricky and we had to use that delta sequence. Here he shows a trick to make the solution near x = 0 be more "computable". The generic idea is to somehow expose the singular part. So instead of dealing with the solution u(x,t), he suggests writing it as u(x,t) = k(x,t) + v(x,t) where k(x,t) has the problem, and v(x,t) is well behaved. He uses this subtraction method for our cold driven half rod problem (but only with h(t) = 1) and gets (7.38) where k = 1 and v = the error function shown. Notice that v is smooth at x = 0 (see 2D erf Maple plot above). This example does not quite demonstrate the idea of exposing the violent part of a solution, but it does demonstrate the method with an actual simple example. [ Well, maybe u(x,t) was 0 at t<0 before it instantly jumps to 1, whereas v(x,t) is continuous at t = 0.? ] [ p 209 the half rod with insulated end, Neumann ∂xg = 0 (image method) ] We now have u(t,x) = !Syntax Error, Idt0 [h(t0){-∂nx0 g(x|x0)} + g(x|x0)∂n0u(tx0,t0) ] (7.14 gen, half rod) but we have a different Green's function g in this case. This problem is similar to the electrostatics problem of a point charge at (x,0,0) next to a planar surface on which ∂xV = 0. We know that in this case the image or mirror charge has the same size and sign as the Green's charge. Such a potential has to be symmetric in x, and so must have ∂xV = 0 at the surface x = 0. And it solves on the right. We appeal to uniqueness to say that the two charges gives the correct solution on the right including on the plane approached from the right. Just so, in our heat problem, we take our original heat pulse and add an image one on the left, same size and sign, so that ∂xg = 0 at x = 0, which is the desired BC. So this means that the solution is g = C(x0) + C(-x0) as shown in (7.39), whereas our problem with g = 0 at the half rod end was given by C(x0) – C(-x0), see p 204 above. The plot now of g is this You can see that, since the left end x=0 is insulated, the temperature there rises at first, then slowly decays after that. If you look at the t slices, we start with our δ(x-2), this spreads out as a Gaussian and raises the left end x = 0 as just stated, but then things decay, but it takes a long time for x=0 to decay because it holds the heat since it is looking at an also-hot abutting section. Unlike on the right, heat cannot flow away very well from the left end of the rod! These are fascinating plots I think, they say a lot. [ p 209 Green's for the half rod with radiative BC ∂ng+θg = 0 at the end.] Our image method does not work for this problem which is a mixture of those previous two problems, but he shows a "trick" to solve the problem. He defines a new Green's function v in terms of the desired g as in p 209 D (very Sneddon-like) and writes up the BV problem for v. But it looks just like 7.24 which was our half-rod with f(x), only here we have a specific f(x). Since the solution to 7.24 was 7.25, he just inserts this special f(x) into 7.24 and gets p 210 A,B, then C with D as usual. In E he has a tricky solution of p 209D for g in terms of v, and then he grinds it all through and writes out the resulting radiative BC Green's Function three different ways, 7.42,43,44 all on page 210. Not rocket science. He claims we shall see this same problem done in a simpler manner soon. Let's do a few steps: v = (∂xg - θg) ∂tv - ∂x2v = ∂t(∂xg - θg) - ∂x2(∂xg - θg) = ∂t∂xg - θ∂tg - ∂x3g + θ∂x2g = ∂x [∂tg - ∂x2g] - θ[∂tg - ∂x2g] = (∂x-θ) [∂tg - ∂x2g] = (∂x-θ) [0] = 0 which is the first equality in E, so v satisfies the homo heat equation. Then the BC at x=0 is this v(x,0|x0,0) = ∂xg(x,0|x0,0) - θg(x,0|x0,0)) = ∂x δ(x-x0) - θ δ(x-x0) and this then verifies all of the system shown in p 209E. But we know the solution to system E and it is 7.25 which we write out in p 210 A. We then do parts on the ∂ξδ(...) term, and then change ∂ξ to either ±∂x depending on which term, and I agree with the four terms B + C with D. Stak then solves the little ODE for g which is driven by v as shown in E. I have shown that applying ∂x to both sides of E gives what we want. Then in F+G he sticks in the result B+C but x is replaced by dummy int variable α. At this point I stop doing the algebra and calculus and I just accept the last three equations stated on page 210. The point is clear: We convert this radiative BC Green's system in g to a Green system in v whose solution we know, then we grind through that known solution to get v. This then drives the ODE for g, and from that we find g. This is the "trick". Now where does he solve this same problem in a simpler manner as claimed bot page 210? Well, he does this with θ = 1 as Example 2 on page 217 using the spatial EF method with continuous λ, and this involves the no-name transform shown there. See below. [ p 211 Green's for an insulated wire ring ] Perimeter is 1 unit, wire is regarded as really a 1D problem (very thin wire, and it has insulated sides, so like a piece of our infinite rod that is just cut and bent into a circle.) The BV problem for this ring is shown in (7.45) which seems very reasonable: this causes in effect the solution and all derivatives to be continuous at the linkage point which is ± 1/2 . If you look at the picture on page 212, you see that we are in fact going to form our wire ring from the piece of straight wire that runs along (-1/2,1/2) on the rod axis. Now consider a different problem, one where we put a positive identical source at x = 0,±1,±2, etc. The solution to THIS problem will be periodic on the infinite rod, and therefore everything in (-1/2,1/2) repeats in every unit interval, and that means that all derivatives at +1/2 will match those at -1/2. Thus, the (-1/2,1/2) segment of this infinite rod satisfies the BV problem specified for our ring! So our ring solution must be simply (7.46). We next come to a sort of technical problem in writing the solution. We start with 7.46 which is just the sum of all the C functions, one for each image source. The problem is that series 7.46 is clumsy to deal with and Stak uses one of his tricks to rewrite this as 7.48 which is much more reasonable. We met this trick back in Volume 1 as he notes, so I am not going to worry about it here, and I just accept the result. The trick is to replace one sum with another sum that is nicer. In particular, you can see that 7.48 is much nicer at large t, and the periodicity is explicit in the cos(2πnx) factor. I expect that as t→∞, result 7.48 should approach a constant around the ring, but that limit is not obvious to me, I would have to do it somehow. Here is a plot using the first 200 terms. A finite term count causes ringing on the t=0 edge as you might expect, but basically we start with a delta function and it just spreads out till it is even everywhere at some positive average value. (this is an insulated ring). You can see this rise in the smaller picture on the right. [ p 212 Green's Function for a finite rod of length l with g=0 at each end] As with our ring, we shall regard this finite rod as the segment of an infinite rod running (-l/2,l/2). In the previous problem where we superposed an infinite number of + image sources, we obtained matching function and derivatives at the segment boundaries, but here we need to make the function actually be zero at all the boundaries. This is very similar to our solution to the 2D vertical strip problem (or width a, not l) in electrostatics which was exercise 6.35 page 166 and where we had this picture (circled charges are negative, others positive) In this electrostatics strip problem, we placed the charges as shown to get g = 0 on BOTH edges of the strip, and to make that work we needed an infinite set of charge pairs, and in fact this made g = 0 on an all the vertical lines obtained from those shown above by shifting those lines ±na. Just so, in our finite rod problem we consider a full infinite rod with a set of point sources exactly as above, and he shows this in a picture on page 212. Since the full rod then has a solution with u = 0 at all the seams, by the usual uniqueness argument, we have solved our finite rod with u=0 on the ends problem. This time, then, our solution is page 213 A which is similar to (7.46) but we have two terms corresponding to our two sets of image sources. We expresses this first in terms of the θ function of 7.46, then claims that some simple trig and fiddling (using 7.48) takes this solution to the form 7.50. In this result, you can see that g in fact vanishes when x = 0 and x = l . This solution is amazingly similar to the thin strip solution which I quote from the Chap 6 meta notes: g(x,y |0,y') = Σn=1∞ (1/nπ) e-(nπ/a)|x| sin(nπy/a) sin(nπy'/a) // = 6.127 B. The Method of Expanding in Spatial Eigenfunctions: h=0, initial f(x) (213) Here 7.51 is our "cold boundary" special case of the heat BV problem. Remember that we have two boundary values really: (1) the initial condition, which here is still f(x) ; (2) the boundary value which is usually h(t,x) but here h = 0. We can as usual try u(x,t) = X(x)T(t) with λ the separation constant. We end up with 7.52 which is our n dimensional spatial EV problem (λ the eigenvalue) where quantization is forced by h = 0 on the spatial boundary, the usual EV situation. For finite region, there will be some discrete spectrum of λi and he orders them at the bottom, allowing as there might be some multiplicities. The spatial EF's as usual are called φi(x), and the atomic form is then φi(x)exp(-λit). The decay sign here is not something we have selected, it arises because the λi are positive. We would never expect expo growth in a non-sourced heat problem! A Smythian form is then 7.54 with coefficients ci. Then the initial condition lets us compute the ci from the boundary condition f(x) as in 7.55, and the grand finale is the solution 7.56 to our problem! So basically we have solved the "cold boundary" heat flow problem with no sources and some initial f(x) by first solving the related spatial EV problem for the cold bounded region, where cold means u = 0 on the boundary. Now how do we explain the nature of the generic solution? Our spatial boundary has u = 0, and that means ∂nu will come out "as it may" and will not be zero. Since the boundary is at u = 0, you can think of the region R as embedded in an infinite reservoir of ice at u = 0 degrees temperature, and basically this is going to eventually suck all the heat out of our region and leave the region at u = 0 throughout. This is then the explanation of the decaying time exponentials in generic solution 7.56. Fascinating. [ 214 Causal Green's for the above.] What is different for Causal Green's looking at 7.51? Comparing to 7.8a, if we set f(x) = δ(x-xo), then our u(x,t) solution becomes the Causal Green's g. Then from 7.55 we are going to have ci = φi(x)* and then 7.56 becomes p 214A where so far we have had t0 = 0. But we have a time-shifting rule 7.9 (p 198, obtained from the BV problem equations), so 7.57 is then our Green's Function, the famous bilinear sum formula we saw earlier in Stak [ where?] . We can then of course use that in our general 7.6 to solve any heat problem in our class of problems where the boundary σ does not change over time. [ Remember: this section is "methods of finding Green's Functions" and so we have here a new method which uses separation of variables and spatial eigenfunctions. ] [ 214 Repeat the above for Neumann case ] This changes 7.52 to say ∂nX = 0 on R. So we now have a different EV problem on the same region, new eigenvalues are called μi and new EF's ψi(x). We know that the smallest of these non-negative numbers is μ0 = 0 as will be shown in a moment. The generic solution for u is p 215A (similar to 7.56). As t→∞, all the time expos decay except that for i = 0 and we find that u(t=∞,x) = the spatial average of f(x). The reason is that ∂nu=0 on σ means region R is fully insulated, so no heat can leave, it just spreads around and in the end we get the average of the initial temperature distribution f(x). So notice that unless <f(x)> = 0, you will always have d0≠0, if we define di as the coefficients for this problem (which Stak did not do). As for the corresponding Causal Green's, it must be 7.57 with φ→ψ and λ→μ, but Stak does not mention this. Comment: From the last two problems, we have sort of indirect evidence that (1) for the first problem, the spectrum will be discrete and entirely positive; (2) for the second problem, same thing with the addition of an EV at 0. If these facts were not true, we would have exponential growth in our heat solutions which is certainly "unphysical". [ 215 R unbounded.] The two spectra just mentioned become continuous, but the same (1) and (2) will be true. The Σi forms will become ∫dλ forms somehow, but Stak does not want to give a general theory for this situation, he will do individual problems soon. [ 215 The general u(x,t) problem by a non-EF method] We already know how to do the general problem, meaning u = h(x,t) on σ instead of u = 0, by taking our Causal Green's 7.57 and jamming it into our general formula 7.14, and now we are allowed to have sources q, initial f, and boundary h. [ I need to rethink 7.14 in the case of ∂nu specification on the boundary. DONE! ] Here, "for variety", Stak will do a direct EF attack on the general problem 7.3 instead of using g. In this direct attack, we get 7.54 but now we have time varying coefficients because we have X = h(x,t) on σ. Thus we have 7.58. [ This reminds me of time dependent perturbation theory where φi would be EF's for no perturbation.] Here, we could for example select our complete set φi to be the same one we used for h = 0, and I think that is what 7.58 is saying. Stak comments that yes, the φi = 0 on the boundary, but u ≠ 0 on the boundary, so expansions like 7.58 are going to have problems near the boundary, which we associate with "non uniform convergence", I have notes on that subject elsewhere. For this reason, Stak steers clear of doing anything term by term without perhaps integrating first. Our first step is this: apply ∫R dx φi(x)* to the PDE 7.3, ∂tu(x,t) – 2u(x,t) = q(x,t) 7.3 so ∫R dx φi(x)* ∂tu(x,t) – ∫R dx φi(x)* 2u(x,t) = ∫R dx φi(x)* q(x,t) Define ∫R dx φi(x) u(x,t) ≡ ui(t), and similarly for q, then we have ∂t ui(t) – ∫R dx φi(x)* 2u(x,t) = qi(t) // which is page 215 C Now a lot happens very quickly here! When we throw over the 2 to the φi* factor, we pick up the usual surface terms as shown p 216 A. This is my Green's #2 ∫V dV [ ψ 2φ – φ 2ψ ] = ∫S dS [ ψ( ∂φ/∂n) – φ( ∂ψ/∂n) ] (6) Green #2 ∫V dV [ φi* 2u – u 2 φi* ] = ∫S dS [φi* ( ∂u/∂n) – u( ∂ φi*/∂n) ] (6) Green #2 The first term on the right vanishes because recall the φi(x) vanish on σ ! Also, we know that 2φ = -λφ so we can treat the thrown over term easily. In the other boundary term, we have u = h. So we end up with 7.59. Assuming you have already solved the h = 0 problem and you know all about the φi(x) and the λi, you see that 7.59 is just a simple first order ODE in variable t for ui(t) where we know everything in the equation. And we know ui(t=0) as shown in p 216 B. So consider ∂tu(t) + λu(t) = r(t) To solve this, change from u(t) to v(t) using u(t) = e-λtv(t) and you get dv/dt = eλtr(t) => v(t) = !Syntax Error, I dt' eλt'r(t') + v(0) so solution is u(t) = e-λtv(0) + !Syntax Error, I dt' eλ(t'-t)r(t') => u(0) = v(0) which we know is fi from B and this last result is exactly p 216 C. So we have solved for ui(t) and we can throw in the rest of 7.58 to get our final solution u(x,t) by summing everything Σi φi and this gives the mess 7.60. To wit ui(t) = e-λit fi + !Syntax Error, I dt' eλi(t'-t)ri(t') from above u(x,t) = Σi φi ui(t) = Σi φi e-λit fi + Σi φi e-λit !Syntax Error, I dt' eλit' ri(t') Then we note that ri(t') = qi(t') – hi(t') and result is u(x,t) = Σi φi e-λit fi + Σi φi e-λit !Syntax Error, I dt' eλit' qi(t') – Σi φi e-λit !Syntax Error, I dt' eλit' hi(t') where of course λi means λi in an exponent. I can improve my notation Ab = e-λt from sym doc u(x,t) = Σi φi e-λt fi + Σi φi e-λt !Syntax Error, I dt' eλt' qi(t') – Σi φi e-λt !Syntax Error, I dt' eλt' hi(t') u(x,t) = Σi φi(x) e-λt fi + Σi φi(x) e-λt !Syntax Error, I dτ eλτ qi(τ) – Σi φi(x) e-λt !Syntax Error, I dτ eλτ hi(τ) and this is a more respectable version of (7.60). So good work Mr Stak, a complete solution of the general problem 7.3 is provided by the EF method where we use the EF's of the homo BC problem! This was an excellent "for variety" exercise, A++. [ Notice that double integrals are required in each of the last two terms in the expression for u(x,t) above, rather ugly. But if you did the Green's Function approach, you would have terms with one integral and one infinite sum, also ugly. But at least you have an answer in each case! ] Question: since φi(x) = 0 on the boundary, how do we get u(x,t) = h(x,t) on the boundary? It must come somehow from that last term. But I don't see any way that can really happen! The last term has this form for some fixed t ( for x very close to the boundary probably only this term survives) u(x,t) = Σi=1∞ φi(x) Hi(t) and we want limx→0 u(x,t) = limx→0 Σi=1∞ φi(x) Hi(t) = limx→0 limn→∞ Σi=1n φi(x) Hi(t) Recall Moore's Theorem (from order interchange folder) which says that both limits must exist and at least one must exist uniformly, and only then can you reverse the order of the limits (which, if you could do it, would say u(x,t) = 0, which is something we don't want!). So consider limn→∞ Σi=1n φi(x) Hi(t) I presume this limit exists certainly for x away from the boundary. But I suspect it is a non-uniform convergence, and I'll bet the other limit is also non-uniform. I don't want to get sidetracked on this issue right now, so will save it for that famous rainy day that will likely never arrive. The bin is getting fuller. C. Examples of the Eigenfunction Expansion Method (216) Example 1: Green's Function for finite rod by EV method u = 0 at both ends. Here we just apply our little theory as shown, and at once we get the same result we got earlier by the fairly complicated infinite image method. That is to say, p 217 A agrees with 7.50. Example 2: The half-rod with radiation from left end (a continuous λ problem) Radiation explained: OK, I finally get it. Recall Sneddon did this as well. Stak pointed out in 7.4 that you can model any heat loss at a boundary where you have Δu by the linear model ∂nu = -θ Δu where θ is some appropriate constant. Certainly we do this for conduction (θ = R value for glass), but he says this often applies to radiation and convection as well. If this happens at a point on some surface, then we would say that ∂nu = -θ(u - u0) where u0 is the ambient temperature outside our surface. Taking this to be "the frozen sea" u0 = 0, we get ∂nu = -θu, which is to say, ∂nu + θu = 0, and THIS is what people keep calling the radiative boundary condition. For our half-rod problem, ∂nu = -∂xu , and we just set θ = 1, and then our BC at the half-rod's end at x = 0 is ∂xu(0,t) - u(u,t) = 0 so we have 7.61. Now, the half-rod is unbounded, and this will be a continuous λ problem. What is the relevance of Lλ = δ(x-ξ) p 217B? Suddenly, in the middle of Chapter 7, Stak rolls out this equation on page 217 (-2 - λ)r(x|ξ;λ) = δ(x-ξ) or Lλ r(x|ξ;λ) = δ(x-ξ) // along with the radiative BC. My question is: how does this "tie in" with the Chapter 7 material we have been studying for a week up to this point? I am completely missing the point. First, what do we know about this in Potential Theory? This operator Lλ is one we studied a lot on Chapter 4, starting for example on page 259 with the string example. Even though the string problem has discrete EV's we call λi, in our operator Lλ the λ is a continuous variable. We can of course compute the EV's of the problem Lλiφi= 0 and these are the φi(x) as shown page 260 A. Our inhomo equation is Lλu = f as in 4.1 with λ continuous. We can expand u and f on the φi getting fi and ui and we think of the φi as defining a spatial transform. Our problem Lλu = f is then solved as in 4.3 and then 4.4 as long as λ is not an eigenvalue. This is all well and good, but it does not explain the relevance of the equation Lλu = f . I shall ignore this relevance issue for the moment and continue onto page 261. We see there the Green's Function problem written now as Lλu = δ and we consider this δ a special case of f, and we just apply 4.4 with f=δ and the result is 4.8 which is the usual g(x|ξ;λ) = Σi φi(x)* φi(ξ)/(λi-λ). Yes, we can integrate this on dλ to obtain an expansion for δ(x-ξ), so that is at least useful. We move now to page 262 where it seems Stak comes up with a simpler expression for g(x|ξ;λ) which is not an infinite sum but is 4.10. Notice that this is only our traditional trapezoidal string g when λ=0, as he shows. Now comes a question: does Stak do anything else of interest with object g(x|ξ;λ) ? Well, he again does the dλ integration with this new and simpler form and obtains again the δ(x-ξ) expression 4.12. So basically, the main purpose of g(x|ξ;λ) seems to be that it generates a transform and the dλ integration gives us the δ(x-ξ) for that transform, and that result is δ(x-ξ) = Σi φi(x)*φi(ξ) . This result of course does not contain the symbol λ anywhere (it contains λi). The dλ integral wraps around the eigenvalues λi (or the cut in the continuous spectrum case). The other part of that transform is the usual orthogonality statement for the φi(x). We then write f(x) = Σifi φi(x) and fi = ∫dx φi(x)*f(x) and that is our transform. So once again, even without thinking about g(x|ξ;λ) , we would know that the φi(x) of Lλiφi= 0 form a complete orthonormal set, so we know δ(x-ξ) = Σi φi(x)*φi(ξ) anyway. But perhaps consideration of the function g(x|ξ;λ) is giving us a PROOF that in fact the φi(x) are really a complete set, which is embodied in that completeness relation δ(x-ξ) = Σi φi(x)*φi(ξ). Apply this to Volume II page 217 for his "second method". Equation B says Lλ r(x|ξ;λ) = δ(x-ξ). It happens that we have the radiative BC here instead of the usual string BC, but I don't think that is the key issue. As we did in the string problem, we can come up with a simple form for the Green's Function r(x|ξ;λ). It is shown in C which is analogous to the "simpler expression" 4.10 for the string problem mentioned in the section above. We can then use this r(x|ξ;λ) to obtain a statement of completeness as in D and E, and yes, this implies the transform of 7.62. I have no problem with any of this. The question is: what is the relevance of all this stuff to the general flow of Chapter 7 ? The relevance must be 7.65 which is an expansion for the function g(x,t|ξ,0), so where is this 7.65 coming from? [Stakgold page 218] He has a typo here which caused me a lot of digression time, he means 7.57 and not 7.47, and 7.57 is this g(x,t |x0,0) = Σi φi(x)*φi(x0) e-λt and OK, the limit for a continuum would be g(x,t |x0,0) = ∫dλ φλ(x)*φλ(x0) e-λt λ = ν2 and THIS then is where 7.65 comes from. So what was the purpose of the Lλ equation shown as p 217 B? We saw above in the discrete spectrum case that the g(x|ξ;λ) object was used so that the dλ integration gave a completeness result and a transform. All this stuff had just λi as a discrete label. We found for example that δ(x-ξ) = Σi φλi(x)*φλi(ξ) . In the continuous spectrum case, we end up with something more like δ(x-ξ) = ∫dλ φλ(x)*φλ(ξ) which is what page 217 E says. Because we are continuous now, the symbol λ DOES appear in this result as the integration variable. So you could argue that the purpose of pondering equation p 217 B is to show that the functions φλ(x) of (7.64) form a complete set in the continuum sense. That is, to show p 217 E. Notice that we never talk about orthogonality, but that must also somehow be true. It is as if completeness alone is all you need to say you have a transform. Oh yes, by using this g(x|ξ;λ) to find δ(x-ξ), we are above to identify the φλ(x)!! This is then a method of finding the eigenfunctions in the first place! Another purpose in computing g(x|ξ;λ) is that the singularity structure of g(x|ξ;λ) determines the spectrum of Lφ = λφ (with BC's). What is the "first method" Stak refers to on page 218 top? This first method is the method Stak used on page 215 where we used time-dependent coefficients. In the continuum λ context, we would write u(x,t) = ∫dλ uλ(t) φλ(x) uλ(t) = ∫dx φλ(x)* u(x,t) // 7.58 new where the φλ(x) are eigenfunctions of -2 φλ(x) = λ φλ(x) with our radiative BC's. (These functions are shown in 7.64 as we now know). If we take our heat equation (∂t-2)u(x,t) = q(x,t) and apply ∫dx φλ(x)* to both sides, we get ∫dx φλ(x)* (∂t-2)u(x,t) = ∫dx φλ(x)* q(x,t) ∂t uλ(t) - ∫dx φλ(x)* 2u(x,t) = qλ(t) (*) which is like page 215C. We then have to use Green's Theorem as shown page 216 A to move 2 to act on φλ , but now things will be different because we have the radiation BC instead of the usual u=0 BC. The Green's Theorem will look just like page 216A, but the surface term requires a different evaluation: ∫dS [ φλ(x)* ∂nu - u∂nφλ(x)* ] We do our usual radiative trick and write φλ* ∂nu - u ∂nφλ = φλ* ∂nu + φλ*u - φλ*u - u ∂nφλ* = φλ*(∂nu +u ) – u (∂n φλ* + φλ*) But we know that (∂nu +u ) = 0 on our surface so our surface term ends up being - ∫dS u (∂n φλ* + φλ*) and we might call u(x,t) = h(x,t) on the surface of R, so this thing is then - ∫dS h(x,t) (∂n φλ*(x) + φλ*(x)) and we then have ∫dx φλ(x)* 2u(x,t) = ∫dx 2φλ(x)* u(x,t) - ∫dS h(x,t) (∂n φλ*(x) + φλ*(x)) = -λ ∫dx φλ(x)* u(x,t) + ∫dS h(x,t) (∂x φλ*(x) - φλ*(x)) = -λ uλ(t) + ∫dS h(x,t) (∂x φλ*(x) - φλ*(x)) so our result (*) from above is now ∂t uλ(t) - ∫dx φλ(x)* 2u(x,t) = qλ(t) ∂t uλ(t) - {-λ uλ(t) + ∫dS h(x,t) (∂x φλ*(x) - φλ*(x))} = qλ(t) ∂t uλ(t) + λ uλ(t) + ∫dS h(x,t) (∂x φλ*(x) - φλ*(x))} = qλ(t) ∂t uλ(t) + λ uλ(t) = qλ(t) - ∫dS h(x,t) (∂x φλ*(x) - φλ*(x))} // 7.59 new and this is our replacement for p 216 7.59. This is a first order thing in uλ(t) which we would solve in the same way he did before on page 216. Now suppose we want to find the Causal g and not u. We replace q→δ(x-x0)δ(t) and then qλ(t) = ∫dx φλ(x)* q(x,t) = ∫dx φλ(x)* δ(x-x0)δ(t) = δ(t) φλ(x0)* so I guess then we have ∂t gλ(t) + λ gλ(t) = δ(t) φλ(x0)* - ∫dS h(x,t) (∂x φλ*(x) - φλ*(x))} where now g(x,t|x0,0) = ∫dλ gλ(t|x0,0) φλ(x) gλ(t|x0,0) = ∫dx g(x,t|x0,0)* u(x,t) // 7.58 new and so then gλ(t|x0,0) is what Stak means on page 218 by "the transform of the Green's Function". Recall from earlier Stak notes on page 216 that we had ∂tu(t) + λu(t) = r(t) => u(x,t) = Σi φi e-λit fi + Σi φi e-λit !Syntax Error, I dt' eλit' ri(t') u(x,t) = ∫dλ φλ(x) e-λt fλ + ∫dλ φλ(x) e-λt !Syntax Error, I dt' eλt' rλ(t') and now we have from above rλ(t) = δ(t) φλ(x0)* - ∫dS h(x,t) (∂x φλ*(x) - φλ*(x))} so the δ(t) will present no problem here, and as before, we get our problem solution u(x,t) as double integrals of the functions f, q and h. I don't need to do this in full detail! I would say that the second method is a better way to solve this problem for g. Comment: When we did the half-rod with u=0 or ∂nu = 0 at x=0 ( p 204 and 209), we were able to use the method of images and the resulting Causal Green's function g in each was very simple, a difference or sum of C functions, and I was able easily to plot these results. But the radiative case has a much messier Green's function which we see in (7.65) and for which we found another result (7.43) which does not involve any integrals and contains the erf function. So I could in fact plot this if I wanted. D. The Transform Method (218) The idea is to replace time t by s, so that ∂t becomes s and we end up with a purely spatial problem where s is just a parameter. The function u(x,t) becomes (x,s) where twiddle means time transform. We have our usual initial condition with f(x), and our usual boundary condition with h(x,t). The usual heat equation becomes ( -2 + s) (x,s) = (x,s) + f(x) // includes the "initial condition" f(x) 7.68 (x,s) = (x,s) for x on σ // this is the "boundary condition" Stak applies this to his first 3D problem in this chapter. Up to now it has all been rod objects. Now at once we are faced with our Lλ equation fretted about so much above, but here it is for a different reason. We set up our G equation like this ( -2 -λ ) G(x|ξ; λ) = δ(x-ξ) 7.69 G(x|ξ; λ) = 0 for x on σ and he adds that G should be square integrable on R Let's assume that we go off and solve for G(x|ξ; λ). The appropriate Green's Function for our s equation above is then G(x|ξ; -s). From our Green's Theorem of regular potential theory we know in general that we can solve the u equation like so ( see "BV problems....doc") (x,s) = ∫V dξ G(x|ξ; -s) [(ξ,s) + f(ξ)] – ∫S dSξ [(ξ,s) ∂nξ G(x|ξ; -s) – G(x|ξ; -s)∂nξ(ξ,s) ] and for our assume boundary condition we have G = 0 on the boundary, so (x,s) = ∫V dξ G(x|ξ; -s) [(ξ,s) + f(ξ)] – ∫S dSξ [(ξ,s) ∂nξG(x|ξ; -s) // agrees with 7.70 The final step is to use the inverse transform 7.71 to get u(x,t). Now we want to do a special case of the above u equation like this ( -2 + s) (x,s) = δ(x-x0) // includes the "initial condition" f(x) 7.68 (x,s) = 0 for x on σ // this is the "boundary condition" which is to say we have set q and h to 0 and f(x) = δ(x-x0). The solution now is (x,s) = G(x|ξ; -s). Back in the time world this would be g(x,t | x0,0), so we must have 7.72 using the inversion formula. Of course we know that G(x|ξ; -s) = Σiφi(x)*φi(ξ)/(λi- [-s]) which is p 219 A, a very general result. This makes poles in the s place for the inversion contour and we end up with B and C, which merely reproduces earlier results. Stak is just doing verification that we have not screwed up somewhere. Example: (p 220) This example has our first 3D and our first 2D spatial examples of heat conduction problems! (3D) We have an infinite space with an empty spherical hole in it. Throughout the volume of this infinite space at t = 0 (outside the hole) we have the initial condition u = 0. But our boundary value on the spherical surface is u = 1. Notice that we have as part of our BC that u = 0 at r = ∞. So we go into spherical coordinates and write 7.73. We now apply our Laplace Transform method just studied in the previous section, so our system is p 220B after t replaced by s. This is the isotropic (n=0) Helmholtz equation with k2 = -s. We know that the Helmholtz solutions are e±ikr/r (see "Separation of Laplace and Helmholtz in spherical coordinates.doc"), and here we have k = ± i so our solutions are then exp(±ir)/ r, as shown in page 221 A. This quickly leads to 221 B (where we just applied our r = ∞ condition of no blowup at least). The claim is that the inverse gives C which evaluates to D (I have not checked these two results yet) where we have a simple result involving erf. As t → ∞, erf(0) → 0 and our limiting solution is 1/r. In heat statics, we would solve our equation without t to get u = 1/r and u = 1 as possible solutions, and I guess 1/r is the one appropriate for u(∞) = 0 as I assumed above. So think of our infinite 3D object as a heat conductor. We have a reservoir of u=1 applied to the sphere boundary, and another reservoir u = 0 at infinity, and u = 1/r is the steady state solution. This fact by itself is totally new to me. As you move away from the sphere or radius 1, u gradually drops off in this manner. It has to be some function of r, and it is 1/r. Probably you can think of the spherical hole as a point charge at the center, etc, which makes 1/r obvious. Stak notes that had we done the problem u = 1 on sphere and u = 1 for initial, then u=1 is the final solution for all time and is also the static solution. As a dynamic 3D problem, u(x,t) is as shown in D which I could plot. ( plot of 221 D) In the foreground as t→0, we get u = 1 at r=0 and u = 0 for r> 1, though the plot does not show this limiting case "all the way". In the back at larger t, we are getting u = 1/r. (2D) Now we just repeat the problem in 2D. We have then an infinite solid with a cylinder cut out of it, and u = 1 on the cylinder inside walls, and u = 0 in the bulk of the infinite solid. The setup is exactly the same except of course the Laplacian is 1/ρ ∂ρ(ρ ∂ρ) in this case, and 1/r2 ∂r(r2 ∂r) in the 3D case. The analog of system 220B is 222A, with only the difference just noted. This is now a 2D Helmholtz again with k2 = -s. This can be found in the Moon and Spencer cylindrical coordinates section page 15, where solutions are in general of the form Jn(ikρ), but here we are isotropic so n = 0, and ikρ = i(-i)ρ = ρ, so the nature of the solution is then J0(ρ). But for our problem we need finite at large ρ, so the solution of interest is then K0(ρ). The constant is selected as in p 222 B so that (ρ=1,s) = 1/s which is the boundary condition in p 222 A. Now the claim is that if you insert B into the inversion formula and deform the contour so it picks up the pole at s = 0 (which gives the "1") as well as the cut discontinuity on (-∞,0), you get result 7.76 where ν = in effect. This is an amazingly messy result for such a simple problem (compare to p 221D for the 3D case). Stak claims the large t limit is 1 now so that our ρ=1 boundary condition overwhelms the u = 0 bulk initialization. The reason is of course that the "heating surface" is a lot "larger" in this problem relative to the object being heated. His final comment is that even his fancy B.10 cannot handle the large t limit of 7.76. 7.4 Uniqueness and Continuous Dependence on the Data [ 223] We certainly expect that if u = 0 on the bottom and sides of the cylinder, that the unique solution is u = 0 everywhere inside. This is the content of Theorem 1 on page 223, which Stak proves two different ways. In the first proof, he uses Green #1 to process the heat equation 7.77 into the form 7.79 which says that the time slope of the R integral of |u|2 is negative. But this integral is 0 at t=0 from the initial condition, so it can either decrease or stay 0. But if it decreases, it goes negative, but the integral is obviously positive, so it must stay 0 and this means u must stay 0. He notes that the proof is not perfect, but good enough for me (ie, no doubt there are pathological counter cases). [224] The Maximum Principle (Theorem 2) says that if u is bounded in both initial and boundary conditions such that |u| ≤ M, then u is similarly bounded everywhere inside the cylinder. In other words, at no point in the interior can |u| exceed its max on the boundaries. I suspect that we may conclude that, barring the case that u = constant inside, the max of u must occur somewhere on the boundary. This is very similar to our conclusion regarding solutions of the equation (Stak does not mention this), and is also reminiscent of a fact about analytic functions in the z plane. Stak first proves the Lemma shown. This lemma says: suppose a function v has a negative spacetime "curvature" as in 7.80 (true, since first order in ∂t, not quite curvature, but I will call it that). Suppose then that v has some spacetime maximum inside our cylinder. If you look near that point (x,t), you find that all derivatives are 0, in particular ∂tu = 0. In addition, all spatial curvatures must be negative, so that then 2u < 0 (as in case). But this makes LHS of 7.80 be positive, which violates our starting assumption. Therefore, given this starting assumption, v cannot have a maximum inside the cylinder. Therefore its maximum must occur in the bottom or side of the cylinder. Now comes the Max Principle proof. Taking u which solves 7.77, construct v as in p 224C. Since 2acting on |r|2 is Σi∂i2 Σjxj2 = Σij 2 δij = Σi 2 = 2n, we get D and of course ε > 0 is assumed. So v meets the condition of our Lemma. And if we surround our assumed bounded region R by a sphere of radius a, then E is also true. Since v = u +εr2, and since the max of u on the boundary is M, we conclude that v ≤ M+εa2 on the boundary. The lemma says v must have its max on the boundary, so we then conclude from our lemma that v ≤ M+εa2 for all points inside the cylinder, which is the RHS of 225 A, while the LHS of 225 A comes from p 224 E. Let ε→0 and 225A says u ≤ M for all points inside the cylinder, and Theorem 2 (Max Principle) is proved. [225] The Minimum Principle is trivially also true therefore, so barring constant, the minimum value of |u| must also occur on the boundary. All as in Potential Theory. So the initial conditions, assuming nothing is infinite in them, bracket the solution everywhere inside, in that A ≤u ≤ B. Corollary 1 says that if you have three solutions of the heat equation (for three sets of BCs) such that u1≤u≤u2 at every point on the boundary (IC and BC), then u1≤u≤u2 also at every point inside. Uniqueness Theorem 3: If you had two (bounded = continuous) solutions of 7.81 with the same BC's, then the difference would be a solution with 0 BC's. But then Theorem 1 says that difference solution must be 0. An alternate proof is to appeal to the min and max principles which forces 0 ≤ Δu ≤ 0 so Δu = 0 everywhere. [226] Theorem 4 is the thing about continuous dependence of the solution on the boundary data. If you make small variations of ε in the boundary data, then you get a small variation in the solution, and this is guaranteed by the min and max principles, nothing to prove. At this point, right in the middle of things, Stak says that it is pretty hard to prove the existence of a solution, although we have proved uniqueness, and he is not even going to attempt a proof of existence. The reason is that it is hard to demonstrate a solution by construction in the general case, as we have done in past existence proofs. But of course in any specific problem, if we can build a solution, then it must exist. Theorem 5 is an extension of the Theorem 1 to unbounded regions: u = 0 on all boundaries means u = 0 inside. In his proof, Stak contrives two solutions v1 and v2 of the q=0 (homo) heat equation where v1 < 0 and v2 > 0 as shown. He thus has v1 ≤ u ≤ v2 on the boundary since u = 0 there. But Corollary 1 says then that v1 ≤ u ≤ v2 everywhere inside (R is the intersection of the unbounded R with a sphere a0) . As temp sphere radius a0 → ∞, this says 0 ≤ u ≤ 0 everywhere inside, so must have u=0 inside. It follows I think that we can also extend uniqueness Theorem 3 to unbounded regions. Yes, that is the corollary on top of page 227. Notice that almost every theorem here requires the restriction that u be a bounded (continuous) solution of the heat equation. This would seem to exclude situations where the q source for example is a delta function, but I suspect we are complete then too by taking a sequence limit. [227] Finally, Stak asks what we know about things like u∞ , the large time limit. If our sources q and boundary sidewall function h both have limits q∞ and h∞, we sort of expect that u∞ would be the solution of a statics problem 7.82, which is just Potential Theory. In other words, we look at some disk cross section of our cylinder at very large time, and ∂tu ~ 0 since things changing so extremely slowly, so we are in the limit. Theorem 6 (without proof) says this expectation is valid, providing that limits like q→q∞ are uniform over R (and the h one uniform on σ). 7.5 Miscellaneous Heat Conduction Equation Topics The Semigroup Connection. This is fascinating to me. In 7.83 you see our usual "free infinite rod" propagator moving u ahead by t units. As Stak notes, and as we all know, things get "spread out" as u moves forward in time. If you regard the propagator as an integral operator, we have ut = Gtu0 . The integral is over R, but time appears in the kernel. It seems pretty reasonable that Gt+t' = Gt Gt' since the zero of time is arbitrary. So the operators Gt almost form a group, but there is no inverse to Gt so it is a semigroup. There is of course an identity G0 where the propagator is just δ as we know. Stak claims you can prove directly Gt+t' = Gt Gt' with the rod propagator as in 7.86. But it is more generally true for any Causal Green's function as shown in 7.87. My Feynman Diagram for 7.78 is this where the double arrow means integration over R, a slice of the cylinder. You add up all the amplitudes passing through intermediate points. Comment on inverse: Imagine a problem where we set f(x) = some arbitrary value for u at t = 0, but we have u = 0 on the boundary (the entire boundary is a heat sucker). The eventual solution to this problem is u = 0 everywhere at t = ∞. It seems pretty obvious that given the t=∞ solution of u=0 that you could not propagate things backwards to recover f(x), because any f(x) would give the same result. The inverse is not unique. Something about the domain and range here. If u(t,x) meets suitable conditions (rather severe ones, see below), then there is a unique solution for the reverse propagator, and I guess in this case Gt would have an inverse and you then have a group instead of a semigroup. Comment: We think of heat flow in thermodynamics as an irreversible process, so it might make some sense that heat conduction problems don't have trivial time reversed solutions. [228] The Backward Problem. Intuitively we can sense that a final state in a heat conduction problem might come from an infinite number of different possible initial states at t = 0 as in our comment above. This jibes with the notion that the propagation operator Gt has no well-defined inverse! For the rod, the inverting problem would require inversion of 7.88 to find a given b at a later time. But it seems to me that we could just diagonalize 7.88 and solve it in the sense of Bk = where Ck = (1/2π) !Syntax Error, Idx exp(-x^2/4) e+ikx = (1/2π) 2exp(-k2) says Maple So then we have Bk = (1/)exp(-k2)Ak. So what is the big problem solving this? Ak = exp(k2) Bk A(x) = !Syntax Error, Idk Ak e-ikx = !Syntax Error, Idk exp(k2-ikx) Bk One issue is that Bk has to be very convergent for large k perhaps. So I guess I miss the point here. Why is this all OK if B(x) is a C∞ function? I found a 1961 paper by Miranker, but it is also unclear on this point. He points out that as you go backwards, u(x,t) has to be analytic in t so is likely not to match some abrupt BC's that you start with at t = 0, maybe that is the point. In any event, Stak claims that the backwards heat flow problem is ill-posed. Aside: Suppose b(x) = constant in 7.88, corresponding to our u = 0 everywhere t=∞ situation. Now 7.88 moves us forward Δt = 1 second. What is Bk in this case? Bk = FT(constant) = δ(k). A(x) = . So this does not seem to present a problem. As I get to page 230, Stak basically does what I did above and says that for given Bk, A(x) is likely to diverge and there is no L2 solution A(x) [ which is his a(x) ]. Then just after 7.91 Stak says that the integral (sum in his case since he does series for finite rod) will in fact converge if b(x) is D∞ . So I am missing something here! [ Well, looking at 7.88, if it were true, we could apply ∂xn to both sides and all this does is create powers under the integral so you get a set of terms of powers against expo decay, and any such integral I guess will converge so all derivatives exist so g(x) is D∞. This is just a rough idea that would have to be checked out. ] I found a very good PDF on Schwartz functions and it does talk about this issue. It says that if you are C∞ on one side of the FT, then you will have expo decay on the other side. So in our example, if b(x) is C∞, then Bk is going to decay very strongly and I guess enough to get the A(x) shown above to converge. OK, so this is just another gaping hole in the PL math knowledge. I guess I will fill it at some point, but not now. I think this set of Schwartz functions was mentioned in Stak earlier under a different name, functions of slow growth perhaps. Maybe Stak will have more to say on this. [ I read the PDF and took some notes, nothing earth shattering.] So you CAN go backwards in time just fine if your starting point in the future b(x) is C∞. However, the solution you get in this case a(x) "does not depend continuously on the data b(x)" and this means the reverse problem is "ill posed". Stak takes a simple example where b(x) = εsin(αx) and a(x) = εeα sin(αx) where we are allowed to pick any α we want. If |b(α) - b(α0)| is small and < ε, |a(α) - a(α0)| can be made arbitrarily large by increasing α. But what we want (for a well-posed problem) is that a small variation in b results in a small variation in the solution a. The discussion is now repeated in terms of the problem of a finite length rod where 7.90 shows the forward propagation from a(x) to b(x), similar to 7.88 for the infinite rod. In this case g is as in p 230 A. Thus, 7.90 is NOT of convolution form but I guess since it is in a separable form, you can still diagonalize the integral equation 7.90 using a Fourier Series. // Well, I studied this and added a section to my group theory document on convolution, go read sections 7 and 8 of that doc. Things diagonalize because we are in the "i basis" which diagonalizes our kernel, and this has nothing to do with group theory or convolution form. Here is Section 8 copied from that doc: _________________________________________________________________________________ Suppose we have this integral equation in heat conduction theory where t > t' : a(x,t) = ∫g(x,t;x',t')b(x',t')dx' g(x,t;x',t') = θ(t-t') Σi φi(x) φi(x')* e-λ(t-t') where the φi and λi come from the problem -2φi(x) = λi φi(x) with φi(x on σ) = 0. That is to say, these are the same φi we had in Section 7 above. We can write our integral equation this way, where we suppress the time labels on things, ax = Σx Gx,x'bx' and obviously our kernel is not diagonal in the x basis. But, we can show that G is diagonal in the "i" basis. Again, we just regard t and t' as parameters which don't change and we could put them in if we wanted, but we suppress them: Gij = ∫dx ∫dx' <i|x><x|G|x'><x'|j> = ∫dx ∫dx' φi(x)*g(x;x')φj(x') = ∫dx ∫dx' φi(x)*{ Σk φk(x)φk(x')* e-λ(t-t') }φj(x') = ∫dx' φi(x')* φj(x') e-λ(t-t') = δij e-λ(t-t') and indeed, in this "i basis", operator G is diagonal with diagonal elements as shown. [ These diagonal elements were of the form 1/λi in the potential version of the above.] Now back to the integral equation which was in the coordinate basis (and we continue to suppress time labels) <x|a> = ∫dx' <x|G|x'><x'|b> We slice off <x| and change the left sides to the i basis, <i|a> = ∫dx' <i|G|x'><x'|b> Then we insert completeness Σj |j><j| = 1 on the RHS after the G <i|a> = ∫dx' <i|G Σj |j><j|x'><x'|b> = Σj Gij∫dx'<j|x'><x'|b> = Σj Gij <j|b> and no surprise we end up with ai(t) = Σj Gijbj = Giibi = e-λ(t-t')bi(t') So here we have diagonalized an integral equation in which the kernel is a heat conduction theory Green's Function. We did not use the "group theory method" since the Green's Function is not of convolution form in spatial coordinates. We simply used the "find a basis in which the integral operator G is diagonal" method. Example: In the discussion of Stak page 229-230, we have our finite length l rod with u = 0 on both ends, and we have a set of φn(x) = sin(nπx/l) so that g(x,t;x',t') = θ(t-t') Σn φn(x) φn(x')* e-λ(t-t') = θ(t-t') (2/l) Σn sin(nπx/l) sin(nπx'/l) e-λ(t-t') where in this problem λn = (nπ/l)2 . Our integral equation b(x,t) = ∫g(x,t;x',t')a(x',t')dx' is stated as 7.90 (notice we have reversed b and a relative to discussion above) where he uses t' = 0 and t = 1 and suppresses the time indices in a and b. We know from above that our diagonalized equation will be this: bn(t) = e-λ(t-t')an(t') or bn(1) = e-λ an(0) or bn = exp[-(nπ/l)2] an or an = exp[ (nπ/l)2] bn which we see agrees with page 230 B and C. The projections are given of course by an = ∫dx φn(x) a(x) = !Syntax Error, Idx sin(nπx/l) a(x) and these are just Fourier Sine Series projections, so we are talking this transform. If we attempt to solve for a(x) given b(x) we get a(x) = Σn an sin(nπx/l) = Σn=1 exp[ (nπ/l)2] bn sin(nπx/l) // page 230 D and Stak comments that it is "unlikely" that this Σn converges. He knows and we know that the only chance this has of converging is if b(x) is infinitely differentiable (and maybe it won't converge even then), see page 230. But assuming Stak's condition E is met by the projections of b(x), there will be a unique solution a(x). _________________________________________________________________________________ Stak wraps up with the data continuity ill-posed issue for this finite length rod problem, all is good. Verify that my red factors should be there. I insert MY a(x) into (7.90) to get b(x) = !Syntax Error, I g(x,1 | ξ,0) a(ξ)dξ = !Syntax Error, I dξ { (2/l) Σn=1 exp(-n2π2/l2) sin(nπx/l)sin(nπξ/l) } { Σm exp[ (mπ/l)2] bn sin(mπξ/l) = (2/l) Σnm exp(-n2π2/l2) sin(nπx/l) exp[ (mπ/l)2] bn * !Syntax Error, I dξ sin(nπξ/l) sin(mπξ/l) Now I look up this required integral I = !Syntax Error, I dξ sin(nπξ/l) sin(mπξ/l) Change variables to θ = πξ/l so that upper endpoint in x becomes π and dθ = (π/l)dξ so I = (l/π) !Syntax Error, Idθ sin(nθ) sin(mθ) I look this up in my transforms.doc Fourier Sine Series Transform section to find that I = (l/π) (π/2) δnm (1-δn0) = (l/2) δnm (1-δn0) Insert this into the above to get b(x) =(2/l) Σnm exp(-n2π2/l2) sin(nπx/l) exp[ (mπ/l)2] bn (l/2) δnm (1-δn0) = Σnm exp(-n2π2/l2) sin(nπx/l) exp[ (mπ/l)2] bn δnm (1-δn0) Both sums begin at 1, so the last factor does nothing, replace it with 1 = Σnm exp(-n2π2/l2) sin(nπx/l) exp[ (mπ/l)2] bn δnm = Σn exp(-n2π2/l2) sin(nπx/l) exp[ (nπ/l)2] bn => b(x) = Σn bn sin(nπx/l) So basically I end up with a(x) = Σn an sin(nπx/l) b(x) = Σn bn sin(nπx/l) Stak is not clear (see 2nd line of text on page 230) exactly how he is defining the an and bn , though presumably both the same way. If he omits the red factor in the first, then I have shown it is missing as well from the second. So things are correct either way. In every one of his equations, things float. So let this one go! The Diffusion Interpretation of the Heat Conduction Equation (230) Recall that electrostatics, static heat flow, no curl fluid flow, and static diffusion are all potential theory worlds. The claim here is that non-static diffusion is like non-static heat flow where u is particle density of the solute in the solvent. Zero temperature u = 0 corresponds to "no particle density", so I guess this is a place where particles are sucked out of your system, just the way u = 0 sucks out heat. So if u = 0 on a boundary, that is an absorbing boundary, it absorbs solute particles. And if ∂nu = 0, no particles cross (since no diffusive force to do so) so an insulating boundary becomes a particle reflecting boundary. I have NEVER done any of this. For some reason, Stak shows how the diffusion equation gets modified if there is a "steady drift" of velocity vx = 2γ in the x direction say (so l = x) , and this is 7.92. This is an "example" of the Fokker-Planck equation says Stak, yet another massive hole in the PL math world, I see it is all tied in with "stochastic processes". But in our simple example here, Stak shows how you can remove the drift and get back to the heat equation by changing from variable u to v. Stak then does a toy example showing how you might start in 1D with a δ source of solute and watch the drift current spread it out. We just solve our normal heat equation for v using our rod propagator C, then convert that back to u. I think I can "see" in the solution 7.95 that our initial puff of solute drifts to the right and spreads out at the same time. The first pair of expos gives a combined 1 when γ(x-x0)-γ2t = 0 which means (x-x0) = γt, so x(t) = x0 + γt so the peak here drifts to the right at speed γ. The other factors describe the spreading out always centered at x0. Asymptotic Formula for The EV's of -2 (231) Remember that even in the heat problem, the EF's φi are just those of regular old potential theory and are thus determined by the boundary σ contained in space R. [232] Exact 7.97 is processed easily into exact 7.98 using <φk|φk> = 1. We consider only a 2D space R like the surface of a heat conducting plate. But since this is the generic problem, we can think of the associated 2D membrane problem. Question: Why is g ≈ C for small t ? // long digression This is a bit of a digression but I think a good one and I can dispel some misunderstandings. (1) In potential theory, I have written up my "Dirichlet method for finding a Green's function". The idea is that you write g(x|ξ) = 1/R + ∫σ dS σ/R = point charge potential + potential due to the induced charge on σ where as usual σ has two unrelated meanings : name of boundary, and surface charge. We go off and focus on the Dirichlet problem to find the potential just of the induced charge, where the Dirichlet prescribed potential is -1/R on the boundary. Then after we solve that, we add 1/R to get g. Notice that the induced term (usually called v) satisfies homo inside R, whereas the 1/R term generates the δ(x-ξ) term. So we can think of the term ∫σ dS σ/R as the homo solution we add to our fundy solution so that our g meets the boundary condition g = 0 on σ. If we think of g0 = 1/R as the "free space propagator" or Green's function or fundy solution, then we write the above as g = 1 g0 + ∫σ dS σ g0 (**) g(x|x0) = g0(x|x0) + ∫σ dSξ [∂nξg(x|x0)] g0(ξ|x0) Now my first clarification is that this point + induced segregation does not align with the Green's Theorem statement which is this (that is, we don't equate first term with first term and surface integral terms with surface integral terms) u(x) = ∫V dξ g(x|ξ)q(ξ) – ∫S dSξ [u(ξ) ∂nξ g(x|ξ) – g(x|ξ)∂nξu(ξ) ] (*) For example, if we set q(ξ) = δ(ξ-x0) as our Green's point charge, and we have g = 0 on the boundary, then u is g, and the above says g(x|x0) = g(x|x0) + 0 + 0. So again, don't align the above two equations with each other, they are different (but both true). And the g in the first term of (*) is the full g, not g0. Now back to that Dirichlet method. Our equation g(x|ξ) = 1/R + ∫σ dS σ/R basically says that the only "thing" that can create potential is charge. And in this Green's problem we have a point charge, and we have induced charge, and there is nothing else!! Then, in this framework, we can say that if we are close to the point charge, the 1/R term is very large, the surfaces are "far away" so their R factors are large, and we can approximately say that g(x|ξ) ≈ g0(x|ξ) when x is close to ξ. (2) Now, how does this play out in heat conduction theory? If we write (∂t-2) g(x,t|x0,t0) = δn(x-x0)δ(t-t0) g = 0 on boundary σ we argue that the analog of a point charge is a spacetime point heat source. Whereas the point charge is a point in space but endures for all time, the spacetime point heat source is localized in both space and time. The argument is then that the only thing that can generate g or u(x,t) in heat theory are such heat sources, and they can either be distributed as q(x,t) within R, or they can appear on a boundary, and there is nothing else. The Poisson equation would be (∂t-2) u(x,t) = q(x,t) for example and you see the distributed sources explicitly here, and the "only other thing" that can affect u are the surface heat sources ∂nu. I would then like to write the analog to (**) above which might be this: ( I set t0 = 0 now) g(x,t|x0,0) = 1 g0(x,t|x0,0) + ∫σ dSξ [∂nξg(x,t|ξ,0)] g0(ξ,t|x0,0) where [∂nξg(x,t|ξ,0)] represents spacetime heat sources or sinks on the boundary. For the g problem we are saying that g = 0 on the boundary, which means the boundary is "frozen" at zero temperature, which means heat wants to flow out through points in the boundary and these points are then negative spacetime heat sources. We then say "this is all there is" in terms of things which can "create g". It is like saying we have the point charge and the induced charge on the boundary. Also, since the first term g0 generates our spacetime δ, the second term is a homo solution of the heat equation, just as we say in potential theory. This second term is the homo solution that adjusts g0 such that g = 0 on the boundary. Now Stak has not talked at all about the notion of the velocity of heat or temperature propagation, but it is implicit in the form of g0 = C where the effective radius of the gaussian ball is ρ = (2) and if we think of this as ρ = vt, then v = 2/ so the spreading velocity is initially infinite, then slows down. Now we know that both g and g0 → δ(x-x0) as t→0+. If we now consider very small time t, then in the above equation there is a small ball (in which the function has significant non-zero value) around x0 from the 1 g0(x,t|x0,0) term, and there are small balls around points on the boundary where we have those driving heat sinks. But for some x0 out in the middle of R, at small t these balls don't overlap, and that is why we can say for small time that g(x,t|x0,0) ≈ 1 g0(x,t|x0,0) in the neighborhood of the point x0 ! If we go to large time, then the balls do overlap, and we cannot make this claim. So it really is a matter of the velocity of propagation of heat flow. Stak only hints at this notion at my arrow marked B on page 232. [233] Everything is fine and we want to replace g by C in 7.98, that is the upshot of all of the above, for small t. The space integral of C (with x = x0 as in 7.98) is just 1/(4πt) times the volume of R which we call AR since we are thinking of a membrane type R. Thus right off the bat this EV sum for small t depends only on the area of our membrane which seems pretty amazing. How do we connect (7.100) with Appendix B? He gives no reference. Appendix B is about t being large, not small and so Appendix B would seem to have no relevance whatsoever! Appendix B does discuss moving the C inversion contour to the left, that is really the only connection! Our 7.100 interest is this J(t) = !Syntax Error, Idλ N(λ)e-tλ our interest is t small J(t) ≈ AR/(4πt2) Let's write this in normal transform notation as in Schaum p 161 J(s) = !Syntax Error, Idt N(t)e-st our interest is s small J(s) ≈ AR/(4πs2) The inverse of this transform is N(t) = (1/2πi) !Syntax Error, Ids est J(s) ≈ (AR/4π) (1/2πi) !Syntax Error, Ids est/s2 N(t) = (AR/4π) (1/2πi) !Syntax Error, Idx' ex't/x'2 The integral here is this: (1/2πi) !Syntax Error, Idx' ex't/(x'-0)2 and if we drag the contour to the left, the double pole gives this contribution (1/2πi) dx' ex't/(x'-0)2 f(x') = ex't where I put things into symbols I use in TK page 15. For m = 2 this integral is 1/(1!) * [D(1)f(x)]|x=0 = [ t ex't] |x=0 = t. So the result then is N(t) = (AR/4π) t or N(λ) = (AR/4π) λ // agrees with 7.101 We can close the contour C to the left on half great circle and get 0 there since s/s2 → 0, so this double pole contribution is all there is. What does 7.101 tell us? The number of eigenvalues below λ is proportional to λ for large λ. This tells me that for large λ, they are linearly spaced. Double λ, double N if λ already large. But Stak's interest is that N(λ) for large λ depends only on the membrane area AR. This whole business really has nothing to do with heat conduction, but he used the heat conduction stuff to derive the result. The result concerns only eigenvalues for a problem. The rest of this section goes as follows: As an improvement in our N(λ) calculation, let's replace the actual boundary σ by a planar boundary that is tangent to the actual boundary at the point of the boundary's closest approach to x0. Perhaps this would be the first term in an expansion for the boundary in some sense. Then we can improve our estimate of the propagator using a negative image source, and this gives us the propagator of p 233D. We then just turn the crank of the previous section and our improved result is 7.102 where the lead term depends on membrane area AR, while the correction term depends on the membrane perimeter LR. Regrettably, Stak gives us no "handle" with which to learn more about this result and Weyl's role in it. But I find it under "Weyl's formula for the eigenvalue distribution".... and there are active papers going on even now in this area. Sometimes N(λ) is called "the counting function" and papers try to put bounds on this thing. I found and saved a long paper giving the history of this Weyl's Law, his work was 1912. Note that boundary volume is often called Ω and boundary area ∂Ω . OK, good for Stak to throw this in. Note added 5.22.11: Suppose the membrane were just a square of edge l. In this case we know λnm = (π/l)2(nx2 + ny2) which we could write as λn = (π/l)2 n2 . How many eigenvalues lie inside a disk of radius λn? The answer is the number of dots in a circular 2D integer grid of radius n = (l/π) dots, which is the same as the area of such a disk, which is π [(l/π) ]2 = π (l/π)2λn = N(λ) = (1/π) AR λn, but we must only count the dots in the first quadrant to avoid overcounting the EV's, so get N(λ) = AR λ/(4π) in agreement with Weyl. I presume the correction term relates to the ambiguity of the dot count near the perimeter of the disk, but I have not looked into that detail. Certainly this ambiguity would be proportional to the perimeter of the disk which in turn is proportional to l which in turn is proportional to the perimeter of the square. A 3D Composite Medium Heat Conduction Problem This problem concerns the interface between the inside of a sphere of radius 1 (u1, material 1) and the outside (u2, material 2). The IC is that inside we have u = 1, and outside we are frozen with u = 0. The heat equations are A and B on page 235, where we write 2 as apropos for a symmetric problem like this one. We need first to wander back to Volume 1 page 327 and review "heat conduction". We assume a medium has a constant thermal conductivity k, and that means at a point on the boundary between two uniform media we are going to have k1∂nu1 = k2∂nu2 at any instant in time. This is just two ways to express the heat flowing through a patch of boundary, and it has to be the same from whichever side you measure it. It does not say the flow of heat through the boundary is 0. And of course we expect u1 = u2 at a point on the boundary, so this explains the BC's C and D on page 235. Next we note from Stak p I.328 that the constant appearing in the heat equation is a = k/c where c is the heat capacity and a is called the thermal diffusivity, since it occupies a diffusion constant position in the heat equation. We then have a1 = k1/c1 and a2 = k2/c2. Equations A and B are Laplace transformed into E and F, and the RHS's here are the -F(0) you see on page 162 Schaum 32.7 (but they are +F(0) when put on the right sides). The solution to either of these equations is a linear combination of atoms exp(±r/a)/r, see notes p 220 above where we recognize this as the solution of isotropic Helmholtz with negative k2. for the outer region we have to take the decaying expo and this explains p 235H. For the inside the only choice is the one shown to get a u1 which does not blow up at the origin (due to that external 1/r factor). I think it is easy to convert the two BC's C and D to the Laplace space (apply operator to both sides) and then you can use these two conditions at r = 1 to determine constants A and B as shown in I and J (I have not checked this). [236] He now sets r = 0 and I am happy with p 236 A, but I have not done the math to show this can be written as in B and C with series f(s) as shown, don't think there is anything unusual here. He claims that when the dust settles you find that for small s, your leading solution term is as in H. In the inversion formula, small s makes the main contribution to the large t result since est is sitting in the integrand. So if we approximate our solution with just the leading small s term and invert, we get a result for the large t behavior of the result. Since we set r = 0, we are just talking about the temperature at the origin which is u1(r=0,t). He has shown that the heat in the inner radius 1 ball will dissipate in such a way the temperature at ball center drops off to 0 at the rate t-3/2 as shown in p 236 I. So this is the main result of this calculation. You start off with a "hot ball" radius 1 that is embedded in "ice", and you sit at the origin and you watch the heat flow out over time. We really have the exact answer for all time t if we can do the inversion on page 235 G with B given by I. As you can see, this is a fairly complicated function of s, and so inversion might be difficult. [237] Stak's final zinger on this problem is shown at the top and I like it. Suppose the two media are the same, and we just started things off with our hot core ball (think of the earth by the way). The total heat in the ball is just Q = (4/3)πR3c = 4πc/3. We could think of a different problem where our initial condition was a big bang with Q/c as a spacetime point source at the origin, which is to say, the heat equation was driven by (4π/3) δ(r)δ(t) [ See p 328 A.11 for this 1/c factor] . It is true that as this bang expands, it will never match the uniform u=1 ball of the original problem, but the total heat is the same. Somehow you must be able to argue that if you wait a long time, the way that heat was arranged at the origin does not really matter in determining the temperature at the origin. This is because 3D space into which the heat has flowed is huge compared to the little initial ball, and it does not know whether the heat came from a ball or a point. Based on this argument, you argue that u is given in this case of "same media" by the free space propagator C in n=3 dimensions. The expo becomes 1 for large time. So looking at 5.140 for n=3, you predict for large t that u will have the form shown in p 237 A, and this does exactly agree with the same-media limit of the original problem with the two different media. Question: Why not just directly do an inverse Laplace Transform on page 236 A ? The reason is that the object B is a very complicated function of s, as shown page 235 I. Aside: If you look just at p 236 H you see that our leading term is ~ s1/2. This by itself is not a "legal" F(s) as I discuss in my Laplace transform notes, but it is legal in the sense of finding large t behavior from this component of the full (and legal) F(s) which we are here approximating. This is discussed in the last part of Appendix B, and the result is consistent with Schaum p 164 32.28 continued to the value n = -1/2. Errata page 237? I agree with 7.104 from 7.103, so write for the denominator, (6π1/2a3/2t3/2 ) = π3/2 π-3/2 43/2 4-3/2 (6π1/2a3/2t3/2 ) = 6π1/2π-3/24-3/2 (43/2π3/2 a3/2t3/2 ) = 6π1/2π-3/24-3/2 (4π at )3/2 = 6 π-12-3 (4π at )3/2 = (6/8π)(4π at )3/2 = (3/4π)(4π at )3/2 Therefore 7.104 becomes 1/denom = (4π/3) / (4π at )3/2 BUT, (5.140) gives the free heat propagator as just the denominator above, so the "point source" must be 4π/3, so I think my errata call is correct. The Stefan Problem (237-238) I found a whole Google book just on this class of problems, where a phase change boundary moves. In 1889 Mr. Stefan (a Slovene) was wondering about the way ice forms on the top of a lake, and the way in which ground freezes. With a uniform cross section, you can think of this as a 1D heat conduction problem, but Stak talks about it in this 3D sense. You have x ≥ 0 filled with ice at u=0 and at t = 0 you apply a metal plate at u = U against this ice (plate in the x = 0 plane). The ice starts to melt and a water layer forms and gets wider with time. We don't worry about density change here, we assume water exactly fills the space of the melted ice. We assume a planar phase interface which moves to the right. The problem then is to find an expression x = ξ(t) which describes the interface location in time, and also to find an expression u(x,t) for temperature in the (ever-thickening) water layer. Temperature is u = 0 to the right of this layer, and we don't care about x < 0 and can assume it is U there if we like. Stak proceeds to set up and then solve this problem. The key fact is this: in a thin layer Δx at the interface, the heat required to melt the layer has to be supplied by the water layer on the left. A feature of this problem is that the boundary moves with time, and one of the BC's reflects this fact. In 7.105 we are saying just this: the heat provided by the water is k ∂nu per unit time so k ∂nu Δt is the heat going into a unit area of the layer. This heat melts the ice and so is equal to νρΔξ where ν is the latent heat of melting per unit mass so this thing is Q/M * M/L3 * L = Q/L2. This results in the problem setup as shown in 7.106 where are region of interest in x is (0,ξ) and ξ is moving (slowly). Stak solves this solution be writing a linear combination of atomic forms which we recognize from earlier work in this chapter, with constants A and B. He finds these constants using the BC's and puts together a complete solution to the problem. The interface position is given by 7.107 where α is a constant you have to get by solving a transcendental equation D (arising from the BCs). Then the temperature in the water region is given by E. Very elegant. Stefan is also known for empirically determining the T4 power law for black body radiation, and his student was Boltzmann who provided the theory to derive the formula, the Stefan-Boltzmann Law! And so ends the heat conduction portion of Chapter 7, apart from the exercises which follow. Think of this as "the heat exam" at the end of this mini course. I will do the exam, then move into the wave equation. The heat part of this chapter was p 194-243, or about 50 pages. The wave part goes to p 331 and so is 90 pages, yeouch! Exercises Exercise 7.4 Solve the ring problem by doing Fourier in the x variable. In the text page 211 this problem was solved by the method of (an infinite number of) images where the result was 7.48. So the first step is to write down Fourier for our ring. Stak wants us to use the complex series form so I quote from my transforms doc ( where n is summed over ALL integers) f(θ) = Σn fn e-inθ // expansion fn = (1/2π) !Syntax Error, Idθ f(θ) e+inθ // projection (1/2π) !Syntax Error, I dθ einθ e-in'θ = δnn' // orthogonality (1/2π) Σn e-inθ' einθ = δ(θ-θ') // completeness Let θ = 2πx so ring runs (-π,π) for x in (-1/2,1/2) and that puts us into the above space. We can see that dθ = 2π dx ∂θ = (1/2π) ∂x ∂θ2 = (1/2π)2 ∂x2 so our homo heat equation is this [∂t - (2π)2∂θ2] u(θ,t) = 0 Expand u(θ,t) as u(θ,t) = Σn un(t)e-inθ ∂θ2 u(θ,t) = Σn un(t)(-in)2e-inθ = – Σn un(t)n2e-inθ Then our homo heat reads Σn ∂tun(t)e-inθ + (2π)2 Σn un(t)n2e-inθ = 0 and since e-inθ form a complete set on our θ interval, we conclude that ∂tun(t) + (2πn)2 un(t) = 0 which is a simple first order ODE the solution of which is un(t) = An exp[-(2πn)2t] where An are constants to be determined. At this point we have this expansion for our solution u(θ,t) = Σn un(t)e-inθ = Σn An exp[-(2πn)2t]e-inθ Now at t = 0 we have some initial temperature prescription u(θ,0) = f(θ) so we write f(θ) = Σn An e-inθ and looking at our transform above, we conclude that An = (1/2π) !Syntax Error, Idθ f(θ) e+inθ So here then is our complete solution to this problem expressed in θ: u(θ,t) = Σn An exp[-(2πn)2t]e-inθ where An = (1/2π) !Syntax Error, Idθ f(θ) e+inθ Now we can re-express this solution in terms of x where θ = 2πx. Let (prime not derivative) u'(x,t) = u(θ,t) f'(x) = f(θ) dθ = 2π dx u(θ,t) = Σn An exp[-(2πn)2t]e-inθ where An = (1/2π) !Syntax Error, Idθ f(θ) e+inθ u'(x,t) = Σn An exp[-(2πn)2t]e-in2πx where An = (1/2π) !Syntax Error, I 2π dx f'(x) e+in2πx and now remove the primes, cancel the 2π, and rename An to be fn , u(x,t) = Σn=-∞∞ fn exp[-(2πn)2t]e-in2πx where fn = !Syntax Error, Idx f(x) e+in2πx and this agrees with Stak except with i → -i. But of course symmetry removes any imaginary part and we can write both our solutions this way: u(θ,t) = Σn=0∞ An εnexp[-(2πn)2t]cos(nθ) where An = (1/π) !Syntax Error, Idθ f(θ) cos(nθ) u(x,t) = Σn=0∞ fn εnexp[-(2πn)2t]cos(2πnx) where fn = 2!Syntax Error, Idx f(x) cos(2πnx) The last we can write as u(x,t) = f0 + 2 Σn=1∞ fn exp[-(2πn)2t]cos(2πnx) where fn = 2!Syntax Error, Idx f(x) cos(2πnx) Now suppose our initial condition were f(x) = δ(x). Then we find that fn = 2*(1/2) = 1 (with fn written as shown here, we pick up half the delta!). Then we have u(x,t) = 1 + 2 Σn=1∞ exp[-(2πn)2t]cos(2πnx) which exactly agrees with our image result 7.48. Remember that our 7.8a form of this problem is the homo equation along with u(x,t) → δ(x) if the 7.8 problem has δ(x)δ(t) on the RHS, which it did for the problem which gave solution 7.48. All done! Exercise 7.5 Redo the Stefan problem with ice at u = V < 0. The two homo heat equations are quite clear, the exact same form in the water and ice regions, in each region we can have heat flow of some sort. The picture is this (imagine white phase change layer is very thin) ( other picture for use later) and there are several new features compared with the text Stefan problem. I have made crude drawings of the temperature in the two regions at some instant of time where the interface layer is at ξ. In the water region we know that u1(0,t) = U. At the layer we know that u1(ξ,t) = u2(ξ,t) = 0 because in this layer the water has come from ice just melted and must be at u=0. To the right there is a region of the ice that is not yet melted, but which is above temperature V < 0. At t = 0 we has u2(x,0) = V because we started out with a solid half space of ice at u=V. So this explains every equation Stak lists off on page 239 except the last one which we shall now consider. At the right edge of our layer, the ice has already been brought up to u = 0 (as shown). Heat flows into our layer from the left, some of that heat is used to convert ice to water, and some of that heat continues on into the ice to raise its temperature. This second heat flow was absent in the text version of the problem. We have heat going to the right at the left edge of the layer = - k1∂xu // notice ∂xu < 0 heat going to the right at the right edge of the layer = -k2∂xu // and again, ∂xu < 0 The difference in these heat flows is what melts the ice so we have - k1∂xu - (-k2∂xu) = νρξ'(t) where the RHS is as in the text problem. In more detail - k1∂xu(ξ-,t) + k2∂xu(ξ+,t) = νρξ'(t) where ∂xu(ξ-,t) ≡ [∂xu(x,t)]x=ξ- etc and we allow for a possible slope discontinuity at the layer in temperature. So now we have all the equations Stak shows to define this problem. Let's blindly follow his text method and see if it works here. If we can find a solution that works, that is it! u1(x,t) = A1+B1erf(x/2) u2(x,t) = A2+B2erf(x/2) If t = 0+, the exponents z1,2 are very large and erf = 1 according to our plot, so we get u1(x,0) = A1+B1 = U ?? u2(x,0) = A2+B2 = #2 At time t = 0, there IS no water. We could imagine a super thin layer of water at x = 0 and it would then have to be at U at time t = 0+. If this is valid, we have A1+B1 = U. The second equation is more reasonable since the ice exists at t = 0+ and is all at V. If x = 0, the exponents z1,2 vanish and erf = 0 according to our plot, so we get u1(0,t) = A1 = #1 u2(0,t) = ?? The first equation is OK since water touches the wall x = 0 at all time. But the second equation does not make sense since the ice does not touch x = 0. If x = ξ we get u1(ξ,t) = A1+B1erf(ξ/2) = 0 (*) BC#3 u2(ξ,t) = A2+B2erf(ξ/2) = 0 BC#4 Perhaps this should have been our first consideration. These tells us that ξ/2 = C1 ξ/2 = C2 but these are stating the same fact, so let's write that fact as in the text problem and say ξ(t) = 2α1 => ξ/(2) = α1 Since A1 = U for sure, we have from (*) that U + B1 erf(α1) = 0 We have so far mentioned all the BC's except BC #5. We have ξ(t) = 2α1 = 2α1 t1/2 ξ'(t) = 2α1 = α1 t-1/2 = α1 νρ ξ'(t) = νρ α1 = α1νρ In passing, we have now verified the RHS of equation p 238 C right, which I did not do while reading the text. Meanwhile, we have [ recall from earlier comment on erf that ∂z erf(z) = (2/) e-z ] u1(x,t) = A1+B1erf[x/(2)] // here z ≡ x/(2) ∂xu1(x,t) = B1 ∂z erf[x/(2)](∂xz) = (2/)B1exp[-x2/(2)2] /(2) = B1exp[-x2/(2)2] / ∂xu1(ξ,t) = B1exp[-ξ2/(2)2] / = B1exp[-α12] / - k1∂xu(ξ-,t) = - B1 k1 exp[-α12] / Again in passing, we have now verified the LHS of equation p 238 C right. Next we have to worry about u2(x,t) = A2+B2erf[x/(2)] // here z ≡ x/(2) ∂xu2(x,t) = B2 ∂z erf[x/(2)](∂xz) = (2/)B2exp[-x2/(2)2] /(2) = B1exp[-x2/(2)2] / ∂xu2(ξ,t) = B2exp[-ξ2/(2)2] / = B2exp[-(a1/a2)ξ2/(2)2] / = B2exp[-(a1/a2)α12] / +k2 ∂xu2(ξ+,t) = + B2 k2 exp[-(a1/a2)α12] / Then our BC #5 reads like so: - k1∂xu(ξ-,t) + k2∂xu(ξ+,t) = νρξ'(t) - B1 k1 exp[-α12] / + B2 k2 exp[-(a1/a2)α12] / = α1νρ Notice that each of these terms has 1/ so we cancel that out to get - B1 k1 exp[-α12] / + B2 k2 exp[-(a1/a2)α12] / = α1νρ and this looks like something we might be able to solve for a value of α1. From one of our results above we have U + B1 erf(α1) = 0 so we have something to install for B1 . But what about B2 ? Well, we do have these conditions from above: A2+B2 = V BC #2 A2+B2erf[ξ/(2)] = 0 BC #4 We can write ξ/(2) = ξ/(2) = α1 Subtracting the two equations above we get B2 { 1 - erf(α1) } = V // and A2 = V - B2 = - B2 erf(α1) B2 = V / { 1 - erf(α1) } // A2 = - V erf(α1) / { 1 - erf(α1) } So our equation to determine α1 is this - B1 k1 exp[-α12] / + B2 k2 exp[-(a1/a2)α12] / = α1νρ B1 = -U/ erf(α1) B2 = V / { 1 - erf(α1) } Let's multiply the first through by to get - B1 k1 exp[-α12] + B2 k2 exp[-(a1/a2)α12] = α1νρ a1 Now insert B1 and B2 to get (**) : U k1 exp[-α12]/erf(α1) + V { 1 - erf(α1) }-1 k2 exp[-(a1/a2)α12] = α1νρ a1 Now in passing, if V = 0 we should recover our text Stefan problem and this says U k1 exp[-α12]/erf(α1) = α1νρ a1 exp[-α12]/erf(α1) = α1νρ a1/ (U k1) and I have now verified p 238 D (for the first time). So despite its ugliness, I assume we could solve (**) for α1 and that the solution is unique. I appeal on physical grounds to its uniqueness, hoping as well that the assume A + B form of the solution I took was general enough. If we start at V = 0 and gradually make it go negative, hopefully our solution α1 moves slowly as well and stays in existence. We can now state the total solution to our problem : first, solve (**) go get the value of α1 as a function of k1,k2,a1,a2,ν,ρ,U,V The location of the melting layer is at ξ(t) = 2α1 which is my pencil plot page 237 The solutions in both regions are these u1(x,t) = A1+B1erf(x/2) u2(x,t) = A2+B2erf(x/2) A1 = U A2 = - V erf(α1) / { 1 - erf(α1) } B1 = -U/ erf(α1) B2 = V / { 1 - erf(α1) } So write them out u1(x,t) = U[ 1 -erf(x/2)/ erf(α1) ] // agrees with p 238 E u2(x,t) = V[-erf(α1) + erf(x/2) ] / { 1 - erf(α1) } We confirm that u2(x,0) = V which is BC#2. We confirm u1(0,t) = U which is BC #1. We confirm u1(ξ,t) = 0 which is BC #3. We confirm u2(ξ,t) = 0 which is BC #4. What happens at large t? We get u1 → U for all x, and u2 → -Verf(α1) / { 1 - erf(α1) } which I don't quite understand. The boundary is way out at ξ(t) = 2α1 and I guess the heat has penetrated the ice and made it be this fraction of -V, so it is still ice. Well, this limit is ill-defined so let's not worry about it. The problem as given in this problem could represent the melting of a pond surface thick sheet of ice by earth-warmed water below it, though you would have to fiddle things a lot, but I think the general idea is valid Exercise 7.6 Freezing of water in a cylindrical tank. I want to find a way to map the solution to the previous exercise onto the solution to this one. Here is a suggestive pair of pictures: Exercise 7.5 Exercise 7.6 where I have made up my own labels for Exercise 7.6. The heavy vertical bar in each picture shows the location of a constant temperature "plate". Here is my proposed mapping between the problems: V = Z U = W but in 7.5 we have U > 0, while in 7.6 we would have U = W < 0 u1 = u3 a1 = a4 k1 = k4 water u2 = u4 a2 = a3 k2 = k3 ice x = h-y If this is correct, then I can for example write the u3 solution this way: u3(y,t) = u1(h-x,t) and recall from 7.5 that u1(h,t) = U[ 1 -erf(x/2)/ erf(α1) ] so we get ( α1 keeps its name and is solution of same equation) u3(y,t) = u1(h-x,t) = U[ 1 -erf((h-y)/2)/ erf(α1) ] = W[ 1 -erf((h-y)/2)/ erf(α1) ] which at least does the right thing when y = h, giving u3 = W. Now our solution 7.5 said that the interface moves at x(t) = 2α1 which I would translate into h-y(t) = 2α1 => y(t) = h - 2α1 So at t = 0, the interface is located at y(0) = h which is correct, we have all water. As t increases, the interface moves toward y = 0. The thing is filled with solid ice at this time y(t) = 0 or h = 2α1 h2 = 4 α12 a4 t t = h2/(4 α12 a4) So our only problem is to compute α1 . We start with our transcendental equation for 7.5 U k1 exp[-α12]/erf(α1) + V { 1 - erf(α1) }-1 k2 exp[-(a1/a2)α12] = α1νρ a1 and we make our replacements as listed above to get W k4 exp[-α12]/erf(α1) + Z { 1 - erf(α1) }-1 k3 exp[-(a4/a3)α12] = α1νρ a3 where now everything is expressed in terms of our 7.6 problem. So the answer is: solve this equation for α1 then t = h2/(4 α12 a4) In his statement of the problem, Stak says to use Z = U and W = V which certainly confuses things, but we can put these into the above equation to get this V k4 exp[-α12]/erf(α1) + U { 1 - erf(α1) }-1 k3 exp[-(a4/a3)α12] = α1νρ a3 We can also write this in terms of our problem 7.5 symbols for ai and ki to get V k4 exp[-α12]/erf(α1) + U { 1 - erf(α1) }-1 k3 exp[-(a1/a2)α12] = α1νρ a1 In this notation, we find that the 7.6 equation is the same as the 7.5 one with V ↔ U. And then our time is t = h2/(4 α12 a1). I suppose if I were more serious, I would just set this problem up from scratch and do it and not relate it to problem 7.5, but it seems a shame to do all that math again. As for his last question, it is unclear whether we are adding a second "plate", or whether we are saying we have only a plate at the base. I think he is adding a second plate, and yes, that is a new problem since now we have three starting temperatures to worry about instead of just Z and W. I think I will skip this part of his question, but it is a reasonable one for him to ask. Exercise 7.7 Obtain the Weyl formula for an n=3 region. We did n = 2 in the text. Looking at page 232 things are the same except where C appears, we have the quantity (4πt)3/2 in place of (4πt) as shown 5.140, as in page 232 C. I think everything else on that page is unchanged. We make this same change (4πt)→ (4πt)3/2 in 7.99 and AR → SR. Page 233 A stays the same. B stays the same. And C is the same as well. The RHS of 7.100 becomes RHS 7.100 = VR/ [(4πt)3/2 t] = VR/ [ (4π)3/2 t5/2] The RHS used to be J(t) ≈ AR/(4πt2) and now it is J(t) ≈ VR/ [ (4π)3/2 t5/2] and J(s) ≈ VR/ [ (4π)3/2 s5/2] Doing exactly what I did above in page 233 notes, we then have N(t) = (1/2πi) !Syntax Error, Ids est J(s) ≈ VR/(4π)3/2 (1/2πi) !Syntax Error, Ids est/s5/2 N(t) = VR/(4π)3/2 (1/2πi) !Syntax Error, Idz eztz-5/2 So this time we have not a order 2 pole but a branch point to worry about. On scratch I think I showed that the discontinuity (top - bottom) across this cut is ( the factor ezt is analytic) Δz-5/2 = -2 |z|-5/2 sin(π5/2) = -2 |z|-5/2 sin(π/2) = -2 |z|-5/2 Δ(eztz-5/2) = ezt (-2 |z|-5/2) But the cut runs x = (-∞,0) and on the cut we have z = x and |z| = x, so we have ∫C dz eztz-5/2 = - !Syntax Error, Idx ext (-2x-5/2) = 2 !Syntax Error, Idx ext x-5/2 = 2 (-1)1/2 !Syntax Error, Idx e-xt x-5/2 = 2i OUCH The integral diverges at x = 0. Does this mean the inverse transform does not exist? No, it just means that I failed to include the "turn" around z = 0 which cancels the infinity of the cut. Suppose I try to just look it up and see what gives. We have f(s) = 1/s5/2 and Schaum p 164 32.28 tells us that F(t) = t3/2/Γ(5/2) = t3/2 / [(3/4)] so we conclude that N(t) = VR/(4π)3/2 {(1/2πi) !Syntax Error, Ids est/s5/2} = VR t3/2/(4π)3/2 / [(3/4)] = VR * (1/4π) * (1/) * (4/3) * (1/ ) * t3/2 = VR * (1/4π) * (1/2π) * (4/3) * t3/2 = VR * (1/8π2) * (4/3) * t3/2 = VR * (1/2π2) * (1/3) * t3/2 = VR * (1/6π2) * t3/2 => N(λ) = (VR/6π2) λ3/2 and this agrees with his given answer to this problem for the leading term. I will skip the part of the exercise to get the next term in the series, I could just mimic what he did on page 234, but I am not that interested in this Weyl's formula stuff right now, it is too far off my path right now. Now if we take a cubic box, I think we get λnmk = (π/l)2(n2 + m2+ k2) for the eigenvalues. We might write this as λn = (π/l)2 |n|2 where now n is a vector. Notice that = (π/l) n. The eigenvalues are integer lattice points inside the sphere of radius |n| and the number of such lattice points is perhaps (4/3)π n3 . Then N(λn) = (4/3)π n3 = (4/3)π [l/π ]3 = (4/3)π λn3/2 l3 /π3 = (4/3π2) VR λn3/2 But, counting each eigenvalue only once means we want 1/8th of the sphere, just the first octant where all the nx, ny and nz are ≥ 0. Then we add a 1/8 factor to get N(λn) = (1/8) (4/3π2) VR λn3/2 = (1/6π2) VR λn3/2 which now agrees with our leading term in the Weyl formula. If n is large, this formula becomes exact. But if n is not so large, then yes, there is some kind of correction term since the lattice spacing is not tiny relative to the sphere radius n. It would no doubt take some work to show that this correction somehow corresponded to the non-leading Weyl term. Exercise 7.8 Obtain the Weyl formula for Neumann and . (a) Well, 7.96 has this new BC, I still call the EF's φi and so ∂nφi(xσ) = 0. The EF's are still orthogonal so 7.98 still valid. Rest of page 232 stays valid. The g of 7.97 has our new ∂ng = 0 BC. I just don't see anything at all different. The only fact of the φi ever used is orthogonality. At small t, we stay away from boundaries so the BC's don't matter! (b) For the radiative BC, ∂nφi(xσ) + θ φi(xσ)= 0 is still a standard issue "homogeneous BC" so I think we have exactly the same conclusion. I am talking here only about the leading term! But as I look at the second term derivation, I see no change there either. Exercise 7.9 Another way to show max and min of harmonic u both lie on the boundary. (a) In the p 224 notes above, I claimed that no point inside R can v exceed the max of v on the boundary, and I commented that this means the max must occur on the boundary. This was for a solution v of the heat conduction equation. But if we take the limit to very large t where we achieve steady state (all those exponentials decay to very smallness), the ∂tv term can be neglected and so our same conclusion applies to a solution as well. Now the question is whether the thing inside | | is positive or negative at the max. (a) If 2v > 0 inside R, then there can be no interior points where the function cups down, so there can be no interior "hilltops". Thus, the max must occur on the boundary. (b) Suppose 2u = 0. Then define v according to v(r) ≡ u(r) + ε |r|2 with ε > 0. then 2v = ε 2r2. But 2r2 = Σi∂i2 (Σjxj2) = Σij ∂i(2xjδij)= Σi ∂i(2xi) = Σi 2 = 2n n # dimensions Thus 2v = ε 2r2 = ε 2n > 0. Therefore by part (a) the max of v occurs on the boundary. Therefore, the max of the quantity u(r) + ε |r|2 must occur on the boundary. But this is true no matter how small ε is, so in the limit this says the max of u also occurs on the boundary. We can restate the whole works by assuming 2v < 0 => min of v is on the boundary. Then do same trick and so min u is on the boundary. Therefore both the min and the max of a solution must be on the boundary. We already showed this in Chapter 6, but here is another proof. Exercise 7.10 Finite rod with dipole source in center. Finite rod, u = 0 at both ends, our source is a unit dipole in the middle. We worried about finite rods on page 212. The Green's function was 7.50. So slap this g into 7.14. The rod's spatial boundary is the two endpoints where u = 0 all the time, so think h = 0 in 7.14. Also, we have no distributed sources. Our initial condition is f(x) = -δ'(x-l/2). [ the reason for the minus sign is that this produces a pulse which is positive on the right side, which is really what you mean by a dipole.] So we have only the middle term in 7.14, to wit u(x,t) = - !Syntax Error, Idy g(x,t | y,0) δ'(y-l/2) If we do parts, the parts will be g(x,t | y,0) (y-l/2) which is of course 0 at both ends just from the δ, so we can throw the derivative over to get u(x,t) = + !Syntax Error, Idy ∂yg(x,t | y,0) δ(y-l/2) = + ∂yg(x,t | y,0)|y=l/2 . So g(x,t | y,0) = (2/l) Σn sin(nπx/l)sin(nπy/l) exp(- n2π2t/l2) –∂yg(x,t | y,0) = (2/l) Σn sin(nπx/l)(nπ/l)cos(nπy/l) exp(- n2π2t/l2) u(x,t) = (2/l) Σn=1∞ sin(nπx/l)(nπ/l) cos(nπ/2) exp(- n2π2t/l2) Now I am happy to say cos(nπ/2) = 1 n = 0 cos(nπ/2) = 0 n = 1 cos(nπ/2) = -1 n = 2 cos(nπ/2) = 0 n = 3 cos(nπ/2) = (-1)n/2 for n = even and = 0 for n = odd, so u(x,t) = (2/l) Σn=0,2,4.. sin(nπx/l)(nπ/l) (-1)n/2 exp(- n2π2t/l2) So let n = 2m and we have u(x,t) = (2/l) Σm=0∞ sin(2mπx/l)(2mπ/l) (-1)m exp(- 4m2π2t/l2) and then replace m by n u(x,t) = (2/l) Σn=0∞ sin(2nπx/l)(2nπ/l) (-1)n exp(- 4n2π2t/l2) = Σn=1∞ (-1)n (4nπ/l2) sin(2nπx/l) exp(- 4n2π2t/l2) where we note that the n=0 term in the sum is 0. I think Stak is off by a factor 1/l . As a check, what do we expect for large t? This is dominated by the n=1 term and we get u(x, large t) = - (4π/l2) sin(2πx/l) exp(- 4π2t/l2) I would argue the sign is right because for this shows a positive wave on the right half of the rod if you draw the picture, and that is the side that had the positive half of the dipole. Easy to see that it is negative at small x and 0 in the middle and both ends, etc, draw a picture. As for the image method, take Figure 7.4 and alternate little dipoles as shown but centered in each interval. This gives the right initial condition for (0,l) rod, and makes sure u = 0 at both ends. Then instead of 7.46. [ ??] So we want to add up the responses of a bunch of dipoles on an infinite rod. But at the top of page 201 we showed that a positive dipole at x = 0 makes ∂x C(x,t |0,0), so I think you would get this result = Σn=1/2,5/2...-3/2,-7/2... ∂x C(x-nl,t |0,0) – Σn=3/2,7/2...-1/2,-5/2... ∂x C(x-nl,t |0,0) So let m = 2n and we have result = Σm=1,5...-3,-7... ∂x C(x-ml/2,t |0,0) – Σm=3,7...-1,-5... ∂x C(x- ml/2,t |0,0) Now in the upper sum let m = 1 + 4k where k = all integers. In the lower sum let m = -1 + 4k. So we get result = Σk ∂x C(x-[1+4k]l/2,t |0,0) – Σk ∂x C(x- [-1+4k]l/,t |0,0) = Σk { ∂x C(x-[1+4k]l/2,t |0,0) – C(x- [-1+4k]l/,t |0,0) } and this then is something I could figure out. I would then apply his same series conversion idea. I see I need to add more notes for page 212, since I went off and did that series conversion trick elsewhere and then I forgot to carry through on page 212. But even if I do no more with this Ex 7.10, I see how to do the image part, and I have done the first part and found a small Stak errors. It is probably not too easy coming up with image problems! Exercise 7.11 Average Temperature is a constant in time. This problem took me a long time because I did not realize you need to use an insulated sphere to have the problem make sense. In order to even talk about "average temperature" we have to restrict to a sphere of radius r which at least has some measureable finite volume. Then in the end we might take r→∞. With such a finite sphere, we do have spatial boundaries to worry about. Suppose we approach this problem first assuming we have a f(x) = u(x,0) = δ(x-x0) as our initial condition in the sense of 7.8a. Can we learn something about average temperature in this case? Certainly at t = 0 we can compute the average temperature. [ for such a problem, u(x,t) = g(x,t|x0,0) ] < u(x,0) > = 1/Vn(r) ∫sphere dx δ(x-x0) = 1/Vn(r) The question then is whether this changes as time moves forward. Notice that our source at x0 is in general not at the center of our sphere. The solution to this problem is then the causal Green's function for a sphere, something I don't think I know. Do we need any kind of boundary conditions on the surface of our sphere? It seems that a good one would be ∂ng = 0 on the boundary which says we have an insulated sphere and we cannot lose heat outwards through the boundary. That should certainly be helpful if we want to have the average temperature in the sphere be a constant. I could go and try to compute this Green's function which could take me a week. If I knew this Green's function I could then compute its average at any instant in time and see if it changes. This is the direct method. I do know that the sphere has some potential theory eigenfunctions we could call ψi(x) and I know that the Green's Function I seek can be written this way (similar to bilinear 7.57 on page 214). g(x,t|x0,0) = Σi ψi(x) ψi(x0) exp (-μit) [ This is one way to write our unknown propagator for the Neumann BC sphere which I just said I don't know in detail, but here it is in the generic bilinear sense.] One thing we learned on that page is that the lowest eigenvalue is μ0 = 0. I think this means that the corresponding eigenfunction is a constant. Let's confirm that in a little digression: Digression: Suppose 2f = 0 in spherical coordinates. [ This is the eigenvalue equation in the case that the eigenvalue is 0.] Smythian interior form is then f(r,θ,φ) = Σnm fnm (r/a)n Pnm(θ) e-imφ ∂rf = Σnm fnm n(r/a)n-1(1/a) Pnm(θ) e-imφ ∂rf|r=a = Σnm fnm n (1/a) Pnm(θ) e-imφ = 0 The spherical harmonics form a complete set in (θ,φ) so we must have fnm n (1/a) = 0 This says that all fnm with n > 0 must vanish, and our only survivor is f00. So our only possible solution to this problem is then f(r,θ,φ) = constant. QED Now, the average temperature is given by [ Here I use u(x,t) = g(x,t|x0,0) being the inside-sphere solution in our special case that we had a localized initial condition f(x) = δ(x-x0) somewhere in the sphere. ] < u(x,t) > = (1/Vn(r)) ∫sphere dx ψi(x) ψi(x0) exp (-μit) = (1/Vn(r)) Σi ψi(x0) exp (-μit) ∫sphere dx ψi(x) Now assume our eigenfunctions are orthonormal. Then we have ψ0(x) = 1/ => 1 = ψ0(x) Then we can write [ here we use orthogonality of the eigenfunctions ] ∫sphere dx ψi(x) = ∫sphere dx ψi(x) ψ0(x) = δi,0 Note that r is a constant in this discussion, the radius of our spherical region R. Therefore, < u(x,t) > = (1/Vn(r)) Σi ψi(x0) exp (-μit) ∫sphere dx ψi(x) = (1/Vn(r)) Σi ψi(x0) exp (-μit) δi,0 = (1/Vn(r)) [ψ0(x0) ]exp (-μ0t) = (1/Vn(r)) * 1 * 1 = (1/Vn(r)) Thus we have shown that for f(x) = u(x,0) = δ(x-x0), our average temperature is a constant over time, and here it is equal to what it started off being. So what we have really shown here so far is this: < g(x,t | x0,0) > = (1/Vn(r)) for x0 lying inside our sphere. We know this fact even though we don't know g(x,t | x0,0). [ Maybe at this point we could say: "By superposition, we could build up any f(x) IC we wanted using weighted δ's, and for each contribution, we get < contribution > = constant by the above argument, therefore this must be true for the general superposed f(x), and thus we have shown that in general < u(x,t) > = constant for any f(x) IC. " ] Now finally let's look at 7.14 for our Neumann choice of BC. We have ( with q = 0) u(t,x) = !Syntax Error, Idnx0 {f(x0) g(x|x0)}|t0=0 + !Syntax Error, Idt0 !Syntax Error, I dSn0 g(x|x0) h(t0,x0) (7.14 Neumann) Putting in our spatial averaging we get < u(x,t) > = !Syntax Error, Idnx0 {f(x0) <g(x|x0)>}|t0=0 + !Syntax Error, Idt0 !Syntax Error, I dSn0 <g(x|x0)> h(t0,x0) The first term gives ∫dxof(x0) (1/Vn(r)). This is independent of t. And since f(x0) is bounded by some M, we know that the integral converges and is ≤ M. So our problem now is dealing with the second term. What does this second term even mean? h in this case is a prescribed boundary heat source distribution, meaning ∂nu(x,t) = h(x,t) on σ. But we don't want heat flowing through our boundary for a general problem for u(x,t) inside sphere r if we want to show that average temperature is constant! It would not be constant if we had loss or gain of heat through the boundary. So we have h = 0 and then this entire second term goes away. We have an insulated sphere! We then find that < u(x,t) > = (1/Vn(r)) ∫dxof(x0) which is our grandiose final result. We know that |< u(x,t) >| = (1/Vn(r)) ∫dxo|f(x0)| ≤ M Thus, as we take r→∞ to get our entire space, we don't have to worry about < u(x,t) > diverging. Not surprisingly, the average temperature is ≤ the maximum it has at any point at t=0 in the case that our boundary is insulated. In the text on page 215 we showed that as t→∞, u approaches this average temperature in an insulated volume. It seems pretty reasonable that the average temperature cannot change inside an insulated sphere, but we don't really have Stakgold physics arguments to show that. It is really energy conservation, but here we are showing it from the heat conduction equation theory. Exercise 7.12 Add a cu linear term into the heat equation and resolve in EF method. Let's do the separation of variables u = X T so we have X∂tT - T2U + cXT = 0 ∂tT/T - 2U/U + c = 0 Suppose we put c into the time problem. Then we write ∂tT/T + c = -λ ∂tT = -(λ+c)T -2U/U = λ -2U = λ U This way, the harder problem (the space one) is still just eigenvalues as in the page 214 text. We take the functions φi which vanish on our boundary. Meanwhile, the T problem is exp(-[λ+c]t). Then if there is a spectrum λi the Smythian form will be u(x,t) = Σi ci φi(x) e-(λi+c)t // the new 7.54 Now at t = 0 we have f(x) = Σi ci φi(x) etc etc and we end up with 7.56 looking like u(x,t) = Σi <f,φi> φi(x) e-(λi+c)t and then f(x) = δ(x-x0) gives us the Green's g(x,t |x0,0) = Σi i(x) φi(x0) e-(λi+c)t which replaces p 214 A So really c just enhances the natural λi ability to reduce temperature or particle density. So one way to write the general solution is 7.14 (Dirichlet version) and we get u(x,t) = the first two terms only of 7.14 since u = h = 0 on the boundary. with g(x,t |x0,0) as stated above. So we allow for some q(x,t) sources, and f(x0) initial condition. Exercise 7.13 Maximum principle for fancy diffusion equation. The k term was explained in the heat flow appendix in heat language. The drift term appeared earlier on page 231 top (Fokker-Plank), the c term we just saw. So here we have a "fancy" diffusion equation which allows for "creation or destruction" of particles, allows for a drift velocity α, and the diffusion coefficient k is a function of position. Eventually I will need to find a book that derives this stuff. That in itself is an interesting question: what kind of "physics" talks about diffusion? It does not fit into any of my usual places, mechanics, E&M, QM, particles, plasma, atomic, nuclear? One place is "kinetic theory of gases" in classical sense. But of course you can have diffusion of holes in a semiconductor. I see there is a lot of chemistry and biology chatter on this subject, and that diffusion is associated with the name Fick, his 1st and 2nd Laws, his year 1855. I think the "drift" thing is a externally imposed velocity of your particles. What happens when you drop dye into a river? OK now is not the time, however. The goal of this problem is to show that, if c = 0 and k > 0, the more complicated drift equation still respects the Maximum Principle we found for the regular heat equation without these new complications. Stak outlines a little proof and leaves a few small steps to the reader. I have not done this problem but it seems reasonable to do. If c > 0 you are destroying particles (reducing temperature) , you would expect the max stuff to still work, since c acts to reduce u. But not the min. And vice versa if c < 0. Exercise 7.14. Show Theorem 1 page 223 for Neumann and for 's. Theorem 1 of Section 7.4 says that if Dirichlet u=0 everywhere on a closed boundary (IC + BC), then continuous heat conduction u = 0 inside. this is similar to the equation fact. Here we are supposed to show this is also true if we assume Neumann or Radiative boundary conditions (θ>0). Exercise 7.15 Show Theorem 1 page 223 for the fancy diffusion equation 7.108 if c ≥ 0. This sounds like more of the same. Put these two problems on the rainy day list. I don't see any big issues here. Good to know the conclusions. Exercise 7.16. The half-rod driven from the left end by various methods. This problem bogged me down for 11 days, see separate doc. I did all of this problem except for the second last sentence which says "Get same result using Fourier cosine stuff". Part (a) states a very special problem where we have e-iωt time on both the q and h BC for the half rod problem. If we look for a separated XT solution with T = e-iωt, we do NOT find ourselves doing the spatial EV problem as I once thought. Rather, we discover a certain steady state solution for the rod that has been driven at x=0 by Aeiωt for all time. That solution is p 242A, and the real part is 242 B (also a solution as I show) which is the actual physical solution: a damped wave travelling to the right, and I made a Maple animation to show this. For this solution, the IC f(x) just comes out a certain sinusoidal way (a photo of the steady state solution at t=0). You cannot also specify f(x) because this problem is not in the same "class" as our generic 7.14 type problem. Here we forced a certain solution class by assuming separation with T = e-iωt. There are no "transients" as part of the solution because this is a steady-state solution, though Stak does not use that term. Part (b) treats the same problem in the regular class we are used to where f(x) = 0 (though he neglects to say this). Using a time transform, I am able to show result 242C where both terms solve the homo heat equation, and the second term is what causes us to achieve u(x,0) = 0. This is a typical situation, you add a homo solution to tune the boundary conditions, here the IC boundary condition. I did find the leading term of the transient integral for large t. Part (c) I call it asks for the same solution via Fourier Cosine spatial transform. I first erroneously did this with a transform and ran into various tricky problems (most of my days were here), and I concluded that you cannot solve the problem this way! I then much later tried the Fourier Cosine method and decided it was "illogical", and so after 11 days I decided I had had enough with this "exercise". Exercise 7.17. Deriving the Weber Transform I am not going to do this because I spent so much time on 7.16. Also, I can see that 7.17 is a fairly standard (though singular at one end) Sturm-Liouville 1D problem. We take the Bessel ODE and look for a Green's Function on the interval (a,∞), which is the new feature here (finite lower endpoint). We need g = 0 at the left end of the interval and g = finite s-norm (limit point) at the right end -- this is what p 242 E says. That solution is found to be p 242 F which involves Bessel functions including a special linear combination of J0 and H0(1) which Stak calls Z0 ( this is a "Smythian form" combination which vanishes at x = a) . The Lλ type L is assumed, which allows us to obtain the spectrum in the usual manner and completeness comes out as in p 243 D. The resulting transform is called The Weber Transform and I imagine one would use it for an exterior problem in cylindrical geometry where u = 0 on the cylinder of radius a. I looked in Google books and found there are three errors on page 243, two of which are the same minus sign, and one of which is the omission of the projection in F (which led me to peruse Google). Exercise 7.18. Practice with the Weber Transform We are supposed to first derive 7.76 by doing the required inverse Laplace Transform, and then we are to redo this page 222 cylinder problem using this new Weber transform. This problem involved an infinite medium with an infinite cylindrical hole taken out where u = 1 ( h = 1 reservoir BC) on this inner boundary and u = 0 in the medium (IC f = 0). There we did a time Laplace transform and got a fancy answer 7.76 using the ILT, whereas here we want to do this thing with a spatial Weber transform. I presume that the starting move is to write p 243 F after expanding u(x,t) this way ( a = 1) u(ρ,t) = !Syntax Error, Idμ μ Z0(ρμ) U(μ,t) / [ J02(μ) + N02(μ)] // expansion U(μ,t) = - !Syntax Error, Idρ ρ Z0(ρμ) u(ρ,t) // projection We then end up with ( from 7.75) ( this is very rough, just trying to get a feel here, I am assuming probably wrongly that I can push the differential operator through the integral, and this very thing has caused big problems elsewhere! ) -ρ-1∂ρ(ρ ∂ρ Z0(ρμ) ) = ρ-1 μ-1 { -[ ρ Z0'(ρμ) ] '} = + ρ-1 μ-1 λ ρ Z0(ρμ) // from 7.75 = (μ2/μ)Z0(ρμ) = μ Z0(ρμ) So then 7.75 becomes !Syntax Error, Idμ μ Z0(ρμ) ∂tU(μ,t) / [ J02(μ) + N02(μ)] = - !Syntax Error, Idμ μ μ Z0(ρμ) U(μ,t) / [ J02(μ) + N02(μ)] and then maybe we are allowed to conclude that ∂tU(μ,t) = -μU(μ,t) => U(μ,t) = Q(μ)e-μt The BC on U(μ,t) I think is U(μ,0) = 0 but that forces Q = 0 which must be wrong and I think this gets back to the "pushing through" mentioned above. So I will just halt this derivation for now and save the rest for the endless rainy day. But you can see how this transform is going to make the time ODE be quite simple and that is why the Weber transform is appropriate to this problem. A better start is probably to apply the transform integral operator to both sides of 7.75 and maybe do some parts. I think then the first term is still ∂tU(μ,t) while the second term would be - !Syntax Error, Idρ ρ Z0(ρμ) { -ρ-1∂ρ(ρ ∂ρ u(ρ,t) ) = !Syntax Error, Idρ Z0(ρμ) {∂ρ(ρ ∂ρ u(ρ,t) ) and then we would do two sequential parts integrations. In the first one, we get no parts because Z0 vanishes at ρ = 1 by design. So we have at that point = - !Syntax Error, Idρ [ ρ ∂ρZ0(ρμ)] (∂ρ u(ρ,t) ) We then attempt a second parts and this time we have some parts: = { [ ρ ∂ρZ0(ρμ)] u(ρ,t)} ρ=1 + !Syntax Error, Idρ ∂ρ [ ρ ∂ρZ0(ρμ)] u(ρ,t) Then perhaps we set u(ρ=1,t) = 1 and replace ( see above) ∂ρ(ρ ∂ρ Z0(ρμ) ) = - μρ Z0(ρμ) to get = [ ∂ρZ0(ρμ)] ρ=1 + !Syntax Error, Idρ [ - μρ Z0(ρμ)] u(ρ,t) = [ ∂ρZ0(ρμ)] ρ=1 -μ !Syntax Error, Idρ ρ Z0(ρμ) u(ρ,t) = [ ∂ρZ0(ρμ)] ρ=1 + μ U(μ,t) So we now have an "improved" version of our earlier equation ∂tU(μ,t) = - [ ∂ρZ0(ρμ)] ρ=1 - μ U(μ,t) and hopefully having this extra constant on the RHS removes our Q paradox from above. I am not going to continue here because it would cost me perhaps 8 hours and I need to "move on". I have just been blindly winging it here and this problem really needs a rigorous presentation in a separate document. End of "heat final exam". This concludes my Stakgold efforts on the heat conduction equation system. we are at page 69 in these raw notes. This covers about 50 pages of the chapter text, and we have another 90 pages to go, so I will open up a new doc.