858 Chapter 19.Partial Differential Equations Two important techniques lead to "rapid"solution of equation (19.4.1)when the sparse matrix is of certain frequently occurring forms.The Fourier transform method is directly applicable when the equations have coefficients that are constant in space.The cyclic reduction method is somewhat more general;its applicability is related to the question of whether the equations are separable(in the sense of "separation of variables").Both methods require the boundaries to coincide with the coordinate lines.Finally,for some problems,there is a powerful combination of these two methods called FACR (Fourier Analysis and Cyclic Reduction).We now consider each method in turn,using equation (19.0.3),with finite-difference representation(19.0.6),as a model example.Generally speaking,the methods in this section are faster,when they apply,than the simpler relaxation methods discussed in $19.5;but they are not necessarily faster than the more complicated multigrid methods discussed in 819.6. Fourier Transform Method RECIPES I The discrete inverse Fourier transform in both z and y is 1 J-1L-1 ujl Umne-2xijm/Je-2xiln/L (19.4.2) 。 Press. n=0n=0 This can be computed using the FFT independently in each dimension,or else all at once via the routine fourn of $12.4 or the routine rlft3 of $12.5.Similarly, IENTIFIC J-1L-1 Pmne-2xijm/Je-2xiln/L 19.4.3) JL n=0n=0 If we substitute expressions (19.4.2)and (19.4.3)in our model problem (19.0.6), we find Numerica 10621 立nn(e2mim/W+e-2mim/J+e2min/L+e-2min/L-4）=mn△2 (19.4.4) 431 or (outside Recipes pmn△2 (19.4.5) 2πm +cos 20一2 L North Thus the strategy for solving equation(19.0.6)by FFT techniques is: .Compute Pmn as the Fourier transform J-1L-1 Pmn= pil e2ximj/Je2xinl/L (19.4.6) j=0=0 .Compute timn from equation (19.4.5).858 Chapter 19. Partial Differential Equations Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). Two important techniques lead to “rapid” solution of equation (19.4.1) when the sparse matrix is of certain frequently occurring forms. The Fourier transform method is directly applicable when the equations have coefficients that are constant in space. The cyclic reduction method is somewhat more general; its applicability is related to the question of whether the equations are separable (in the sense of “separation of variables”). Both methods require the boundaries to coincide with the coordinate lines. Finally, for some problems, there is a powerful combination of these two methods called FACR (Fourier Analysis and Cyclic Reduction). We now consider each method in turn, using equation (19.0.3), with finite-difference representation (19.0.6), as a model example. Generally speaking, the methods in this section are faster, when they apply, than the simpler relaxation methods discussed in §19.5; but they are not necessarily faster than the more complicated multigrid methods discussed in §19.6. Fourier Transform Method The discrete inverse Fourier transform in both x and y is ujl = 1 JL J
−1 m=0 L
−1 n=0 umne−2πijm/J e−2πiln/L (19.4.2) This can be computed using the FFT independently in each dimension, or else all at once via the routine fourn of §12.4 or the routine rlft3 of §12.5. Similarly, ρjl = 1 JL J
−1 m=0 L
−1 n=0 ρmne−2πijm/J e−2πiln/L (19.4.3) If we substitute expressions (19.4.2) and (19.4.3) in our model problem (19.0.6), we find umn e2πim/J + e−2πim/J + e2πin/L + e−2πin/L − 4  = ρmn∆2 (19.4.4) or umn = ρmn∆2 2  cos 2πm J + cos 2πn L − 2  (19.4.5) Thus the strategy for solving equation (19.0.6) by FFT techniques is: • Compute ρmn as the Fourier transform ρmn = J
−1 j=0 L
−1 l=0 ρjl e2πimj/J e2πinl/L (19.4.6) • Compute umn from equation (19.4.5).