- Email: [email protected]
The Poisson equation on the unit disk: a multigrid solver using polar coordinates
The AM
Transmitted by U. Trottenherg
A mukigri~ S&S for the Poison equation on the unit disk is presented. It is shown that the problems due to the use of polar coordinates, namely the singularity at the origin and the an&otropy of the operator, constitute no obstade to an efficient realization of a multigid program. Our numerical results attest to the efficiency of the solver, also in comparison with well-known packages.
INTRODUCTION The aim of this work was to develop an efficient mukigrid solver for the Poisson equation on the unit disk, with the use of polar coordinates and Dirichlet boundary conditions. The use of polar coordinates allows a simple treatment of the boundary conditions, even for those of Neumann type. On the other han& it results in an anisotropic operator, which also has a sing&&y at the origin. Our results show that it is possible to overcome thez difficulties by using line and local relaxation, thus preserving the multigrid efficiency of the solver. We have ccmpared the program developed with the d ~~01, which represent different weMnown packages MGol, msmm, approaches to solving this problem. The HSHPAKprogram [l] is a dlirect *On leave from the University of Sac?Pa&
Brazil.
APPLIED MATHEMATKS AND COMPuTATlON 25: 123-135 (1988) @ Elsevier Science Publishing Co., Inc., 1988 52 Vanderbilt Ave., New ‘York, NY 10017
123
SAUL0 R. M. BARROS
be competitive, being at least as efficient as the themofspecialcoordinates,in*ted~ ixAe~ofmuhigridsoaverS. we refer to Stiik and Trot&&erg [7&
-&a=
iId-
f
*(
f,
t)
((X,y)(X~+ye,l},
= (l-r)cosf
1
(1-r)sint'
maps the zxstaqle R:= ((t,r)(Ogt<2r, Tbis e theMifdisL.~~lmedoperatoris -p:=-
A rectanguhr
-1
1
l-r
[ (1-t)
p -+ar ate
p I
(1 -‘)z
a
Ogr
. )I
grid for the rectangie R is then constructed as fon0Ws:
c := ((ti,rj)lf=ih,, the mesh fj = jh,, tit& O
and ht = 2w/N,. The numbers Nr and Nt are chosen, for structural reasons, as powers of two.
Mdtigrid for Polar Coordinates
125
The operator-A* is discretized through centi differences at a point P = &, fj) = (ih,, jh,), resulting in the following fivepoint star:
-9p(l--
.
jh;:)
1 hF(l-
-1
jh,)
l-
-1
jh,
1- jhr
-#(l-
jhr+$)
.
where 9 = h,/h, and 2 is the sum of the four neighboring coefficients. The problem is solved on the rectangular grid G. Coarser grids with t&e same structure are used in the multigrid method. On the line with r coordinate equal to 0, we have DiricbIet boundary conditions, and on the ones with t coordinate equal to 0 and 2~, periodic boundary conditions. The line with 7 coordinate equal to 1 mrresponds to the origin, where the transformed operator is not defined. We discretize the equation at this point using the usual fivqoint star for the Laphce operator. More precisely, we use an average of rotations of this star, in order that aMpoints on the inner circle (which Correspond to points with r coordinate 1 - h, on the rec+&ugk) m involved in the discretization of the equation at the origin. In all grids we -have iL; 3 4, in order to permit this discIetization. Throughout this work we will occasionaBy not dissguish between points on the disk and the mrresponding ones on the rectangle. This simplifies the writing without kading to mnfusion. For example, we will say the inner circle, instead of saying the line with r coordinate 1 - hr. The anisotropy of the transformed Laplace operator depends on 9. In some cases, the “strong mnnections” of the operator Eieon circle lines, while in others they lie either on circle or on radial lines, depending on the m@om of the domain. This behavior plays an important role in the mnstruction of the solver.
2.
THE CONSTlNJCTlON OF THEZMULTIGRID SOLVER
(a) Stnrcture of Gridsad
@emtom
grid is defined as G := ((ti, r$I tt = ih,, r = ,&, *th 0 G i < N,, 06 j
SAUL0 R. hf. RARROS
to use line rekuation in amnection with standard
Eclfientiaifw~~~theaVprapriaterelaxationbtocarry~tthe locgllFourier analyst (see, e.g., [2])-The analysis which we present here was -ted also in the work of J. linden [3]. (In that work the P&son ~~inpolar~is~~~on~annuhts.~~~~~nonfhe full disk is then treated through a composite mesh discretizatio~~ methti.) The equation on the azmhs with internal radius R, and external radius Rs < 1 is &a cazsidsed here (z&c with Dirichlet bcmdary conditions). ‘JL’he case R,= 0 and R,=P represents the full cbk. The ,rmoothing factors
Multi@
fm Poiiu coordiruttes
127
obtaind are the following: Zebm row relarati
(circle relaxation):
km=
if
R,=O.
Zebm column relazQtson(radial daxati011):
In the form&s3 q = h,/h, = 2laN,/N, = 29T(z”rWM9.These smoothing factors correspond to the case in which two smoothing steps are performed (we perform one step befole and one after the aarseqid correction). For the full disk, we obtain
kR=
0.34
if
M,aM,+3,
0.51
if
M,=M,+2,
if
M,
. ~9.82
,LL~*~ = P (independently of the size of q). This shows that if M, 3 M, +3, we a~ exp& good results using row (circle) relaxation, while performing only column (radial) relaxation never leads to good convergence factors. The reason for that is that the radial relaxation is never efficient on relaxing in a neighborhood of the origix3,because in this part of the domain the “strong connections” are in circle lines (only the size of this neighborhood depends on q, not the existence of it). On the other hand, when M, < M, + 2, circle relaxation is no longer efficient gear the boundary, because the “strong connections” lie in r~&al lines in this part of the domain. It is then necessary to use alternating-direction relaxation in order to obtain good factors, bcause the anisotropic b&avior of the operator is not the same in the whole domain. aal
J
SAUU) R M. BABROS
TABLE
wt*WI
rE!batb
(1,4)
1
FACI'OBSFOB x=i PBOBLBlU ON Convergence factor
~YooIuPulsDco~
Dim. = 0.8
ANGUS’
0.9
0.99
0.999
0.04
0.06
0.17 0.097
0.20 0.12
-
0.077
0.27 0.051
0.41 0.16
Row
(73 Ah. zebra
AN
“Upper entry: V cycles (1,l); lower entry: W cycles (1, 1).
iUuliigid fm Polw Comdhates
129
injection in this case produces slightly worse results. On the other hand, the use of the average when AAf,< M, +2 even causes divergence in some cases.) 3.
TESTING AND IMPROVING THE SOLVER THROUGH LOGAL RELAXATION
The firstresults concerning the equation on the full disk are displayed in e 2- The numerical tests presented here, as well as the other ones in this work, were carried out in an IBM l 3983B. We can see that with the use of Vcycles, factors of about 0.2 are obtained. This is not fully satisfactory for such a problem. The use of W cycles changes the factors to a reasonable 0.12. We should however, consider the additional amount of work (about 3096) in order to achieve these factors with W cycles. For each relaxation sweep line relaxations are performed, at least in one directiofi. Each relaxation of a radial line includes the solution of a tridiagonaI system, which differs from line to line, just by the right-hand side. It is thus necessary to decompose the matrix only once per sweep. The circle relaxat%n involves linear systems that are no longer tridiagona& since the corresponding matrix has two more nonzero entries, due to the periodic boundary conditions. The systems change from line to line, whicli makes this relaxation very expensive (these systems are solved using a method suggested in Pereyra [4]). Observing the error between the discrete sohrtion and the iterations from our method, we have noted a great dominance of the error at the origin. Starting with approximations with very small errors near the origin led to initial sweeps with very good conveq~nce rates. These rates were no longer good when the errors near the ori@ shed the same order of the global error. The same behavior was observed with W cycles. These observations and some comments from Brandt [2] guided us to try ebcQZ r&x&m. The local relaxation was done only at the inner circle and at the origin. We started each smoothing step by relaxing the inner circle line and then the TABLE 2 NUMERICALLY COMPUTED CONVEEENCE
FAC7l’ORS’
Convergence factor M,
M,=4
5
6
7
7
6.167 0.110
am 0.116
0.200 0.m
0.201 0.116
8
0.176 0.114
0.190 0.113
0.216 0.115
0.214 0.115
“Upper entry: V cycles (1,l); lower entry: W CYcleS(l,l)-
SAUL0 R. M. BARROS
amountofegbiawork
siuceonlyonelineofthegp%dis betterresultsTheusea# Vcycles weseeinTable3.(!SeeahoTables4and5
INam
PkshmhgGy?al~witbWcydesalsoimprovedtbe
factozs IF-
example, with (I&, l&)=(7,7) we obtained a ourgpalwasmainlytoachieve
49 THEFINALSOLVER
The use of local re#amth prtduad very goad convergencefacto=, but the~~~in~tbese~was,iBsomecases,stinlarge. tried to reduce the amount of work, CQnvergeNex%te!G. If A+$&na, +3, Qnly ~~isperformedandtheworlrdoneisnotso~e.Inthiiscase we don’t change anything. If M&W,+2, alternating zebra is used and iixxe work is retpi& From the smoothing analysis, we see that, in these cases, the circle rehaticm is efficient on relaxinglines near the origin,but not near the IIUIJII~~IY, when the radial relaxationdoes the work. This analysis suggests that circle relaxation is important only for circle lines inside the
Multigridfor Polar Coor&nates
131 TABLE4
v
CYCLES (&I)
(IMPROVED VERSION)
Convergencefactor [CPU time per cycle (set)] M,
M,-4
5
6
7
7
0.058 0.081
0.088 0.18
0.078 0.27
0.080 0.48
8
0.060 0.17
0.065 0.28
0.062 0.67
0.082 1.10
circle of radius P = l/9. On the other hand, not doiug circle relaxation outside this circle saves considerable work (depending on 9), since this relaxation is expensive. We have then tried a relaxation, composed of the steps: (1) localrelaxation, (2) circle relaxation inside the circle of radius l/9, (3) radialrelaxation (on the full disk). With ?his new relaxation strategy, we improved the efficiency of the soher, producing only small changes in the convergence factors. ‘I’abks 4 and 5 present these new factors and times. In order to solve the differential equation we have used the multigrid cyc!es in the fu&m&igrid mode (see [2], [7J). It is a very efficient way to obtain a solution with the discretization accuracy (rather than using the cycles in an iterative manner). Results concerning the application of our solver in the full-mu&rid mode are presented in the next section.
TABLE5 w
CYCLES (&I)
(IMPItOVED VERSION)*
Convergencefactor [CPU time per cycle (se@]
n-i,
M,=4
5
6
7
7
O.PlO 0.077
0.116 0.23
0.116 0.37
0.116 0*67
8
0.114 0.18
0.113 0.32
0.115 0.92
0.115 1.54
“Without local relaxation.
SAUL0 R. M. BARROS
132
of our program (see Section 4) ~~01, and FX~HPAIL We have Vcydes with one rehxation before and one after The bfGol program was performed in the same ,
systemobtainedthrou&*ediscretizationwi*~~ same system solved by our “geometric” multigrid In the tabk we refer to our program as MGRXAR. Convergence andCPUtimesareshowninTable6. I& Table 7 we present results concerning the x&ion of the Poisson eqqt& with m sohtion g=sin(3x+y). Our program, as well n&l, was~inthefnlE-multigrid~Wepresentalsoresults~ve~ tbe~obwcyclies(~~lacal~n)inaurprogram(~~to~ ~WrThe~~~~isadirect~ffforthePoissanequation and uses polar coo&m&s Since ~~01 is not suitable for the fuknukigrid mode and since its use as an iterative solver would got produce competitive xesd& it is not inin this comparison. The results furnished here the error in the uniform norm, in relation to the ~<#ses:~ecaseM~=~~,inw~~wehave thesamenlmnberOa~~ineachdirection,andthecaseM,=M,+2,in which the mesh s&q in polar coordinates, are approximately the same in with the FI~AK program are the exact discretization errors in polar ccnrdinates (up to roundoff error). Through them we see that
TAB=6 COMPARISON OF CONVERGENCE FACTORS
Convergencefactor
[CPUtimepercycle(sec)] Program (&J&)=(7,4) (7S) (796)(797)(696)(6J) 0.16 0.16 0.1150.17 0.16 AM01 0.36 0.75 0.20 0.37 0.68 0.17 0.16 MGOl 0.40 0.12 MGPOLAR
0.0780.0890.0840.090 0.18 0.27 0.48 0.12 0.24 0.08ii 0.058 0.
in unifolm “In seeon&.
“Error
MGFOLARW
BESULTS FORTHEEQUATION WlTIi -AK
EXACT SOLUTION g =
norm.
lAm 4.24 033 1.25 5.21
0.403-a3 0.101 -03 Oas-03 0.115-03 0287-04
1.26 5.34 0.39 1.60 6.79
OAOS -03 0.101 -03 0.458-03 0.115-03 0.287 -04
1.21 5.22 0:31 1.27 5.60
tim& 0.42 1.46 5.38 -
CPU
MGol
0.5uO-03 0.132-03 0.329-04 -
Error”
sin(3X + $f)
CPUtim& Error” CPUtiId Error” Error” cPutim& O-156 - 02 0.26 0.159 - 02 0.32 0.162-02 0.29
0386 -03 OS63 -04 0.491-0s 0.3.25 -03 (991) 0.313-04
(7,7) (WI (7,5) (680
04s WJ fs,s)
MGFOLAR
COIUP-
TABLE 7
E
SAUlD R M. BARROS
REFERENGES Adams, P. Swantn#rbes and R Sweet, FISEF+E A packageof fortran subprogpms iix the sbhion obsepuabbellipticpartialdifferential equations (Vexsion3, June1919),NationalCenterfbr Atmos~ericResearch Boulder,CO 80307. A. Bran&, Mu&grid Techiques: 1984 Guide with Applications to Fluid Dynamics, GMIMtudie No. 85, GMD, St. k-tin, 1984. JoLhden, M-m fiir die Poisso&lei&mg in Kreisund Ftinggebiet under varwendnng beak Ibordinaten, Diplomrbeit, Inst. f&r Angewandte Maths , univ. Bonn, 1981.
j.
Multi@ fm PO&Wcoordinates
135
4 V. Pereym, Iteaateddefened correctionsfor no&near bamdary value problems, hrarma. M&. l&111-125 ( 5 J.R~pndK.~~o~~~offinitedifferenceandfiniteelement eqgathsbyalgeb&mukgrid(AMG),in MulrigridMethodi~jbrZ~and ZX~Equatfona(D. J.PaddonandH.Ho~Eds.),Inst.of~cs and its AppkWns Gakence Ser.,New Sex.,No. 3, clarendon Press, Oxford, 1985, m 169-al2 6 K. St&n, A multi-@d pqgam to solve AU-c(x,y)U=f(x,y) (on Q), U= g(x, II)(on Ws 0~1NoxhreamgularBounded Domains 62,IMA Report!!E!ma& GMD, Si. Am 198% 7 K. Stiiben and U. ‘I@Mu&rid metbods: Fundamental algorithms, dpaabtempna;lysisd~,in~NotesinMcrdzemotics,vol. 960(W. Haddwsdhand U. Trcaakg, Eds.), Springer,1982, pp. l-176.
Transmitted by U. Trottenherg
A mukigri~ S&S for the Poison equation on the unit disk is presented. It is shown that the problems due to the use of polar coordinates, namely the singularity at the origin and the an&otropy of the operator, constitute no obstade to an efficient realization of a multigid program. Our numerical results attest to the efficiency of the solver, also in comparison with well-known packages.
INTRODUCTION The aim of this work was to develop an efficient mukigrid solver for the Poisson equation on the unit disk, with the use of polar coordinates and Dirichlet boundary conditions. The use of polar coordinates allows a simple treatment of the boundary conditions, even for those of Neumann type. On the other han& it results in an anisotropic operator, which also has a sing&&y at the origin. Our results show that it is possible to overcome thez difficulties by using line and local relaxation, thus preserving the multigrid efficiency of the solver. We have ccmpared the program developed with the d ~~01, which represent different weMnown packages MGol, msmm, approaches to solving this problem. The HSHPAKprogram [l] is a dlirect *On leave from the University of Sac?Pa&
Brazil.
APPLIED MATHEMATKS AND COMPuTATlON 25: 123-135 (1988) @ Elsevier Science Publishing Co., Inc., 1988 52 Vanderbilt Ave., New ‘York, NY 10017
123
SAUL0 R. M. BARROS
be competitive, being at least as efficient as the themofspecialcoordinates,in*ted~ ixAe~ofmuhigridsoaverS. we refer to Stiik and Trot&&erg [7&
-&a=
iId-
f
*(
f,
t)
((X,y)(X~+ye,l},
= (l-r)cosf
1
(1-r)sint'
maps the zxstaqle R:= ((t,r)(Ogt<2r, Tbis e theMifdisL.~~lmedoperatoris -p:=-
A rectanguhr
-1
1
l-r
[ (1-t)
p -+ar ate
p I
(1 -‘)z
a
Ogr
. )I
grid for the rectangie R is then constructed as fon0Ws:
c := ((ti,rj)lf=ih,, the mesh fj = jh,, tit& O
and ht = 2w/N,. The numbers Nr and Nt are chosen, for structural reasons, as powers of two.
Mdtigrid for Polar Coordinates
125
The operator-A* is discretized through centi differences at a point P = &, fj) = (ih,, jh,), resulting in the following fivepoint star:
-9p(l--
.
jh;:)
1 hF(l-
-1
jh,)
l-
-1
jh,
1- jhr
-#(l-
jhr+$)
.
where 9 = h,/h, and 2 is the sum of the four neighboring coefficients. The problem is solved on the rectangular grid G. Coarser grids with t&e same structure are used in the multigrid method. On the line with r coordinate equal to 0, we have DiricbIet boundary conditions, and on the ones with t coordinate equal to 0 and 2~, periodic boundary conditions. The line with 7 coordinate equal to 1 mrresponds to the origin, where the transformed operator is not defined. We discretize the equation at this point using the usual fivqoint star for the Laphce operator. More precisely, we use an average of rotations of this star, in order that aMpoints on the inner circle (which Correspond to points with r coordinate 1 - h, on the rec+&ugk) m involved in the discretization of the equation at the origin. In all grids we -have iL; 3 4, in order to permit this discIetization. Throughout this work we will occasionaBy not dissguish between points on the disk and the mrresponding ones on the rectangle. This simplifies the writing without kading to mnfusion. For example, we will say the inner circle, instead of saying the line with r coordinate 1 - hr. The anisotropy of the transformed Laplace operator depends on 9. In some cases, the “strong mnnections” of the operator Eieon circle lines, while in others they lie either on circle or on radial lines, depending on the m@om of the domain. This behavior plays an important role in the mnstruction of the solver.
2.
THE CONSTlNJCTlON OF THEZMULTIGRID SOLVER
(a) Stnrcture of Gridsad
@emtom
grid is defined as G := ((ti, r$I tt = ih,, r = ,&, *th 0 G i < N,, 06 j
SAUL0 R. hf. RARROS
to use line rekuation in amnection with standard
Eclfientiaifw~~~theaVprapriaterelaxationbtocarry~tthe locgllFourier analyst (see, e.g., [2])-The analysis which we present here was -ted also in the work of J. linden [3]. (In that work the P&son ~~inpolar~is~~~on~annuhts.~~~~~nonfhe full disk is then treated through a composite mesh discretizatio~~ methti.) The equation on the azmhs with internal radius R, and external radius Rs < 1 is &a cazsidsed here (z&c with Dirichlet bcmdary conditions). ‘JL’he case R,= 0 and R,=P represents the full cbk. The ,rmoothing factors
Multi@
fm Poiiu coordiruttes
127
obtaind are the following: Zebm row relarati
(circle relaxation):
km=
if
R,=O.
Zebm column relazQtson(radial daxati011):
In the form&s3 q = h,/h, = 2laN,/N, = 29T(z”rWM9.These smoothing factors correspond to the case in which two smoothing steps are performed (we perform one step befole and one after the aarseqid correction). For the full disk, we obtain
kR=
0.34
if
M,aM,+3,
0.51
if
M,=M,+2,
if
M,
. ~9.82
,LL~*~ = P (independently of the size of q). This shows that if M, 3 M, +3, we a~ exp& good results using row (circle) relaxation, while performing only column (radial) relaxation never leads to good convergence factors. The reason for that is that the radial relaxation is never efficient on relaxing in a neighborhood of the origix3,because in this part of the domain the “strong connections” are in circle lines (only the size of this neighborhood depends on q, not the existence of it). On the other hand, when M, < M, + 2, circle relaxation is no longer efficient gear the boundary, because the “strong connections” lie in r~&al lines in this part of the domain. It is then necessary to use alternating-direction relaxation in order to obtain good factors, bcause the anisotropic b&avior of the operator is not the same in the whole domain. aal
J
SAUU) R M. BABROS
TABLE
wt*WI
rE!batb
(1,4)
1
FACI'OBSFOB x=i PBOBLBlU ON Convergence factor
~YooIuPulsDco~
Dim. = 0.8
ANGUS’
0.9
0.99
0.999
0.04
0.06
0.17 0.097
0.20 0.12
-
0.077
0.27 0.051
0.41 0.16
Row
(73 Ah. zebra
AN
“Upper entry: V cycles (1,l); lower entry: W cycles (1, 1).
iUuliigid fm Polw Comdhates
129
injection in this case produces slightly worse results. On the other hand, the use of the average when AAf,< M, +2 even causes divergence in some cases.) 3.
TESTING AND IMPROVING THE SOLVER THROUGH LOGAL RELAXATION
The firstresults concerning the equation on the full disk are displayed in e 2- The numerical tests presented here, as well as the other ones in this work, were carried out in an IBM l 3983B. We can see that with the use of Vcycles, factors of about 0.2 are obtained. This is not fully satisfactory for such a problem. The use of W cycles changes the factors to a reasonable 0.12. We should however, consider the additional amount of work (about 3096) in order to achieve these factors with W cycles. For each relaxation sweep line relaxations are performed, at least in one directiofi. Each relaxation of a radial line includes the solution of a tridiagonaI system, which differs from line to line, just by the right-hand side. It is thus necessary to decompose the matrix only once per sweep. The circle relaxat%n involves linear systems that are no longer tridiagona& since the corresponding matrix has two more nonzero entries, due to the periodic boundary conditions. The systems change from line to line, whicli makes this relaxation very expensive (these systems are solved using a method suggested in Pereyra [4]). Observing the error between the discrete sohrtion and the iterations from our method, we have noted a great dominance of the error at the origin. Starting with approximations with very small errors near the origin led to initial sweeps with very good conveq~nce rates. These rates were no longer good when the errors near the ori@ shed the same order of the global error. The same behavior was observed with W cycles. These observations and some comments from Brandt [2] guided us to try ebcQZ r&x&m. The local relaxation was done only at the inner circle and at the origin. We started each smoothing step by relaxing the inner circle line and then the TABLE 2 NUMERICALLY COMPUTED CONVEEENCE
FAC7l’ORS’
Convergence factor M,
M,=4
5
6
7
7
6.167 0.110
am 0.116
0.200 0.m
0.201 0.116
8
0.176 0.114
0.190 0.113
0.216 0.115
0.214 0.115
“Upper entry: V cycles (1,l); lower entry: W CYcleS(l,l)-
SAUL0 R. M. BARROS
amountofegbiawork
siuceonlyonelineofthegp%dis betterresultsTheusea# Vcycles weseeinTable3.(!SeeahoTables4and5
INam
PkshmhgGy?al~witbWcydesalsoimprovedtbe
factozs IF-
example, with (I&, l&)=(7,7) we obtained a ourgpalwasmainlytoachieve
49 THEFINALSOLVER
The use of local re#amth prtduad very goad convergencefacto=, but the~~~in~tbese~was,iBsomecases,stinlarge. tried to reduce the amount of work, CQnvergeNex%te!G. If A+$&na, +3, Qnly ~~isperformedandtheworlrdoneisnotso~e.Inthiiscase we don’t change anything. If M&W,+2, alternating zebra is used and iixxe work is retpi& From the smoothing analysis, we see that, in these cases, the circle rehaticm is efficient on relaxinglines near the origin,but not near the IIUIJII~~IY, when the radial relaxationdoes the work. This analysis suggests that circle relaxation is important only for circle lines inside the
Multigridfor Polar Coor&nates
131 TABLE4
v
CYCLES (&I)
(IMPROVED VERSION)
Convergencefactor [CPU time per cycle (set)] M,
M,-4
5
6
7
7
0.058 0.081
0.088 0.18
0.078 0.27
0.080 0.48
8
0.060 0.17
0.065 0.28
0.062 0.67
0.082 1.10
circle of radius P = l/9. On the other hand, not doiug circle relaxation outside this circle saves considerable work (depending on 9), since this relaxation is expensive. We have then tried a relaxation, composed of the steps: (1) localrelaxation, (2) circle relaxation inside the circle of radius l/9, (3) radialrelaxation (on the full disk). With ?his new relaxation strategy, we improved the efficiency of the soher, producing only small changes in the convergence factors. ‘I’abks 4 and 5 present these new factors and times. In order to solve the differential equation we have used the multigrid cyc!es in the fu&m&igrid mode (see [2], [7J). It is a very efficient way to obtain a solution with the discretization accuracy (rather than using the cycles in an iterative manner). Results concerning the application of our solver in the full-mu&rid mode are presented in the next section.
TABLE5 w
CYCLES (&I)
(IMPItOVED VERSION)*
Convergencefactor [CPU time per cycle (se@]
n-i,
M,=4
5
6
7
7
O.PlO 0.077
0.116 0.23
0.116 0.37
0.116 0*67
8
0.114 0.18
0.113 0.32
0.115 0.92
0.115 1.54
“Without local relaxation.
SAUL0 R. M. BARROS
132
of our program (see Section 4) ~~01, and FX~HPAIL We have Vcydes with one rehxation before and one after The bfGol program was performed in the same ,
systemobtainedthrou&*ediscretizationwi*~~ same system solved by our “geometric” multigrid In the tabk we refer to our program as MGRXAR. Convergence andCPUtimesareshowninTable6. I& Table 7 we present results concerning the x&ion of the Poisson eqqt& with m sohtion g=sin(3x+y). Our program, as well n&l, was~inthefnlE-multigrid~Wepresentalsoresults~ve~ tbe~obwcyclies(~~lacal~n)inaurprogram(~~to~ ~WrThe~~~~isadirect~ffforthePoissanequation and uses polar coo&m&s Since ~~01 is not suitable for the fuknukigrid mode and since its use as an iterative solver would got produce competitive xesd& it is not inin this comparison. The results furnished here the error in the uniform norm, in relation to the ~<#ses:~ecaseM~=~~,inw~~wehave thesamenlmnberOa~~ineachdirection,andthecaseM,=M,+2,in which the mesh s&q in polar coordinates, are approximately the same in with the FI~AK program are the exact discretization errors in polar ccnrdinates (up to roundoff error). Through them we see that
TAB=6 COMPARISON OF CONVERGENCE FACTORS
Convergencefactor
[CPUtimepercycle(sec)] Program (&J&)=(7,4) (7S) (796)(797)(696)(6J) 0.16 0.16 0.1150.17 0.16 AM01 0.36 0.75 0.20 0.37 0.68 0.17 0.16 MGOl 0.40 0.12 MGPOLAR
0.0780.0890.0840.090 0.18 0.27 0.48 0.12 0.24 0.08ii 0.058 0.
in unifolm “In seeon&.
“Error
MGFOLARW
BESULTS FORTHEEQUATION WlTIi -AK
EXACT SOLUTION g =
norm.
lAm 4.24 033 1.25 5.21
0.403-a3 0.101 -03 Oas-03 0.115-03 0287-04
1.26 5.34 0.39 1.60 6.79
OAOS -03 0.101 -03 0.458-03 0.115-03 0.287 -04
1.21 5.22 0:31 1.27 5.60
tim& 0.42 1.46 5.38 -
CPU
MGol
0.5uO-03 0.132-03 0.329-04 -
Error”
sin(3X + $f)
CPUtim& Error” CPUtiId Error” Error” cPutim& O-156 - 02 0.26 0.159 - 02 0.32 0.162-02 0.29
0386 -03 OS63 -04 0.491-0s 0.3.25 -03 (991) 0.313-04
(7,7) (WI (7,5) (680
04s WJ fs,s)
MGFOLAR
COIUP-
TABLE 7
E
SAUlD R M. BARROS
REFERENGES Adams, P. Swantn#rbes and R Sweet, FISEF+E A packageof fortran subprogpms iix the sbhion obsepuabbellipticpartialdifferential equations (Vexsion3, June1919),NationalCenterfbr Atmos~ericResearch Boulder,CO 80307. A. Bran&, Mu&grid Techiques: 1984 Guide with Applications to Fluid Dynamics, GMIMtudie No. 85, GMD, St. k-tin, 1984. JoLhden, M-m fiir die Poisso&lei&mg in Kreisund Ftinggebiet under varwendnng beak Ibordinaten, Diplomrbeit, Inst. f&r Angewandte Maths , univ. Bonn, 1981.
j.
Multi@ fm PO&Wcoordinates
135
4 V. Pereym, Iteaateddefened correctionsfor no&near bamdary value problems, hrarma. M&. l&111-125 ( 5 J.R~pndK.~~o~~~offinitedifferenceandfiniteelement eqgathsbyalgeb&mukgrid(AMG),in MulrigridMethodi~jbrZ~and ZX~Equatfona(D. J.PaddonandH.Ho~Eds.),Inst.of~cs and its AppkWns Gakence Ser.,New Sex.,No. 3, clarendon Press, Oxford, 1985, m 169-al2 6 K. St&n, A multi-@d pqgam to solve AU-c(x,y)U=f(x,y) (on Q), U= g(x, II)(on Ws 0~1NoxhreamgularBounded Domains 62,IMA Report!!E!ma& GMD, Si. Am 198% 7 K. Stiiben and U. ‘I@Mu&rid metbods: Fundamental algorithms, dpaabtempna;lysisd~,in~NotesinMcrdzemotics,vol. 960(W. Haddwsdhand U. Trcaakg, Eds.), Springer,1982, pp. l-176.