The Heat and Wave Equations in 2D and 3D

The heat and wave equations in 2D and 3D 18.303 Linear Partial Differential Equations Matthew J. Hancock Fall 2006 1 2D and 3D Heat Equation Ref: Myint-U & Debnath §2.3 – 2.5 [Nov 2, 2006] Consider an arbitrary 3D subregion V of R3 (V ⊆ R3 ), with temperature u (x, t) defined at all points x = (x, y, z) ∈ V . We generalize the ideas of 1-D heat flux to find an equation governing u. The heat energy in the subregion is defined as � � heat energy = cρu dV V Recall that conservation of energy implies rate of change heat energy into V from heat energy generated = + of heat energy boundaries per unit time in solid per unit time We desire the heat flux through the boundary S of the subregion V , which is the normal component of the heat flux vector φ, φ · n̂, where n̂ is the outward unit normal at the boundary S. Hats on vectors denote a unit vector, |n̂| = 1 (length 1). If the heat flux vector φ is directed inward, then φ · n̂ < 0 and the outward flow of heat is negative. To compute the total heat energy flowing across the boundaries, �� we sum φ · n̂ over the entire closed surface S, denoted by a double integral dS. S Therefore, the conservation of energy principle becomes � � � � � � � � d cρu dV = − φ · n̂ dS + Q dV (1) dt V S V 1 1.1 Divergence Theorem (a.k.a. Gauss’s Theorem) For any volume V with closed smooth surface S, � � � � � ∇ · A dV = A · n̂ dS V S where A is any function that is smooth (i.e. continuously differentiable) for x ∈ V . Note that � � ∂ ∂ ∂ ∂ ∂ ∂ , , = êx + êy + êz ∇= ∂x ∂y ∂z ∂x ∂y ∂z where êx , êy , êz are the unit vectors in the x, y, z directions, respectively. The divergence of a vector valued function F = (Fx , Fy , Fz ) is ∇·F= ∂Fx ∂Fy ∂Fz + + . ∂x ∂y ∂z The Laplacian of a scalar function F is ∇2 F = ∂2F ∂2F ∂2F + + . ∂x2 ∂y 2 ∂z 2 Applying the Divergence theorem to (1) gives � � � � � � � � � d cρu dV = − ∇ · φ dV + Q dV dt V V V Since V is independent of time, the integrals can be combined as � � � � � ∂u cρ + ∇ · φ − Q dV = 0 ∂t V Since V is an arbitrary subregion of R3 and the integrand is assumed continuous, the integrand must be everywhere zero, cρ 1.2 ∂u +∇·φ−Q=0 ∂t (2) Fourier’s Law of Heat Conduction The 3D generalization of Fourier’s Law of Heat Conduction is φ = −K0 ∇u (3) where K0 is called the thermal diffusivity. Substituting (3) into (2) gives ∂u Q = κ∇2 u + ∂t cρ 2 (4) where κ = K0 / (cρ). This is the 3D Heat Equation. Normalizing as for the 1D case, x̃ = x , l t̃ = κ t, l2 Eq. (4) becomes (dropping tildes) the non-dimensional Heat Equation, ∂u = ∇2 u + q, ∂t (5) where q = l2 Q/ (κcρ) = l2 Q/K0 . 2 2D and 3D Wave equation The 1D wave equation can be generalized to a 2D or 3D wave equation, in scaled coordinates, utt = ∇2 u (6) This models vibrations on a 2D membrane, reflection and refraction of electromagnetic (light) and acoustic (sound) waves in air, fluid, or other medium. 3 Separation of variables in 2D and 3D Ref: Guenther & Lee §10.2, Myint-U & Debnath §4.10, 4.11 We consider simple subregions D ⊆ R3 . We assume the boundary conditions are zero, u = 0 on ∂D, where ∂D denotes the closed surface of D (assumed smooth). The 3D Heat Problem is ut = ∇2 u, u (x, t) = 0, u (x, 0) = f (x) , x ∈ D, t > 0, x ∈ ∂D, (7) x ∈ D. The 3D wave problem is utt = ∇2 u, u (x, t) = 0, u (x, 0) = f (x) , x ∈ D, x ∈ ∂D, t > 0, (8) x ∈ D, ut (x, 0) = g (x) , x ∈ D. We separate variables as u (x, t) = X (x) T (t) 3 (9) The 3D Heat Equation implies T′ ∇2 X = = −λ = const T X (10) where λ = const since the l.h.s. depends solely on t and the middle X ′′ /X depends solely on x. The 3D wave equation becomes ∇2 X T ′′ = = −λ = const T X (11) On the boundaries, X (x) = 0, x ∈ ∂D The Sturm-Liouville Problem for X (x) is ∇2 X + λX = 0, X (x) = 0, 4 x∈D (12) x ∈ ∂D Solution for T (t) Suppose that the Sturm-Liouville problem (12) has eigen-solution Xn (x) and eigen­ value λn , where Xn (x) is non-trivial. Then for the 3D Heat Problem, the problem for T (t) is, from (10), T′ = −λ (13) T with solution Tn (t) = cn e−λn t (14) and the corresponding solution to the PDE and BCs is un (x, t) = Xn (x) Tn (t) = Xn (x) cn e−λn t . For the 3D Wave Problem, the problem for T (t) is, from (11), with solution T ′′ = −λ T (15) �� � �� � λn t + βn sin λn t Tn (t) = αn cos (16) and the corresponding normal mode is un (x, t) = Xn (x) Tn (t). 4 5 Uniqueness of the 3D Heat Problem Ref: Guenther & Lee §10.3 We now prove that the solution of the 3D Heat Problem ut = ∇2 u, x∈D u (x, t) = 0, u (x, 0) = f (x) , x ∈ ∂D x∈D is unique. Let u1 , u2 be two solutions. Define v = u1 − u2 . Then v satisfies vt = ∇2 v, v (x, t) = 0, v (x, 0) = 0, Let V (t) = x ∈ ∂D x∈D � � � D 2 x∈D v 2 dV ≥ 0 V (t) ≥ 0 since the integrand v (x, t) ≥ 0 for all (x, t). Differentiating in time gives � � � dV (t) = 2vvt dV dt D Substituting for vt from the PDE yields � � � dV (t) = 2v∇2 vdV dt D By result (26) derived below, � � � � � dV (t) = 2 v∇v · n̂dS − 2 |∇v|2 dV dt ∂D D But on ∂D, v = 0, so that the first integral on the r.h.s. vanishes. Thus � � � dV (t) = −2 |∇v|2 dV ≤ 0 dt D Also, at t = 0, V (0) = � � � (17) (18) v 2 (x, 0) dV = 0 D Thus V (0) = 0, V (t) ≥ 0 and dV /dt ≤ 0, i.e. V (t) is a non-negative, non-increasing function that starts at zero. Thus V (t) must be zero for all time t, so that v (x, t) must be identically zero throughout the volume D for all time, implying the two solutions are the same, u1 = u2 . Thus the solution to the 3D heat problem is unique. For insulated BCs, ∇v = 0 on ∂D, and hence v∇v · n̂ = 0 on ∂D. Thus we can still derive Eq. (18) from (17), and the uniqueness proof still holds. Thus the 3D Heat Problem with Type II homogeneous BCs also has a unique solution. 5 6 Sturm-Liouville problem Ref: Guenther & Lee §10.2, Myint-U & Debnath §7.1 – 7.3 Both the 3D Heat Equation and the 3D Wave Equation lead to the Sturm-Liouville problem ∇2 X + λX = 0, x ∈ D, X (x) = 0, 6.1 (19) x ∈ ∂D. Green’s Formula and the Solvability Condition Ref: Guenther & Lee §8.3, Myint-U & Debnath §10.10 (Exercise 1) For the Type I BCs assumed here (u (x, t) = 0, for x ∈ ∂D), we now show that all eigenvalues are positive. To do so, we need a result that combines some vector calculus with the Divergence Theorem. From vector calculus, for any scalar function G and vector valued function F, ∇ · (GF) = G∇ · F + F · ∇G (20) Using the divergence theorem, � � � � � (GF) · n̂dS = ∇ · (GF) dV S Substituting (20) into (21) gives � � � � � � � � (GF) · n̂dS = G∇ · FdV + F · ∇GdV S (21) V V (22) V Result (22) applied to G = v1 and F = ∇v2 gives � � � � � � � (v1 ∇v2 ) · n̂dS = v1 ∇2 v2 + ∇v2 · ∇v1 dV S (23) V Result (22) applied to G = v2 and F = ∇v1 gives � � � � � � � (v2 ∇v1 ) · n̂dS = v2 ∇2 v1 + ∇v1 · ∇v2 dV S (24) V Subtracting (23) and (24) gives Green’s Formula (also known as Green’s second iden­ tity): � � � � � � � (v1 ∇v2 − v2 ∇v1 ) · n̂dS = v1 ∇2 v2 − v2 ∇2 v1 dV (25) S V This is also known as a Solvability Condition, since the values of v1 and v2 on the boundary of D must be consistent with the values of v1 and v2 on the interior of D. Result (25) holds for any smooth function v defined on a volume V with closed smooth surface S. 6 6.2 Positive, real eigenvalues (for Type I BCs) Ref: Myint-U & Debnath §7.2 Choosing G = v and F = ∇v, for some function v, in (22) gives � � � � � � 2 � v∇v · n̂dS = v∇ v + ∇v · ∇v dV S V = � � � V � � v∇2 v + |∇v|2 dV. (26) Result (26) holds for any smooth function v defined on a volume V with closed smooth surface S. We now apply result (26) to a solution X (x) of the Sturm-Liouville problem (19). Letting v = X (x), S = ∂D and V = D, Eq. (26) becomes � � � � � � � � 2 X∇X · n̂dS = X∇ XdV + |∇X|2 dV (27) ∂D D Since X (x) = 0 for x ∈ ∂D, � � ∂D Also, from the PDE in (19), � � � D X∇X · n̂dS = 0 2 D X∇ XdV = −λ (28) � � � X 2 dV Substituting (28) and (29) into (27) gives � � � � � � 2 0 = −λ X dV + |∇X|2 dV D (29) D (30) D For non-trivial solutions, X �= 0 at some points in D and hence by continuity of X, ��� X 2 dV > 0. Thus (30) can be rearranged, D ��� |∇X|2 dV D λ= � � � ≥0 X 2 dV D (31) Since X is real, the the eigenvalue λ is also real. If ∇X = 0 for all points in D, then integrating and imposing the BC X (x) = 0 for x ∈ ∂D gives X = 0 for all x ∈ D , i.e. the trivial solution. Thus ∇X is nonzero ��� at some points in D, and hence by continuity of ∇X, |∇X|2 dV > 0. Thus, D from (31), λ > 0. 7 6.3 Orthogonality of eigen-solutions to Sturm-Liouville problem Ref: Myint-U & Debnath §7.2 Suppose v1 , v2 are two eigenfunctions with eigenvalues λ1 , λ2 of the 3D SturmLiouville problem ∇2 v + λv = 0, x∈D v = 0, x ∈ ∂D We make use of Green’s Formula (25) with V = D, S = ∂D, which holds for any functions v1 , v2 defined on a volume V with smooth closed connected surface S. Since v1 = 0 = v2 on ∂D and ∇2 v1 = −λ1 v1 and ∇2 v2 = −λ2 v2 , then Eq. (25) becomes � � � 0 = (λ1 − λ2 ) v1 v2 dV D Thus if λ1 �= λ2 , � � � v1 v2 dV = 0, (32) D and the eigenfunctions v1 , v2 are orthogonal. 7 Heat and Wave problems on a 2D rectangle, ho­ mogeneous BCs Ref: Guenther & Lee §10.2, Myint-U & Debnath §9.4 – 9.6 [Nov 7, 2006] 7.1 Sturm-Liouville Problem on a 2D rectangle We now consider the special case where the subregion D is a rectangle D = {(x, y) : 0 ≤ x ≤ x0 , 0 ≤ y ≤ y0 } The Sturm-Liouville Problem (19) becomes ∂ 2v ∂ 2v + + λv = 0, (x, y) ∈ D, ∂x2 ∂y 2 v (0, y) = v (x0 , y) = 0, 0 ≤ y ≤ y0 , v (x, 0) = v (x, y0 ) = 0, 8 0 ≤ x ≤ x0 . (33) (34) (35) Note that in the PDE (33), λ is positive a constant (we showed above that λ had to be both constant and positive). We employ separation of variables again, this time in x and y: substituting v (x, y) = X (x) Y (y) into the PDE (33) and dividing by X (x) Y (y) gives Y ′′ X ′′ +λ=− Y X Since the l.h.s. depends only on y and the r.h.s. only depends on x, both sides must equal a constant, say µ, Y ′′ X ′′ +λ=− =µ (36) Y X The BCs (34) and (35) imply X (0) Y (y) = X (x0 ) Y (y) = 0, 0 ≤ y ≤ y0 , X (x) Y (0) = X (x) Y (y0 ) = 0, 0 ≤ x ≤ x0 . To have a non-trivial solution, Y (y) must be nonzero for some y ∈ [0, y0 ] and X (x) must be nonzero for some x ∈ [0, x0 ], so that to satisfy the previous 2 equations, we must have X (0) = X (x0 ) = Y (0) = Y (y0 ) = 0 (37) The problem for X (x) is the 1D Sturm-Liouville problem X ′′ + µX = 0, 0 ≤ x ≤ x0 (38) X (0) = X (x0 ) = 0 We solved this problem in the chapter on the 1D Heat Equation. We found that for non-trivial solutions, µ had to be positive and the solution is � � �2 � mπ mπx , µm = , m = 1, 2, 3, ... (39) Xm (x) = am sin x0 x0 The problem for Y (y) is Y ′′ + νY = 0, 0 ≤ y ≤ y0 , (40) Y (0) = Y (y0 ) = 0 where ν = λ − µ. The solutions are the same as those for (38), with ν replacing µ: � � � �2 nπy nπ Yn (y) = bn sin , νn = , n = 1, 2, 3, ... (41) y0 y0 The eigen-solution of the 2D Sturm-Liouville problem (33) – (35) is � � � � nπy mπx sin , m, n = 1, 2, 3, ... (42) vmn (x, y) = Xm (x) Yn (y) = cmn sin x0 y0 with eigenvalue � 2 � m n2 2 λmn = µm + νn = π + 2 . x20 y0 9 7.2 Solution to heat equation on 2D rectangle The heat problem on the 2D rectangle is the special case of (7), ∂2u ∂2u + , (x, y) ∈ D, ∂x2 ∂y 2 u (x, y, t) = 0, (x, y) ∈ ∂D, ut = u (x, y, 0) = f (x, y) , t > 0, (x, y) ∈ D, where D is the rectangle D = {(x, y) : 0 ≤ x ≤ x0 , 0 ≤ y ≤ y0 }. We reverse the separation of variables (9) and substitute solutions (14) and (42) to the T (t) problem (13) and the Sturm-Liouville problem (33) – (35), respectively, to obtain � � � � mπx nπy −λmn t umn (x, y, t) = Amn sin sin e x0 y0 „ « � � � � 2 2 nπy −π2 mx2 + ny2 t mπx 0 0 sin e = Amn sin x0 y0 To satisfy the initial condition, we sum over all m, n to obtain the solution, in general form, � � � � ∞ � ∞ ∞ � ∞ � � nπy −λmn t mπx sin e (43) u (x, y, t) = umn (x, y, t) = Amn sin x0 y0 m=1 n=1 m=1 n=1 Setting t = 0 and imposing the initial condition u (x, y, 0) = f (x, y) gives � � � � ∞ � ∞ ∞ � ∞ � � mπx nπy f (x, y) = u (x, y, 0) = Amn vmn (x, y) = Amn sin sin x y0 0 m=1 n=1 m=1 n=1 � mπx x0 � � nπy y0 � where vmn (x, y) = sin sin are the eigenfunctions of the 2D Sturm­ Liouville problem on a rectangle, (33) – (35). Multiplying both sides by vmn̂ ˆ (x, y) (m̂, n̂ = 1, 2, 3, ...) and integrating over the rectangle D gives � � � � ∞ � ∞ � f (x, y) vmn̂ Amn vmn (x, y) vmn̂ (44) ˆ (x, y) dA = ˆ (x, y) dA D D m=1 n=1 where dA = dxdy. Note that � � � vmn (x, y) vmn̂ ˆ (x, y) dA = D � � � mπx ˆ mπx sin dx sin x0 x0 0 � � � � � y0 nπy nπy ˆ × sin sin dy y0 y0 0 � x0 y0 /4, m = m̂ and n = n̂ = 0, otherwise 10 x0 � Thus (44) becomes � � f (x, y) vm̂n̂ (x, y) dA = D Am̂n̂ x0 y0 4 Since m̂, n̂ are dummy variables, we replace m̂ by m and n̂ by n, and rearrange to obtain � � 4 f (x, y) vmn̂ Amn = ˆ (x, y) dA x0 y0 � � � � � x0D� y0 4 mπx nπy = f (x, y) sin sin dydx (45) x0 y0 0 x0 y0 0 7.3 Solution to wave equation on 2D rectangle, homogeneous BCs The solution to the wave equation on the 2D rectangle follows similarly. The general 3D wave problem (8) becomes ∂ 2u ∂ 2u + , (x, y) ∈ D, ∂x2 ∂y 2 u (x, y, t) = 0, (x, y) ∈ ∂D, utt = u (x, y, 0) = f (x, y) , t > 0, (x, y) ∈ D, ut (x, y, 0) = g (x, y) , (x, y) ∈ D, where D is the rectangle D = {(x, y) : 0 ≤ x ≤ x0 , 0 ≤ y ≤ y0 }. We reverse the separation of variables (9) and substitute solutions (16) and (42) to the T (t) problem (15) and the Sturm Liouville problem (33) – (35), respectively, to obtain � � � � �� � �� �� nπy � mπx sin αnm cos umn (x, y, t) = sin λnm t + βnm sin λnm t x0 y0 � � � � nπy mπx sin = sin x0 y0 � � � � � � �� m 2 n2 m 2 n2 × αnm cos πt + 2 + βnm sin πt + 2 x20 y0 x20 y0 To satisfy the initial condition, we sum over all m, n to obtain the solution, in general form, ∞ � ∞ � u (x, y, t) = umn (x, y, t) (46) m=1 n=1 Setting t = 0 and imposing the initial conditions u (x, y, 0) = f (x, y) , ut (x, y, 0) = g (x, y) 11 gives f (x, y) = u (x, y, 0) = ∞ � ∞ � αmn vmn (x, y) m=1 n=1 = ∞ � ∞ � αmn sin m=1 n=1 � mπx x0 � sin � nπy y0 � ∞ � ∞ � � λmn βmn vmn (x, y) g (x, y) = ut (x, y, 0) = m=1 n=1 = ∞ � ∞ � � λmn βmn sin m=1 n=1 � � mπx x0 � sin � nπy y0 � � � � nπy where vmn (x, y) = sin mπx sin are the eigenfunctions of the 2D Sturm Li­ x0 y0 ouville problem on a rectangle, (33) – (35). As above, multiplying both sides by vmn̂ ˆ (x, y) (m̂, n̂ = 1, 2, 3, ...) and integrating over the rectangle D gives � � 4 f (x, y) vmn (x, y) dxdy αmn = x0 y0 D � � 4 √ βmn = g (x, y) vmn (x, y) dxdy x0 y0 λmn D 8 Heat and Wave equations on a 2D circle, homo­ geneous BCs Ref: Guenther & Lee §10.2, Myint-U & Debnath §9.4, 9.13 (exercises) We now consider the special case where the subregion D is the unit circle (we may assume the circle has radius 1 by choosing the length scale l for the spatial coordinates as the original radius): � � D = (x, y) : x2 + y 2 ≤ 1 The Sturm-Liouville Problem (19) becomes ∂ 2v ∂ 2v + + λv = 0, ∂x2 ∂y 2 v (x, y) = 0, (x, y) ∈ D, (47) x2 + y 2 = 1, (48) where we already know λ is positive and real. It is natural to introduce polar coor­ dinates via the transformation x = r cos θ, y = r sin θ, w (r, θ, t) = u (x, y, t) for 0 ≤ r ≤ 1, −π ≤ θ < π. 12 You can verify that ∂2v ∂2v 1 ∂ ∇v= + 2 = 2 ∂x ∂y r ∂r 2 � ∂w r ∂r The PDE becomes � � ∂w 1 ∂ 2w 1 ∂ r + 2 2 + λw = 0, r ∂r ∂r r ∂θ � + 1 ∂2w r2 ∂θ2 0 ≤ r ≤ 1, −π ≤ θ < π (49) The BC (48) requires w (1, θ) = 0, −π ≤ θ < π. (50) We use separation of variables by substituting w (r, θ) = R (r) H (θ) (51) into the PDE (49) and multiplying by r2 / (R (r) H (θ)) and then rearranging to obtain � � d dR 1 d2 H 1 r r + λr2 = − 2 dr dr R (r) dθ H (θ) Again, since the l.h.s. depends only on r and the r.h.s. on θ, both must be equal to a constant µ, � � dR 1 d2 H 1 d r + λr2 = − 2 =µ (52) r dr dr R (r) dθ H (θ) The BC (50) becomes w (1, θ) = R (1) H (θ) = 0 which, in order to obtain non-trivial solutions (H (θ) = � 0 for some θ), implies R (1) = 0 (53) In the original (x, y) coordinates, it is assumed that v (x, y) is smooth (i.e. contin­ uously differentiable) over the circle. When we change to polar coordinates, we need to introduce an extra condition to guarantee the smoothness of v (x, y), namely, that w (r, −π) = w (r, π) , wθ (r, −π) = wθ (r, π) . (54) dH dH (−π) = (π) . dθ dθ (55) Substituting (51) gives H (−π) = H (π) , The solution v (x, y) is also bounded on the circle, which implies R (r) must be bounded for 0 ≤ r ≤ 1. 13 The problem for H (θ) is d2 H + µH (θ) = 0; dθ2 H (−π) = H (π) , dH dH (−π) = (π) . dθ dθ (56) You can show that for µ < 0, we only get the trivial solution H (θ) = 0. For µ = 0, we have H (θ) = const, which works. For µ > 0, non-trivial solutions are found only when µ = m2 , Hm (θ) = am cos (mθ) + bm sin (mθ) Thus, in general, we may assume λ = m2 , for m = 0, 1, 2, 3, ... The equation for R (r) in (52) becomes � � d dR 1 r r + λr2 = µ = m2 , m = 0, 1, 2, 3, ... dr dr R (r) Rearranging gives r2 � d2 Rm dRm � 2 +r + λr − m2 Rm = 0; 2 dr dr Rm (1) = 0, |Rm (0)| < ∞ (57) We know already that λ > 0, so we can let √ ¯ m (s) = Rm (r) R s = λr, so that (57) becomes s2 � d2 R̄m dR̄m � 2 +s + s − m2 R̄m = 0; 2 ds ds R̄m �√ � λ = 0, � � �R̄m (0)� < ∞ (58) The ODE is called Bessel’s Equation which, for each m = 0, 1, 2, ... has two linearly independent solutions, Jm (s) and Ym (s), called the Bessel functions of the first and second kinds, respectively, of order m. The function Jm (s) is bounded at s = 0; the function Ym (s) is unbounded at s = 0. The general solution to the ODE is R̄m (s) = cm1 Jm (s) + cm2 Ym (s) where cmn are constants of integration. Our bound­ � � edness criterion �R̄m (0)� < ∞ at s = 0 implies c2m = 0. Thus �√ � λr . R̄m (s) = cm Jm (s) , Rm (r) = cm Jm [Nov 9, 2006] The Bessel Function Jm (s) of the first kind of order m has power series Jm (s) = � k=0 (−1)k s2k+m . k! (k + m)!22k+m (59) Jm (s) can be expressed in many ways, see Handbook of Mathematical Functions by Abramowitz and Stegun, for tables, plots, and equations. The fact that Jm (s) is 14 expressed as a power series is not a drawback. It is like sin (x) and cos (x), which are also associated with power series. Note that the power series (59) converges absolutely for all s ≥ 0 and converges uniformly on any closed set s ∈ [0, L]. To see this, note that each term in the sum satisfies � � � (−1)k s2k+m � L2k+m � � ≤ � � � k! (k + m)!22k+m � k! (k + m)!22k+m Note that the sum of numbers � k=0 L2k+m k! (k + m)!22k+m converges by the Ratio Test, since the ratio of successive terms in the sum is � � � � �2 � L2(k+1)+m 2 2 � (k+1)!(k+1+m)!22(k+1)+m � L L L � �= � � (k + 1) (k + 1 + m) 4 ≤ (k + 1)2 4 = 2 (k + 1) L2k+m � � 2k+m k!(k+m)!2 Thus for k > N = ⌈L/2⌉, � � � � � �2 L2(k+1)+m � (k+1)!(k+1+m)!2 2(k+1)+m � L � � <1 � � < 2 (N + 1) L2k+m � � k!(k+m)!22k+m Since the upper bound is less than one and is independent of the summation index k, then by the Ratio test, the sum converges absolutely. By the Weirstrass M-Test, the infinite sum in (59) converges uniformly on [0, L]. Since L is arbitrary, the infinite sum in (59) converges uniformly on any closed subinterval [0, L] of the real axis. Each Bessel function Jm (s) has an infinite number of zeros (roots) for s > 0. Let jm,n be the n’th zero of the function Jm (s). Note that 1 2 3 J0 (s) j0,1 = 2.4048 j0,2 = 5.5001 j0,3 = 8.6537 J1 (s) j1,1 = 3.852 j1,2 = 7.016 j1,3 = 10.173 The second BC requires �√ � �√ � λ = Jm λ =0 √ This has an infinite number of solutions, namely λ = jm,n for n = 1, 2, 3, .... Thus the eigenvalues are 2 λmn = jm,n , m, n = 1, 2, 3, ... ¯m Rm (1) = R with corresponding eigenfunctions Jm (rjm,n ). The separable solutions are thus � J0 (rj0,n ) n = 1, 2, 3, ... vmn (x, y) = wmn (r, θ) = Jm (rjm,n ) (αmn cos (mθ) + +βmn sin (mθ)) m, n = 1, 2, 3, ... (60) 15 8.1 Solution to heat equation on the 2D circle, homogeneous BCs The heat problem on the 2D circle is the special case of (7), ∂2u ∂2u + , (x, y) ∈ D, ∂x2 ∂y 2 u (x, y, t) = 0, (x, y) ∈ ∂D, ut = u (x, y, 0) = f (x, y) , t > 0, (x, y) ∈ D, where D is the circle D = {(x, y) : x2 + y 2 ≤ 1}. We reverse the separation of vari­ ables (9) and substitute solutions (14) and (42) to the T (t) problem (13) and the Sturm Liouville problem (47) – (48), respectively, to obtain 2 umn (x, y, t) = vmn (x, y) e−λmn t = vmn (x, y) e−Jm,n t where vmn is given in (60). To satisfy the initial condition, we sum over all m, n to obtain the solution, in general form, u (x, y, t) = = ∞ � ∞ � m=1 n=1 ∞ � ∞ � umn (x, y, t) = ∞ � ∞ � vmn (x, y) e−λmn t m=1 n=1 Jm (rjm,n ) (αmn cos (mθ) + βmn sin (mθ)) e−λmn t (61) m=1 n=1 � where r = x2 + y 2 and tan θ = y/x. Setting t = 0 and imposing the initial condition u (x, y, 0) = f (x, y) gives f (x, y) = u (x, y, 0) = ∞ � ∞ � vmn (x, y) m=1 n=1 We can use orthogonality relations to find αmn and βmn . 9 The Heat Problem on a square with inhomoge­ neous BC We now consider the case of the heat problem on the 2D unit square D = {(x, y) : 0 ≤ x, y ≤ 1} , 16 (62) where a hot spot exists on the left side, ut = ∇2 u, (x, y) ∈ D � u0 /ε, {x = 0, |y − y0 | < ε/2} u (x, y, t) = 0, otherwise on ∂D (63) u (x, y, 0) = f (x, y) where the hot spot is confined to the left side: 0 ≤ y0 − ε/2 ≤ y ≤ y0 + ε/2 ≤ 1. 9.1 Equilibrium solution and Laplace’s Eq. on a rectangle Ref: Guenther & Lee §8.1, 8.2, Myint-U & Debnath §6.6, 8.7 As in the 1D case, we first find the equilibrium solution uE (x, y), which satisfies the PDE and the BCs, ∇2 uE = 0, � uE (x, y) = (x, y) ∈ D u0 /ε {x = 0, |y − y0 | < ε/2} 0 otherwise on ∂D The PDE for uE is called Laplace’s equation. Laplace’s equation is an example of an elliptic PDE. The wave equation is an example of a hyperbolic PDE. The heat equation is a parabolic PDE. These are the three types of second order (i.e. involving double derivatives) PDEs: elliptic, hyperbolic and parabolic. We proceed via separation of variables: uE (x, y) = X (x) Y (y), so that the PDE becomes Y ′′ X ′′ = = −λ − X Y where λ is constant since the l.h.s. depends only on x and the middle only on y. The BCs are Y (0) = Y (1) = 0, X (1) = 0 and X (0) Y (y) = � u0 /ε |y − y0 | < ε/2 0 otherwise We first solve for Y (y), since we have 2 easy BCs: Y ′′ + λY = 0; Y (0) = Y (1) = 0 The non-trivial solutions, as we have found before, are Yn = sin (nπy) with λn = n2 π 2 , for each n = 1, 2, 3, ... Now we consider X (x): X ′′ − n2 π 2 X = 0 17 and hence X (x) = c1 enπx + c2 e−nπx An equivalent and more convenient way to write this is X (x) = c3 sinh nπ (1 − x) + c4 cosh nπ (1 − x) Imposing the BC at x = 1 gives X (1) = c4 = 0 and hence X (x) = c3 sinh nπ (1 − x) Thus the equilibrium solution to this point is uE (x, y) = ∞ � n=1 An sinh (nπ (1 − x)) sin (nπy) You can check that this satisfies the BCs on x = 1 and y = 0, 1. Also, from the BC on x = 0, we have � ∞ � u0 /ε |y − y0 | < ε/2 uE (0, y) = An sinh (nπ) sin (nπy) = 0 otherwise n=1 Multiplying both sides by sin (mπy) an integrating in y gives ∞ � n=1 An sinh (nπ) � 1 sin (nπy) sin (mπy) dy = 0 � y0 +ε/2 y0 −ε/2 u0 sin (mπy) dy ε From the orthogonality of sin’s, we have � y0 +ε/2 1 u0 Am sinh (mπ) = sin (mπy) dy 2 ε y0 −ε/2 Thus, Am 2u0 = ε sinh (mπ) � y0 +ε/2 sin (mπy) dy y0 −ε/2 � �y +ε/2 2u0 cos (mπy) 0 = − ε sinh (mπ) mπ y0 −ε/2 2u0 = (cos (mπ (y0 − ε/2)) − cos (mπ (y0 + ε/2))) εmπ sinh (mπ) � � 4u0 sin (mπy0 ) sin mπε 2 = εmπ sinh (mπ) 18 Thus � � ∞ 4u0 � sin (nπy0 ) sin nπε 2 uE (x, y) = sinh (nπ (1 − x)) sin (nπy) επ n=1 n sinh (nπ) To solve the transient problem, we proceed as in 1-D by defining the function v (x, y, t) = u (x, y, t) − uE (x, y) so that v (x, y, t) satisfies v t = ∇2 v v = 0 on ∂D v (x, y, 0) = f (x, y) − uE (x, y) 9.2 First term approximation To approximate the equilibrium solution uE (x, y), note that sinh nπ (1 − x) enπ(1−x) − e−nπ(1−x) = sinh nπ enπ − e−nπ For sufficiently large n, we have enπ(1−x) sinh nπ (1 − x) ≈ = e−nπx sinh nπ enπ Thus the terms decrease in magnitude (x > 0) and hence uE (x, y) can be approxi­ mated the first term in the series, � � 4u0 sin (πy0 ) sin πε 2 uE (x, y) ≈ sinh (π (1 − x)) sin (πy) επ sinh (π) A plot of sinh (π (1 − x)) sin (πy) is given below. The temperature in the center of the square is approximately � � � � �π � �π � 4u0 sin (πy0 ) sin πε 1 1 2 , ≈ sinh sin uE 2 2 επ sinh (π) 2 2 9.3 Easy way to find steady-state temperature at center For y0 = 1/2 and ε = 1, we have uE � 1 1 , 2 2 � � � 4u0 sinh π2 u0 ≈ ≈ . π sinh (π) 4 19 z 10 7.5 5 1 0.75 2.5 0.5 0.25 y 0 0.25 0.5 0.75 1 x Figure 1: Plot of sinh (π(1 − x)) sin(πy). It turns out there is a much easier way to derive this last result. Consider a plate � � with BCs u = u0 on one side, and u = 0 on the other 3 sides. Let α = uE 21 , 21 . � � Rotating the plate by 90o will not alter uE 21 , 21 , since this is the center of the plate. Let uEsum be the sums of the solutions corresponding to the BC u = u0 on each of the four different sides. Then by linearity, uEsum = u0 on all sides and hence uEsum = u0 across the plate. Thus � � 1 1 = 4α u0 = uEsum , 2 2 � � Hence α = uE 12 , 21 = u0 /4. 9.4 Placement of hot spot for hottest steady-state center Note that uE � 1 1 , 2 2 � � � �π � 4u0 sin (πy0 ) sin πε 2 ≈ sinh επ sinh (π) 2 �π� � πε � 4u0 sinh 2 sin (πy0 ) sin 2 = π sinh (π) ε � � �� u0 2 sin πε 2 ≈ sin (πy0 ) πε 4 π 2 2 u0 sin (πy0 ) 4 π u0 sin (πy0 ) = 2π ≈ 20 for small ε. Thus the steady-state center temperature is hottest when the hot spot is placed in the center of the side, i.e. y0 = 1/2. 9.5 Solution to inhomogeneous heat problem on square We now use the standard trick and solve the inhomogeneous Heat Problem (63), ut = ∇2 u, (x, y) ∈ D � u0 /ε, {x = 0, |y − y0 | < ε/2} u (x, y, t) = 0, otherwise on ∂D (64) u (x, y, 0) = f (x, y) using the equilibrium solution uE . We define the transient part of the solution as v (x, y, t) = u (x, y, t) − uE (x, y) where u (x, y, t) is the solution to (64). The problem for v (x, y, t) is therefore vt = ∇2 v, v (x, y, t) = 0, (x, y) ∈ D (x, y) ∈ ∂D v (x, y, 0) = f (x, y) − uE (x, y) This is the Heat Problem with homogeneous PDE and BCs. We found the solution to this problem above in (43), (provided we use x0 = 1 = y0 in (43); don’t mix up the side length y0 in (43) with the location y0 of the hot spot in the homogeneous problem): ∞ � ∞ � 2 2 2 v (x, y, t) = Amn sin (mπx) sin (nπy) e−π (m +n )t m=1 n=1 where Amn = 4 � 1 0 � 1 0 (f (x, y) − uE (x, y)) sin (mπx) sin (nπy) dydx Thus we have found the full solution u (x, y, t). 10 Heat problem on a circle with inhomogeneous BC Consider the heat problem ut = ∇2 u, u (x, y, t) = g (x, y) , u (x, y, 0) = f (x, y) 21 (x, y) ∈ D (x, y) ∈ ∂D (65) where D = {(x, y) : x2 + y 2 ≤ 1} is a circle (disc) of radius 1. 10.1 Equilibrium solution and Laplace’s Eq. on a circle Ref: Guenther & Lee §8.1, 8.2, Myint-U & Debnath §8.4, 8.5 To solve the problem, we must first introduce the steady-state u = uE (x, y) which satisfies the PDE and BCs, ∇2 uE = 0, uE (x, y) = g (x, y) , (x, y) ∈ D, (x, y) ∈ ∂D. As before, switch to polar coordinates via x = r cos θ, y = r sin θ, wE (r, θ) = uE (x, y) for 0 ≤ r ≤ 1, −π ≤ θ < π. The problem for uE becomes � � ∂wE 1 ∂ 2 wE 1 ∂ r + 2 = 0, r ∂r ∂r r ∂θ2 w (r, −π) = w (r, π) , 0 ≤ r ≤ 1, −π ≤ θ < π wθ (r, −π) = wθ (r, π) , |w (0, θ)| < ∞ w (1, θ) = ĝ (θ) , (66) (67) (68) −π ≤ θ < π (69) where ĝ (θ) = g (x, y) for (x, y) ∈ ∂D and θ = arctan (y/x). We separate variables w (r, θ) = R (r) H (θ) and the PDE becomes r d R (r) dr � � dR d2 H 1 r =− 2 dr dθ H (θ) Since the l.h.s. depends only on r and the r.h.s. on θ, both must be equal to a constant µ, � � dR d2 H 1 r d r =− 2 =µ R (r) dr dr dθ H (θ) The problem for H (θ) is, as before, d2 H + µH (θ) = 0; dθ2 H (−π) = H (π) , 22 dH dH (−π) = (π) , dθ dθ with eigen-solutions Hm (θ) = am cos (mθ) + bm sin (mθ) , m = 0, 1, 2, ... The problem for R (r) is � � d dR d2 R dR 0=r r − m2 R = r 2 2 + r − m2 R dr dr dr dr Try R (r) = rα to obtain the auxiliary equation α (α − 1) + α − m2 = 0 whose solutions are α = ±m. Thus for each m, the solution is Rm (r) = c1 rm + c2 r−m . For m > 0, r−m blows up as r → 0. Our boundedness criterion (68) implies c2 = 0. Hence Rm (r) = cm rm where cm are constants to be found by imposing the BC (69). The separable solutions satisfying the PDE (66) and conditions (67) to (68), are wm (r, θ) = rm (Am cos (mθ) + Bm sin (mθ)) , m = 0, 1, 2, ... The full solution is the infinite sum of these over m, uE (x, y) = wE (r, θ) = ∞ � rm (Am cos (mθ) + Bm sin (mθ)) (70) m=0 We still need to find the Am , Bm . Imposing the BC (69) gives ∞ � (Am cos (mθ) + Bm sin (mθ)) = ĝ (θ) (71) m=0 The orthogonality relations are � π cos (mθ) sin (nθ) dθ = 0 −π � � � π� cos (mθ) cos (nθ) π, m = n dθ = , sin (mθ) sin (nθ) 0, m = � n −π (m > 0) . Multiplying (71) by sin nθ or cos nθ and applying these orthogonality relations gives � π 1 ĝ (θ) dθ A0 = 2π −π � 1 π Am = ĝ (θ) cos (mθ) dθ (72) π −π � 1 π Bm = ĝ (θ) sin (mθ) dθ π −π 23 y u=0 D 1 u=u0/(ST0 T0 x Figure 2: Setup for hot spot problem on circle. 10.1.1 Hot spot on boundary Suppose ĝ (θ) = � u0 θ0 +π 0 −π ≤ θ ≤ θ0 otherwise which models a hot spot on the boundary. The Fourier coefficients in (72) are thus A0 = u0 , 2π Bm = − Am = u0 sin (mθ0 ) mπ (θ0 + π) u0 (cos (mθ0 ) − (−1)m ) mπ (θ0 + π) Thus the steady-state solution is u0 � m uE = + r 2π m=1 ∞ 10.1.2 � � u0 (cos (mθ0 ) − (−1)m ) u0 sin (mθ0 ) cos (mθ) − sin (mθ) mπ (θ0 + π) mπ (θ0 + π) Interpretation The convergence of the infinite series is rapid if r ≪ 1. If r ≈ 1, many terms are required for accuracy. The center temperature (r = 0) at equilibrium (steady-state) is � π 1 u0 = ĝ (θ) dθ uE (0, 0) = wE (0, θ) = 2π 2π −π 24 i.e., the mean temperature of the circumference. This is a special case of the Mean Value Property of solutions to Laplace’s Equation ∇2 u = 0. We now consider plots of uE (x, y) for some interesting cases. We draw the level curves (isotherms) uE = const as solid lines. Recall from vector calculus that the gradient of uE , denoted by ∇uE , is perpendicular to the level curves. Recall also from the physics that the flux of heat is proportional to ∇uE . Thus heat flows along the lines parallel to ∇uE . Note that the heat flows even though the temperature is in steady-state. It is just that the temperature itself at any given point does not change. We call these lines the “heat flow lines” or the “orthogonal trajectories”, and draw these as dashed lines in the figure below. Note that lines of symmetry correspond to (heat) flow lines. To see this, let nl be the normal to a line of symmetry. Then the flux at a point on the line is ∇u · nl . Rotate the image about the line of symmetry. The arrow for the normal to the line of symmetry is now pointing in the opposite direction, i.e. −nl , and the flux is −∇u · nl . But since the solution is the same, the flux across the line must still be ∇u · nl . Thus ∇u · nl = −∇u · nl which implies ∇u·nl = 0. Thus there is no flux across lines of symmetry. Equivalently, ∇u is perpendicular to the normal to the lines of symmetry, and hence ∇u is parallel to the lines of symmetry. Thus the lines of symmetry are flow lines. Identifying the lines of symmetry help draw the level curves, which are perpendicular to the flow lines. Also, lines of symmetry can be thought of as an insulating boundary, since ∇u · nl = 0. (i) θ0 = 0. Then ∞ 2u0 � 2n−1 sin ((2n − 1) θ) u0 − 2 r uE = 2π π m=1 2n − 1 Use the BCs for the boundary. Note that the solution is symmetric with respect to the y-axis (i.e. even in x). The solution is discontinuous at {y = 0, x = ±1}, or {r = 1, θ = 0, π}. See plot. (ii) −π < θ0 < −π/2. The sum for uE is messy, so we use intuition. We start with the boundary conditions and use continuity in the interior of the plate to obtain a qualitative idea of the level curves and heat flow lines. See plot. (iii) θ0 → −π + (a heat spot). Again, use intuition to obtain a qualitative sketch of the level curves and heat flow lines. Note that the temperature at the hot point is infinite. See plot. The heat flux lines must go from a point of hot temperature to a point of low temperature. Since the temperature along boundary is zero except at 25 u=0 u=0 u u=u0/S u=0 u of u0 S T0 T 0 o S T0 S T0  Figure 3: Plots of steady-state temperature due to a hot segment on boundary. θ = −π, then all the heat flux lines must leave the hot spot and extend to various points on the boundary. 10.2 Solution to inhomogeneous heat problem on circle We now use the standard trick and solve the inhomogeneous Heat Problem (65), ut = ∇2 u, u (x, y, t) = g (x, y) , u (x, y, 0) = f (x, y) (x, y) ∈ D (x, y) ∈ ∂D (73) using the equilibrium solution uE . We define the transient part of the solution as v (x, y, t) = u (x, y, t) − uE (x, y) where u (x, y, t) is the solution to (73). The problem for v (x, y, t) is therefore vt = ∇2 v, v (x, y, t) = 0, (x, y) ∈ D (x, y) ∈ ∂D v (x, y, 0) = f (x, y) − uE (x, y) This is the Heat Problem with homogeneous PDE and BCs. We found the solution to this problem above in (61), v (x, y, t) = ∞ � ∞ � Jm (rjm,n ) (αmn cos (mθ) + βmn sin (mθ)) e−λmn t m=1 n=1 where αmn and βmn are found from orthogonality relations and f (x, y) − uE (x, y). Thus we have found the full solution u (x, y, t). 26 11 Mean Value Property for Laplace’s Eq. Ref: Guenther & Lee §8.4, Myint-U & Debnath §8.4 Theorem [Mean Value Property] Suppose v (x, y) satisfies Laplace’s equation in a 2D domain D, ∇2 v = 0, (x, y) ∈ D. (74) Then at any point (x0 , y0 ) in D, v equals the mean value of the temperature around any circle centered at (x0 , y0 ) and contained in D, � π 1 v (x0 , y0 ) = v (x0 + R cos θ, y0 + R sin θ) dθ. (75) 2π −π Note that the curve {(x0 + R cos θ, y0 + R sin θ) : −π ≤ θ < π} traces the circle of radius R centered at (x0 , y0 ). Proof: To prove the Mean Value Property, we first consider Laplace’s equation (74) on the unit circle centered at the origin (x, y) = (0, 0). We already solved this problem, above, when we solved for the steady-state temperature uE that took the value ĝ (θ) on the boundary. The solution is Eqs. (70) and (72). Setting r = 0 in (70) gives the center value � π 1 ĝ (θ) dθ (76) uE (0, 0) = A0 = 2π −π On the boundary of the circle (of radius 1), (x, y) = (cos θ, sin θ) and the BC implies that uE = ĝ (θ) on that boundary. Thus, ĝ (θ) = uE (cos θ, sin θ) and (76) becomes � π 1 uE (0, 0) = A0 = uE (cos θ, sin θ) dθ. (77) 2π −π We now consider the region � � B(x0 ,y0 ) (R) = (x, y) : (x − x0 )2 + (y − y0 )2 ≤ R2 , which is a disc of radius R centered at (x, y) = (x0 , y0 ). Since B(x0 ,y0 ) (R) ⊆ D, Laplace’s equation (74) holds on this disc. Thus ∇2 v = 0, (x, y) ∈ B(x0 ,y0 ) (R) . (78) We make the change of variable x̂ = x − x0 , R ŷ = y − y0 , R uE (x̂, ŷ) = v (x, y) (79) to map the circle B(x0 ,y0 ) (R) into the unit disc {(x̂, ŷ) : x̂2 + ŷ 2 ≤ 1}. Laplace’s equa­ tion (78) becomes � � ˆ 2 uE = 0, ∇ (x̂, ŷ) ∈ (x̂, ŷ) : x̂2 + ŷ 2 ≤ 1 27 ˆ 2 = (∂ 2 /∂x̂2 , ∂ 2 /∂ŷ 2 ). The solution is given by Eqs. (70) and (72), and we where ∇ found the center value above in Eq. (77). Reversing the change of variable (79) in Eq. (77) gives � π 1 v (x0 , y0 ) = v (x0 + R cos θ, y0 + R sin θ) dθ 2π −π as required. � For the heat equation, the Mean Value Property implies the equilibrium tempera­ ture at any point (x0 , y0 ) in D equals the mean value of the temperature around any disc centered at (x0 , y0 ) and contained in D. 12 Maximum Principle for Laplace’s Eq. Ref: Guenther & Lee §8.4, Myint-U & Debnath §8.2 Theorem [Maximum Principle] Suppose v (x, y) satisfies Laplace’s equation in a 2D domain D, ∇2 v = 0, (x, y) ∈ D. Then the function v takes its maximum and minimum on the boundary of D, ∂D. Proof: Let (x0 , y0 ) be an interior point in D where v takes it’s maximum. In particular, (x, y) is not on the boundary ∂D of D. The Mean Value Property implies that for any disc B(x0 ,y0 ) (R) of radius R > 0 centered at (x0 , y0 ), � π 1 v (x0 + R cos θ, y0 + R sin θ) dθ v (x0 , y0 ) = 2π −π = Average of v on boundary of circle But the average is always between the minimum and maximum. Thus v (x0 , y0 ) must be between the minimum and maximum value of v (x, y) on the boundary of the disc. But since v (x0 , y0 ) is the maximum of v on D, then the entire boundary of the circle must have value v (x0 , y0 ). Since R > 0 is arbitrary, this holds for all values in the disc B(x0 ,y0 ) (R). Keep increasing R until the disc hits the boundary of D. Then v (x0 , y0 ) is a value of v along the boundary of D. A similar argument holds for the minimum of v. � Corollary If v = 0 everywhere on the boundary, then v must be zero at everywhere in D. For the homogeneous heat equation (i.e. no sinks/sources), this implies the equi­ librium temperature attains its maximum/minimum on the boundary. 28 If the boundary is completely insulated, ∇v = 0 on ∂D, then the equilibrium temperature is constant. To see this, we apply result (26) to the case where v is the solution to Laplace’s equation with ∇v = 0 on ∂D, � � � � � � � � � 2 2� 0= v∇v · n̂dS = v∇ v + |∇v| dV = |∇v|2 dV ∂D D D Thus, by continuity, |∇v| = 0 and ∇v = 0 everywhere in D and hence v = const in D. 13 Eigenvalues on different domains In this section we revisit the Sturm-Liouville problem ∇2 v + λv = 0, v = 0, x∈D x ∈ ∂D and consider the effect of the shape of the domain on the eigenvalues. Definition Rayleigh Quotient ��� ∇v · ∇vdV � �D � R (v) = v 2 dV D (80) Theorem Given a domain D ⊆ R3 and any function v that is piecewise smooth on D, non-zero at some points on the interior of D, and zero on all of ∂D, then the smallest eigenvalue of the Laplacian on D satisfies λ ≤ R (v) and R (h) = λ if and only if h (x) is an eigen-solution of the Sturm Liouville problem on D. Sketch Proof: We use result (26) derived for any smooth function v defined on a volume V with closed smooth surface S. � � � � � � 2 � v∇v · n̂dS = v∇ v + ∇v · ∇v dV ∂D D In the statement of the theorem, we assumed that v = 0 on ∂D, and hence � � � � � � ∇v · ∇vdV = − v∇2 vdV D (81) D Let {φn } be an orthonormal basis of eigenfunctions on D, i.e. all the functions φn which satisfy ∇2 φn + λn φn = 0, φn = 0, 29 x∈D x ∈ ∂D and � � � φn φm dV = D � 1, m = n 0, m �= n We can expand v in the eigenfunctions, v (x) = ∞ � An φn (x) n=1 where the An are constants. Assuming we can differentiate the sum termwise, we have ∞ ∞ � � 2 2 ∇v= An ∇ φn = − λn An φn (82) n=1 n=1 ��� 2 The orthonormality property (i.e., the orthogonality property with φ dV = 1) D n implies � � � � � � ∞ � ∞ ∞ � � 2 v dV = An Am φn φm dV = A2n (83) D D n=1 m=1 n=1 and, from (82), � � � D 2 v∇ vdV = − ∞ � ∞ � λn An Am � � � D m=1 n=1 φn φm dV = − ∞ � λn A2n (84) n=1 Substituting (84) into (81) gives � � � D ∇v · ∇vdV = − � � � 2 v∇ vdV = D ∞ � λn A2n (85) n=1 Substituting (83) and (85) into (80) gives ��� �∞ 2 ∇v · ∇vdV n=1 λn An D� � � � = R (v) = ∞ 2 v 2 dV n=1 An D We assume the eigenfunctions are arranged in increasing order. In particular, λn ≥ λ1 . Thus �∞ �∞ 2 λ1 A2n n=1 An n=1 � = λ = λ1 . R (v) ≥ �∞ 1 ∞ 2 2 n=1 An n=1 An If v is an eigenfunction with eigenvalue λ1 , then R (v) = λ1 A21 = λ1 A21 If v is not an eigenfunction corresponding to λ1 , then there exists an n > 1 such that λn > λ1 and An = � 0, so that R (v) > λ1 . � 30 � and D in R2 satisfy Theorem If two domains D D � D � but D = � i.e., D ⊂ D � D, � , λ1 and then the smallest eigenvalues of the Sturm-Liouville problems on D and D λ̂1 , respectively, satisfy �1 < λ1 λ � that contains the sub-domain D is associated with a In other words, the domain D smaller eigenvalue. Proof: Note that the Sturm-Liouville problems are ∇2 v + λv = 0, (x, y) ∈ D v = 0, (x, y) ∈ ∂D �v� = 0, ∇2 v� + λ � (x, y) ∈ D � (x, y) ∈ ∂D v� = 0, Let v1 be the eigenfunction corresponding to λ1 on D. Then, as we have proven before, λ1 = R (v1 ) , (86) where R (v) is the Rayleigh Quotient. Extend the function v1 continuously from D � to obtain a function v̂1 on D � which satisfies to D � v1 (x, y) , (x, y) ∈ D v̂1 = � 0 (x, y) ∈ D, (x, y) ∈ /D The extension is continuous, since v1 is zero on the boundary of D. Applying the � and function v̂1 (which satisfies all the requirements previous theorem to the region D of the theorem) gives λ̂1 ≤ R (v̂1 ) . Equality happens only if v̂1 is the eigenfunction corresponding to λ̂1 . Useful fact [stated without proof]: the eigenfunction(s) corresponding to the � are nonzero in the interior of D �. smallest eigenvalue λ̂1 on D �, From the useful fact, v̂1 cannot be an eigenfunction corresponding to λ̂1 on D � (outside D). Thus, as the previous theorem states, since it is zero in the interior of D λ̂1 < R (v̂1 ) . 31 (87) � in the Rayleigh quotient reduce to Since v̂1 = 0 outside D, the integrals over D integrals over D, where v̂1 = v1 , and hence R (v̂1 ) = R (v1 ) . (88) Combining (86), (87), and (88) gives the result, λ̂1 < λ1 . � Example: Consider two regions, D1 is a rectangle of length x0 , height y0 and D2 is a circle of radius R. Recall that the smallest eigenvalue on the rectangle D1 is � � 1 1 2 + λ11 = π x20 y02 The smallest eigenvalue on the circle of radius 1 is λ01 = J02,1 where the first zero √ of the Bessel function J0 (s) of the first kind is J0,1 = 2.4048. Since λ multiplied r in the Bessel function, then�for a�circle of�� radius�R, we’d rescale by the change of √ λ r̂ where r̂ goes from 0 to 1. Thus variable r̂ = r/R, so that Jm λr = Jm R2 on the circle of radius R, the smallest eigenvalue is � �2 J0,1 λ01 = , J0,1 = 2.4048. R Suppose the rectangle is actually a square of side length 2R. Then � �2 4.934 J0,1 5.7831 π2 = , λ01 = = λ11 = 2 2 2R R R R2 Thus, λ11 < λ01 , which confirms the second theorem, since D2 ⊂ D1 , i.e., the circle is contained inside the square. Now consider the function � r �2 v (r) = 1 − R You can show that ∇v · ∇v = and R (v) = �� ∇v · ∇vdA �� = v 2 dA D2 D2 � dv dr �2 , � π � R � dv �2 rdrdθ dr −π 0 �π �R v 2 rdrdθ −π 0 = 6 > λ01 R2 This confirms the first theorem, since v (r) is smooth on D2 , v (R) = 0 (zero on the boundary of D2 ), and v is nonzero in the interior. 32 13.1 Faber-Kahn inequality Thinking about the heat problem on a 2D plate, what shape of plate will cool the slowest? It is a geometrical fact that of all shapes of equal area, the circle (disc) has the smallest circumference. Thus, on physical grounds, we expect the circle to cool the slowest. Faber and Kahn proved this in the 1920s. Faber-Kahn inequality For all domains D ⊂ R2 of equal area, the disc has the smallest first eigenvalue λ1 . √ Example. Consider the circle of radius 1 and the square of side length π. Then both the square and circle have the same area. The first eigenvalues for the square and circle are, respectively, � � 1 1 2 λ1SQ = π + = 2π = 6.28, λ1CIRC = (J0,1 )2 = 5.7831 π π and hence λ1SQ > λ1CIRC , as the Faber-Kahn inequality states. 14 Nodal lines Consider the Sturm-Liouville problem ∇2 v + λv = 0, v = 0, x∈D (89) x ∈ ∂D Nodal lines are the curves where the eigenfunctions of the Sturm-Liouville problem are zero. For the solution to the vibrating membrane problem, the normal modes unm (x, y, t) are zero on the nodal lines, for all time. These are like nodes on the 1D string. Here we consider the nodal lines for the square and the disc (circle). 14.1 Nodal lines for the square For the square, the eigenfunctions and eigenvalues are a special case of those we found for the rectangle, with side length x0 = y0 = a, � mπx � � mπy � � π2 � vmn (x, y) = sin sin , λmn = 2 m2 + n2 , n, m = 1, 2, 3, ... a a a The nodal lines are the lines on which vmn (x, y) = 0, and are x= ka , m y= la , n 1 ≤ k ≤ m − 1, 1≤l ≤n−1 for m, n ≥ 2. Note that v11 (x, y) has no nodal lines on the interior - it is only zero on boundary ∂D. Since λmn = λnm , then the function fnm = Avmn + Bvnm is also an 33 eigenfunction with eigenvalue λmn , for any constants A, B. The nodal lines for fnm can be quite interesting. Examples: we draw the nodal lines on the interior and also the lines around the boundary, where vnm = 0. (i) m = 1, n = 1. � πy � � πx � sin v11 = sin a a This is positive on the interior and zero on the boundary. Thus the nodal lines are simply the square boundary ∂D. (ii) m = 1, n = 2. � � � πx � 2πy v12 = sin sin a a The nodal lines are the boundary ∂D and the horizontal line y = a/2. (iii) m = 1, n = 3 � � � πx � 3πy sin v13 = sin a a The nodal lines are the boundary ∂D and the horizontal lines y = a/3, 2a/3. (iv) m = 3, n = 1 � � � πy � 3πx v31 = sin sin a a The nodal lines are the boundary ∂D and the vertical lines x = a/3, 2a/3. (v) Consider v13 − v31 . Since λ31 = λ13 = 10π 2 /a2 , this is a solution to SL problem (89) on D with λ = λ13 . To find the nodal lines, we use the identity � � sin 3θ = (sin θ) 3 − 4 sin2 θ to write � πy � � � πy �� � πx � sin 3 − 4 sin2 v13 = sin a a a � πx � � πy � � � πx �� v31 = sin sin 3 − 4 sin2 a a a � πy � � � πx � � πy �� � πx � sin sin2 − sin2 v13 − v31 = 4 sin a a a a The nodal lines are the boundary of the square, ∂D, and lines such that � � � � � � πx � � πy �� � � πx � � πy �� 2 πx 2 πy 0 = sin − sin = sin − sin sin + sin a a a a a a i.e., sin Thus � πx � a � πy � = sin ± a πy πx =± + kπ, a a 34 any integer k. Hence the nodal lines are y = ±x + ka for all integers k. We’re only concerned with the nodal lines that intersect the interior of the square plate: y = x, y = −x + a Thus the nodal lines of v13 − v31 are the sides and diagonals of the square. Let DT be the isosceles right triangle whose hypotenuse lies on the bottom hor­ izontal side of the square. The function v13 − v31 is zero on the boundary ∂DT , positive on the interior of DT , and is thus the eigenfunction corresponding to the first (smallest) eigenvalue λ = λ13 of the Sturm-Liouville problem (89) on the triangle DT . In this case, we found the eigenfunction without using separation of variables, which would have been complicated on the triangle. (vi) With v13 , v31 given above, adding gives � πy � � 3 � � � � ��� � πx � 2 πx 2 πy sin − sin + sin v13 + v31 = 4 sin a a 2 a a The nodal lines for v13 + v31 are thus the square boundary and the closed nodal line defined by � πx � � πy � 3 sin2 + sin2 = . a a 2 Let Dc be the region contained within this closed nodal line. The function − (v13 + v31 ) is zero on the boundary ∂Dc , positive on the interior of Dc , and thus − (v13 + v31 ) is the eigenfunction corresponding to the first (smallest) eigenvalue λ13 of the SturmLiouville problem (89) on Dc . (vii) Find the first eigenvalue on the right triangle � � √ DT 2 = 0 ≤ y ≤ 3x, 0 ≤ x ≤ 1 . Note that separation of variables is ugly, because you’d have to impose the BC �√ � X (x) Y 3x = 0 We proceed by placing the triangle inside a rectangle of horizontal and vertical side √ lengths 1 and 3, respectively. The sides of the rectangle coincide with the perpen­ dicular sides of the triangle. Thus, all eigenfunctions vmn for the rectangle are already zero on two sides of the triangle. However, any particular eigenfunction vmn will not be zero on the triangle’s hypotenuse, since all the nodal lines of vmn are horizontal or vertical. Thus we need to add multiple eigenfunctions. However, to satisfy the Sturm-Liouville problem, all the eigenfunctions must be associated with the same 35 √ eigenvalue. For the rectangle with side lengths 3 and 1, the eigenvalues are given by � 2 � � m π2 � 2 n2 2 λmn = π = + 3m + n2 1 3 3 We create a table of 3m2 + n2 and look for repeated values: n, m 1 2 3 4 5 1 4 7 12 19 28 2 13 16 21 28 37 3 4 28 31 36 43 The smallest repeated value of 3m2 + n2 is 28: λ31 = λ24 = λ15 = 28π 2 /3 The linear combination of the corresponding eigenfunctions is itself an eigenfunction of the SL problem (89) on the rectangle D, v = Av31 + Bv24 + Cv15 . We wish v to also be zero along the hypotenuse of DT 2 , and seek A, B, C which satisfy Av31 + Bv24 + Cv15 = 0 √ on y = 3x. If we can find such constants A, B, C then we have found an eigenfunc­ tion of the SL problem (89) on the triangle DT 2 . If this eigenfunction is positive on the interior of DT 2 then we have also found the smallest eigenvalue (in this case λ31 ). If we cannot find such A, B, C, then we seek the next largest repeated value of 3m2 + n2 and continue as above. If we can find (A, B, C) such that v = 0 on √ y = 3x, but v is zero on the interior of DT 2 , then λ31 is not the smallest eigenvalue and we must pursue another method of finding the smallest eigenvalue. The Rayleigh quotient can help us identify an upper bound. 14.2 Nodal lines for the disc (circle) For the disc of radius 1, we found the eigenfunctions and eigenvalues to be vmnS = Jm (rjm,n ) sin mθ, vmnC = Jm (rjm,n ) cos mθ 36 with 2 λmn = π 2 jm,n , n, m = 1, 2, 3, ... The nodal lines are the lines on which vmnC = 0 or vmnS = 0. Examples. (i) m = 0, n = 1. v01 = J0 (rj0,1 ) The nodal lines are the boundary of the disc (circle of radius 1). (ii) m = 0, n = 2. v02 = J0 (rj0,2 ) The nodal lines are two concentric circles, one of radius r = 1, the other or radius r = J0,1 /J0,2 < 1. (iii) m = 1 , n = 1 and sine. v11S = J0 (rj1,1 ) sin θ The nodal lines are the boundary (circle of radius 1) and the line θ = −π, 0, π (horizontal diameter). 15 Steady-state temperature in a 3D cylinder Suppose a 3D cylinder of radius a and height L has temperature u (r, θ, z, t). We assume the axis of the cylinder is on the z-axis and (r, θ, z) are cylindrical coordinates. Initially, the temperature is u (r, θ, z, 0). The ends are kept at a temperature of u = 0 and sides kept at u (a, θ, z, t) = g (θ, z). The steady-state temperature uE (r, θ, z) in the 3D cylinder is given by ∇2 uE = 0, −π ≤ θ < π, 0 ≤ r ≤ a, 0 ≤ z ≤ L, uE (r, θ, 0) = uE (r, θ, L) = 0, uE (a, θ, z) = g (θ, z) , −π ≤ θ < π, 0 ≤ r ≤ a, −π ≤ θ < π, 0 ≤ z ≤ L. In cylindrical coordinates (r, θ, z), the Laplacian operator becomes � � ∂v 1 ∂ 2v ∂ 2v 1 ∂ 2 r + 2 2+ 2 ∇v= r ∂r ∂r r ∂θ ∂z We separate variables as v (r, θ, z) = R (r) H (θ) Z (z) 37 (90) (91) so that (90) becomes 1 d rR (r) dr � � dR (r) 1 d2 H(θ) 1 d2 Z (z) r + 2 + = 0. dr r H(θ) dθ2 Z (z) dz 2 Rearranging gives 1 d rR (r) dr � � dR (r) 1 d2 H(θ) 1 d2 Z (z) r + 2 = − =λ dr r H(θ) dθ2 Z (z) dz 2 (92) where λ is constant since the l.h.s. depends only on r, θ while the middle depends only on z. The function Z (z) satisfies d2 Z + λZ = 0 dz 2 The BCs at z = 0, L imply 0 = u (r, θ, 0) = R (r) H (θ) Z (0) 0 = u (r, θ, L) = R (r) H (θ) Z (L) To obtain non-trivial solutions, we must have Z (0) = 0 = Z (L) . As we’ve shown many times before, the solution for Z (z) is, up to a multiplicative constant, � nπz � � nπ �2 Zn (z) = sin , λn = , n = 1, 2, 3, ... L L Multiplying Eq. (92) by r2 gives � � dR (r) 1 d2 H(θ) r d r − λn r 2 = − =µ (93) R (r) dr dr H(θ) dθ2 where µ is constant since the l.h.s. depends only on r and the middle only on θ. We assume g (θ, z) is smooth and 2π-periodic in θ. Hence H ′ (−π) = H ′ (π), H ′ (−π) = H ′ (π) We solved this problem above and found that µ = m2 for m = 0, 1, 2, ..., H0 (θ) = const and Hm (θ) = Am cos mθ + Bm sin mθ, m = 1, 2, 3, ... Multiplying (93) by R (r) gives � � � � d dR (r) r r − λn r2 + m2 R (r) = 0 dr dr 38 � √ � This is the Modified Bessel Equation with linearly independent solutions Im r λn � √ � and Km r λn called the modified Bessel functions of order m of the first and second kinds, respectively. Thus � � � � � � Rmn (r) = c1m Im r λn + c2m Km r λn Since the Im ’s are regular (bounded) at r = 0, while the Km ’s are singular (blow up), and since Rmn (r) must be bounded, we must have c2m = 0, or � nπr � � � � Rmn (r) = c1m Im r λn = c1m Im L Thus the general solution is, by combining constants, uE (r, θ, z) = ∞ � ∞ � m=0 n=1 Im � nπr � L sin � nπz � L [Amn cos (mθ) + Bm sin (mθ)] In theory, we can now impose the condition u (a, θ, z) = g (θ, z) and find Amn , Bmn using orthogonality of sin, cos. 39