19.5 Relaxation Methods for Boundary Value Problems 863 FACR Method The best way to solve equations of the form(19.4.28),including the constant coefficient problem(19.0.3),is a combination of Fourier analysis and cyclic reduction, the FACR method [3-6].If at the rth stage of CR we Fourier analyze the equations of the form (19.4.32)along y,that is,with respect to the suppressed vector index,we will have a tridiagonal system in the z-direction for each y-Fourier mode: 砖-2r+鸿+鸿+2=△2gk (19.4.35) 81 Here is the eigenvalue ofT()corresponding to the kth Fourier mode.For the equation (19.0.3),equation (19.4.5)shows that will involve terms like cos(2k/L)-2 raised to a power.Solve the tridiagonal systems for at the levels j=2,2 x 2r,4 x 2,...,J-2r.Fourier synthesize to get the y-values on these z-lines.Then fill in the intermediate z-lines as in the original CR algorithm. The trick is to choose the number of levels of CR so as to minimize the total number of arithmetic operations.One can show that for a typical case of a 128x128 mesh,the optimal level is r =2;asymptotically,r-log2(log2J). A rough estimate of running times for these algorithms for equation(19.0.3) is as follows:The FFT method(in both x and y)and the CR method are roughly comparable.FACR with r=0(that is,FFT in one dimension and solve the tridiagonal equations by the usual algorithm in the other dimension)gives about a factor of two gain in speed.The optimal FACR with r=2 gives another factor ss是0 of two gain in speed. CITED REFERENCES AND FURTHER READING: 6 Swartzrauber,P.N.1977,S/AM Review,vol.19,pp.490-501.[1] Buzbee,B.L,Golub,G.H.,and Nielson,C.W.1970,SIAM Journal on Numerical Analysis,vol.7, pp.627-656:see als0op.ct.vol.11,pp.753-763.2] Hockney,R.W.1965,Journal of the Association for Computing Machinery,vol.12,pp.95-113.[3] Hockney,R.W.1970,in Methods of Computational Physics,vol.9(New York:Academic Press). pp.135-211.[4] Hockney,R.W.,and Eastwood,J.W.1981.Computer Simulation Using Particles (New York: Recipes McGraw-Hill),Chapter 6.[5] Numerical Recipes 10621 43106 Temperton,C.1980,Journal of Computationa/Physics,vol.34,pp.314-329.[6] (outside North 19.5 Relaxation Methods for Boundary Value Problems As we mentioned in $19.0,relaxation methods involve splitting the sparse matrix that arises from finite differencing and then iterating until a solution is found There is another way of thinking about relaxation methods that is somewhat more physical.Suppose we wish to solve the elliptic equation Cu=p (19.5.1)19.5 Relaxation Methods for Boundary Value Problems 863 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). FACR Method The best way to solve equations of the form (19.4.28), including the constant coefficient problem (19.0.3),is a combination of Fourier analysis and cyclic reduction, the FACR method [3-6]. If at the rth stage of CR we Fourier analyze the equations of the form (19.4.32) along y, that is, with respect to the suppressed vector index, we will have a tridiagonal system in the x-direction for each y-Fourier mode: uk j−2r + λ(r) k uk j + uk j+2r = ∆2g(r)k j (19.4.35) Here λ(r) k is the eigenvalue of T(r) corresponding to the kth Fourier mode. For the equation (19.0.3), equation (19.4.5) shows that λ(r) k will involve terms like cos(2πk/L) − 2 raised to a power. Solve the tridiagonal systems for uk j at the levels j = 2r, 2 × 2r, 4 × 2r, ..., J − 2r. Fourier synthesize to get the y-values on these x-lines. Then fill in the intermediate x-lines as in the original CR algorithm. The trick is to choose the number of levels of CR so as to minimize the total number of arithmetic operations. One can show that for a typical case of a 128×128 mesh, the optimal level is r = 2; asymptotically, r → log2(log2 J). A rough estimate of running times for these algorithms for equation (19.0.3) is as follows: The FFT method (in both x and y) and the CR method are roughly comparable. FACR with r = 0 (that is, FFT in one dimension and solve the tridiagonal equations by the usual algorithm in the other dimension) gives about a factor of two gain in speed. The optimal FACR with r = 2 gives another factor of two gain in speed. CITED REFERENCES AND FURTHER READING: Swartzrauber, P.N. 1977, SIAM Review, vol. 19, pp. 490–501. [1] Buzbee, B.L, Golub, G.H., and Nielson, C.W. 1970, SIAM Journal on Numerical Analysis, vol. 7, pp. 627–656; see also op. cit. vol. 11, pp. 753–763. [2] Hockney, R.W. 1965, Journal of the Association for Computing Machinery, vol. 12, pp. 95–113. [3] Hockney, R.W. 1970, in Methods of Computational Physics, vol. 9 (New York: Academic Press), pp. 135–211. [4] Hockney, R.W., and Eastwood, J.W. 1981, Computer Simulation Using Particles (New York: McGraw-Hill), Chapter 6. [5] Temperton, C. 1980, Journal of Computational Physics, vol. 34, pp. 314–329. [6] 19.5 Relaxation Methods for Boundary Value Problems As we mentioned in §19.0, relaxation methods involve splitting the sparse matrix that arises from finite differencing and then iterating until a solution is found. There is another way of thinking about relaxation methods that is somewhat more physical. Suppose we wish to solve the elliptic equation Lu = ρ (19.5.1)