Schwarz Methods: To Symmetrize or Not to Symmetrize

0.39 MB·English

by  Michael Holst

SCHWARZ METHODS: TO SYMMETRIZE OR NOT TO SYMMETRIZE
MICHAELHOLSTANDSTEFANVANDEWALLE
ABSTRACT. A preconditioning theory is presented which establishes sufficient condi-
tions for multiplicative and additive Schwarz algorithms to yield self-adjoint positive
definite preconditioners. It allows for the analysis and use of non-variational and non-
convergent linear methods as preconditioners for conjugate gradient methods, and it is
appliedtodomaindecompositionandmultigrid. Itisillustratedwhysymmetrizingmay
0 be a bad idea for linear methods. It is conjectured that enforcing minimal symmetry
1 achievesthebestresultswhencombinedwithconjugategradientacceleration.Also,itis
0
shownthatabsenceofsymmetryinthelinearpreconditionerisadvantageouswhenthe
2
linear method is accelerated by using the Bi-CGstab method. Numerical examples are

n
presentedfortwotestproblemswhichillustratethetheoryandconjectures.
a
J

8

]
A CONTENTS
N 1. Introduction 1
. 2. Krylovaccelerationoflineariterativemethods 2
h
2.1. Backgroundmaterialandnotation 2
t
a 2.2. Linearmethods 3
m
2.3. KrylovaccelerationofSPDlinearmethods 5
[ 2.4. Krylovaccelerationofnonsymmetriclinearmethods 6

3. MultiplicativeSchwarzmethods 6
1
v 3.1. Aproductoperator 7
2 3.2. Multiplicativedomaindecomposition 7
6 3.3. Multiplicativemultigrid 10
3
4. AdditiveSchwarzmethods 13
1
4.1. Asumoperator 13
.
1
4.2. Additivedomaindecomposition 14
0
4.3. Additivemultigrid 16
0
1 5. Tosymmetrizeornottosymmetrize 17
: 6. Numericalresults 17
v
6.1. Example1 17
i
X 6.2. Example2 20
r 7. Concludingremarks 22
a
8. Acknowledgments 23
References 23
1. INTRODUCTION
Domain decomposition (DD) and multigrid (MG) methods have been studied exten-
sively in recent years, both from a theoretical and numerical point of view. DD methods
were first proposed in 1869 by H. A. Schwarz as a theoretical tool in the study of ellip-
tic problems on non-rectangular domains [22]. More recently, DD methods have been
Date:August7,1995.
Keywordsandphrases. multigrid,domaindecomposition,Krylovmethods,Schwarzmethods,conju-
gategradients,Bi-CGstab.
ThisworkwassupportedinpartbytheNSFunderCooperativeAgreementNo.CCR-9120008.
1
2 M.HOLSTANDS.VANDEWALLE
reexamined for use as practical computational tools in the (parallel) solution of general
elliptic equations on complex domains [16]. MG methods were discovered much more
recently [10]. They have been extensively developed both theoretically and practically
sincethelateseventies[6,11],andtheyhaveproventobeextremelyefficientforsolving
very broad classes of partial differential equations. Recent insights in the product nature
ofcertainMGmethodshaveledtoaunifiedtheoryofMGandDDmethods,collectively
referredtoasSchwarzmethods[5,9,27].
Inthispaper,weconsideradditiveandmultiplicativeSchwarzmethodsandtheiraccel-
eration with Krylov methods, for the numerical solution of self-adjoint positive definite
(SPD) operator equations arising from the discretization of elliptic partial differential
equations. The standard theory of conjugate gradient acceleration of linear methods re-
quires that a certain operator associated with the linear method – the preconditioner –
be symmetric and positive definite. Often, however, as in the case of Schwarz-based
preconditioners,thepreconditionerisknownonlyimplicitly,andsymmetryandpositive
definiteness are not easily verified. Here, we try to construct natural sets of sufficient
conditions that are easily verified and do not require the explicit formulation of the pre-
conditioner. Moreprecisely,wederiveconditionsfortheconstituentcomponentsofMG
andDDalgorithms(smoother,subdomainsolver,transferoperators,etc.),thatguarantee
symmetry and positive definiteness of the preconditioning operator which is (explicitly
orimplicitly)definedbytheresultingSchwarzmethod.
We examine the implications of these conditions for various formulations of the stan-
dardDDandMGalgorithms. Thetheorywedevelophelpstoexplaintheoftenobserved
behavior of a poor or even divergent MG or DD method which becomes an excellent
preconditioner when accelerated by a conjugate gradient method. We also investigate
the role of symmetry in linear methods and preconditioners. Both analysis and numer-
ical evidence suggest that linear methods should not be symmetrized when used alone,
and only minimally symmetrized when accelerated by conjugate gradients, in order to
achieve the best possible convergence results. In fact, the best results are often obtained
when a very nonsymmetric linear iteration is used in combination with a nonsymmetric
systemsolversuchasBi-CGstab,eventhoughtheoriginalproblemisSPD.
Theoutlineofthepaperisasfollows. Webeginin§2byreviewingbasiclinearmeth-
ods for SPD linear operator equations, and examine Krylov acceleration strategies. In
§3 and §4, we analyze multiplicative and additive Schwarz preconditioners. We de-
velop a theory that establishes sufficient conditions for the multiplicative and additive
algorithms to yield SPD preconditioners. This theory is used to establish sufficient con-
ditions for multiplicative and additive DD and MG methods, and allows for analysis of
non-variational and even non-convergent linear methods as preconditioners. A simple
lemma, given in §5, illustrates why symmetrizing may be a bad idea for linear methods.
In §6, results of numerical experiments obtained with finite-element-based DD and MG
methodsappliedtosomenon-trivialtestproblemsarereported.
2. KRYLOV ACCELERATION OF LINEAR ITERATIVE METHODS
In this section, we review some background material on self-adjoint linear opera-
tors, linear methods, and conjugate gradient acceleration. More thorough reviews can
befoundin[12,18].
2.1. Background material and notation. Let H be a real finite-dimensional Hilbert
spaceequippedwiththeinner-product(·,·)inducingthenorm(cid:107)·(cid:107) = (·,·)1/2. H canbe
SCHWARZMETHODS:TOSYMMETRIZEORNOTTOSYMMETRIZE 3
thought of as, for example, the Euclidean space Rn, or as an appropriate finite element
space.
The adjoint of a linear operator A ∈ L(H,H) with respect to (·,·) is the unique
operator AT satisfying (Au,v) = (u,ATv), ∀u,v ∈ H. An operator A is called self-
adjoint or symmetric if A = AT; a self-adjoint operator A is called positive definite or
simply positive, if (Au,u) > 0, ∀u ∈ H, u (cid:54)= 0. If A is self-adjoint positive definite
(SPD) with respect to (·,·), then the bilinear form (Au,v) defines another inner-product
onH,whichwedenoteas(·,·) . Itinducesthenorm(cid:107)·(cid:107) = (·,·)1/2.
A A A
The adjoint of an operator M ∈ L(H,H) with respect to (·,·) , the A-adjoint, is the
A
uniqueoperatorM∗satisfying (Mu,v) = (u,M∗v) , ∀u,v ∈ H. Fromthisdefinition
A A
itfollowsthat
M∗ = A−1MTA. (2.1)
M iscalledA-self-adjointifM = M∗,andA-positiveif(Mu,u) > 0, ∀u ∈ H,u (cid:54)= 0.
A
If N ∈ L(H ,H ), then the adjoint of N, denoted as NT ∈ L(H ,H ), is defined as
1 2 2 1
theuniqueoperatorrelatingtheinner-productsinH andH asfollows:
1 2
(Nu,v) = (u,NTv) , ∀u ∈ H , ∀v ∈ H . (2.2)
H2 H1 1 2
Since it is usually clear from the arguments which inner-product is involved, we shall
often drop the subscripts on inner-products (and norms) throughout the paper, except
whennecessarytoavoidconfusion.
We denote the spectrum of an operator M as σ(M). The spectral theory for self-
adjointlinearoperatorsstatesthattheeigenvaluesoftheself-adjointoperatorM arereal
andlieintheclosedinterval[λ (M),λ (M)]definedbytheRayleighquotients:
min max
(Mu,u) (Mu,u)
λ (M) = min , λ (M) = max . (2.3)
min max
u(cid:54)=0 (u,u) u(cid:54)=0 (u,u)
Similarly, if an operator M is A-self-adjoint, then its eigenvalues are real and lie in the
interval defined by the Rayleigh quotients generated by the A-inner-product. A well-
known property is that if M is self-adjoint, then the spectral radius of M, denoted as
ρ(M), satisfies ρ(M) = (cid:107)M(cid:107). This property can also be shown to hold in the A-norm
forA-self-adjointoperators(or,moregenerally,forA-normaloperators[1]).
Lemma2.1. IfAisSPDandM isA-self-adjoint,then ρ(M) = (cid:107)M(cid:107) .
A
2.2. Linear methods. Given the equation Au = f, where A ∈ L(H,H) is SPD, con-
sider the preconditioned equation BAu = Bf, with B ∈ L(H,H). The operator B,
the preconditioner, is usually chosen so that a Krylov or Richardson method applied to
the preconditioned system has some desired convergence properties. A simple linear
iterativemethodemployingtheoperatorB takestheform
un+1 = un −BAun +Bf = (I −BA)un +Bf, (2.4)
where the convergence behavior of (2.4) is determined by the properties of the so-called
errorpropagationoperator,
E = I −BA. (2.5)
The spectral radius of the error propagator E is called the convergence factor for the
linearmethod,whereasthenormisreferredtoasthecontractionnumber. Werecalltwo
well-knownlemmas;seeforexample[17]or[20].
Lemma 2.2. For arbitrary f and u0, the condition ρ(E) < 1 is necessary and sufficient
forconvergenceofthelinearmethod(2.4).
4 M.HOLSTANDS.VANDEWALLE
Lemma 2.3. The condition (cid:107)E(cid:107) < 1, or the condition (cid:107)E(cid:107) < 1, is sufficient for
A
convergenceofthelinearmethod(2.4).
We now state a series of simple lemmas that we shall use repeatedly in the following
sections. Theirshortproofsareaddedforthereader’sconvenience.
Lemma2.4. IfAisSPD,thenBAisA-self-adjointifandonlyifB isself-adjoint.
Proof. Note that: (ABAu,v) = (BAu,Av) = (Au,BTAv). The lemma follows since
BA = BTAifandonlyifB = BT. (cid:3)
Lemma2.5. IfAisSPD,thenE isA-self-adjointifandonlyifB isself-adjoint.
Proof. Notethat: (AEu,v)=(Au,v)−(ABAu,v)=(Au,v)−(Au,(BA)∗v)=(Au,(I−
(BA)∗)v). Therefore, E∗ = E if and only if BA = (BA)∗. By Lemma 2.4, this holds if
andonlyifB isself-adjoint. (cid:3)
Lemma2.6. IfAandB areSPD,thenBAisA-SPD.
Proof. ByLemma2.4,BAisA-self-adjoint. Also,wehave(ABAu,u) = (BAu,Au) =
(B1/2Au,B1/2Au) > 0, ∀u (cid:54)= 0.Hence,BAisA-positive,andtheresultfollows. (cid:3)
Lemma2.7. IfAisSPDandB isself-adjoint,then(cid:107)E(cid:107) = ρ(E).
A
Proof. ByLemma2.5,E isA-self-adjoint. ByLemma2.1theresultfollows. (cid:3)
Lemma2.8. IfE∗ istheA-adjointofE,then(cid:107)E(cid:107)2 = (cid:107)EE∗(cid:107) .
A A
Proof. TheprooffollowsthatofafamiliarresultfortheEuclidean2-norm[12]. (cid:3)
Lemma2.9. IfAandB areSPD,andE isA-non-negative,then(cid:107)E(cid:107) < 1.
A
Proof. ByLemma2.5,E isA-self-adjoint. AsE isA-non-negative,itholds(Eu,u) ≥
A
0, or (BAu,u) ≤ (u,u) . By Lemma 2.6, BA is A-SPD, and we have that 0 <
A A
(BAu,u) ≤ (u,u) , ∀u (cid:54)= 0,which,by(2.3),impliesthat0 < λ ≤ 1, ∀λ ∈ σ(BA).
A A i i
Thus,ρ(E) = 1−min λ < 1.Finally,byLemma2.7,wehave(cid:107)E(cid:107) = ρ(E). (cid:3)
i i A
Wewillalsohaveuseforthefollowingtwosimplelemmas.
Lemma2.10. IfAisSPDandB isself-adjoint,andE issuchthat:
−C (u,u) ≤ (Eu,u) ≤ C (u,u) , ∀u ∈ H,
1 A A 2 A
forC ≥ 0andC ≥ 0,thenρ(E) = (cid:107)E(cid:107) ≤ max{C ,C }.
1 2 A 1 2
Proof. By Lemma 2.5, E is A-self-adjoint, and by (2.3) λ (E) and λ (E) are
min max
boundedby−C andC ,respectively. TheresultthenfollowsbyLemma2.7. (cid:3)
1 2
Lemma2.11. IfAandB areSPD,thenLemma2.10holdsforsomeC < 1.
2
Proof. ByLemma2.6,BAisA-SPD,whichimpliesthattheeigenvaluesofBAarereal
and positive. Hence, we must have that λ (E) = 1 − λ (BA) < 1, ∀i. Since C in
i i 2
Lemma2.10boundsthelargestpositiveeigenvalueofE,wehavethatC < 1. (cid:3)
2
SCHWARZMETHODS:TOSYMMETRIZEORNOTTOSYMMETRIZE 5
2.3. KrylovaccelerationofSPDlinearmethods. Theconjugategradientmethodwas
developed by Hestenes and Stiefel [13] as a method for solving linear systems Au = f
with SPD operators A. In order to improve convergence, it is common to precondition
thelinearsystembyanSPDpreconditioningoperatorB ≈ A−1,inwhichcasethegener-
alizedorpreconditionedconjugategradientmethodresults([8]). Ourgoalinthissection
is to briefly review some relationships between the contraction number of a basic linear
preconditionerandthatoftheresultingpreconditionedconjugategradientalgorithm.
Westartwiththewell-knownconjugategradientcontractionbound([12]):
(cid:32) (cid:33)i+1
2
(cid:107)ei+1(cid:107) ≤ 2 1− (cid:107)e0(cid:107) = 2δi+1 (cid:107)e0(cid:107) . (2.6)
A (cid:112) A cg A
1+ κ (BA)
A
The ratio of extreme eigenvalues of BA appearing in the derivation of the bound gives
rise to the generalized condition number κ (BA) appearing above. This ratio is often
A
mistakenly called the (spectral) condition number κ(BA); in fact, since BA is not self-
adjoint, this ratio is not in general equal to the usual condition number (this point is
discussed in great detail in [1]). However, the ratio does yield a condition number in the
A-norm. ThefollowinglemmaisaspecialcaseofCorollary4.2in[1].
Lemma2.12. IfAandB areSPD,then
λ (BA)
κ (BA) = (cid:107)BA(cid:107) (cid:107)(BA)−1(cid:107) = max . (2.7)
A A A
λ (BA)
min
Remark 2.13. Often a linear method requires a parameter α in order to be convergent,
leading to an error propagator of the form E = I − αBA. Equation (2.7) shows that
the A-condition number does not depend on the particular choice of α. Hence, one can
usetheconjugategradientmethodasanacceleratorforthemethodwithoutaparameter,
avoidingthepossiblycostlyestimationofagoodα.
The following result gives a bound on the condition number of the operator BA in
terms of the extreme eigenvalues of the error propagator E = I −BA; such bounds are
oftenusedintheanalysisoflinearpreconditioners(cf. Proposition5.1in[26]). Wegive
ashortproofofthisresultforcompleteness.
Lemma2.14. IfAandB areSPD,andE issuchthat:
−C (u,u) ≤ (Eu,u) ≤ C (u,u) , ∀u ∈ H, (2.8)
1 A A 2 A
forC ≥ 0andC ≥ 0,thentheabovemustholdwithC < 1,anditfollowsthat:
1 2 2
1+C
1
κ (BA) ≤ .
A
1−C
2
Proof. First, since A and B are SPD, by Lemma 2.11 we have that C < 1. Since
2
(Eu,u) = (u,u) −(BAu,u) ,itisclearthat
A A A
(1−C )(u,u) ≤ (BAu,u) ≤ (1+C )(u,u) , ∀u ∈ H.
2 A A 1 A
ByLemma2.6,BAisA-SPD.Itseigenvaluesarerealandpositive,andlieintheinterval
definedbytheRayleighquotientsgeneratedbytheA-inner-product. Hence,thatinterval
isgivenby[(1−C ),(1+C )],andbyLemma2.12theresultfollows. (cid:3)
2 1
Remark 2.15. Even if a linear method is not convergent, it may still be a good precon-
ditioner. If it is the case that C << 1, and if C > 1 does not become too large, then
2 1
κ (BA) will be small and the conjugate gradient method will converge rapidly, even
A
thoughthelinearmethoddiverges.
6 M.HOLSTANDS.VANDEWALLE
IfonlyaboundonthenormoftheerrorpropagatorE = I−BAisavailable,thenthe
following result can be used to bound the condition number of BA. This result is used
forexamplein[27].
Corollary2.16. IfAandB areSPD,and(cid:107)I −BA(cid:107) ≤ δ < 1,then
A
1+δ
κ (BA) ≤ . (2.9)
A
1−δ
Proof. ThisfollowsimmediatelyfromLemma2.14withδ = max{C ,C }. (cid:3)
1 2
The next result connects the contraction number of the preconditioner to the contrac-
tionnumberofthepreconditionedconjugategradientmethod. Itshowsthattheconjugate
gradientmethodalwaysacceleratesalinearmethod(iftheconditionsofthelemmahold).
Lemma2.17. IfAandB areSPD,and(cid:107)I −BA(cid:107) ≤ δ < 1,thenδ < δ.
A cg
Proof. Anabbreviatedproofappearsin[27],amoredetailedproofin[14]. (cid:3)
2.4. Krylov acceleration of nonsymmetric linear methods. The convergence theory
of the conjugate gradient iteration requires that the preconditioned operator BA be A-
self-adjoint (see [2] for more general conditions), which from Lemma 2.4 requires that
B be self-adjoint. If a Schwarz method is employed which produces a nonsymmetric
operator B, then although A is SPD, the theory of the previous section does not ap-
ply,andanonsymmetricsolversuchasconjugategradientsonthenormalequations[2],
GMRES [21], CGS [23], or Bi-CGstab [25] must be used for the now non-A-SPD pre-
conditionedsystem,BAu = Bf.
The conjugate gradient method for SPD problems has several nice properties (good
convergence rate, efficient three-term recursion, and minimization of the A-norm of the
error at each step), some of which must be given up in order to generalize the method to
nonsymmetric problems. For example, while GMRES attempts to maintain a minimiza-
tion property and a good convergence rate, the three-term recursion must be sacrificed.
Conjugate gradients on the normal equations maintains a minimization property as well
as the efficient three-term recursion, but sacrifices convergence speed (the effective con-
ditionnumberisthesquareoftheoriginalsystem). MethodssuchasCGSandBi-CGstab
sacrifice the minimization property, but maintain good convergence speed and the effi-
cientthree-termrecursion. Forthesereasons,methodssuchasCGSandBi-CGstabhave
becomethemethodsofchoiceinmanyapplicationsthatgiverisetononsymmetricprob-
lems. Bi-CGstab has been shown to be more attractive than CGS in many situations due
to the more regular convergence behavior [25]. In addition, Bi-CGstab does not require
the application of the adjoint of the preconditioning operator, which can be difficult to
implementinthecaseofsomeSchwarzmethods.
In§6,weshallusethepreconditionedBi-CGstabalgorithmtoacceleratenonsymmet-
ric Schwarz methods. In a sequence of numerical experiments, we shall compare the
effectivenessofthisapproachwithunacceleratedsymmetricandnonsymmetricSchwarz
methods,andwithsymmetricSchwarzmethodsacceleratedwithconjugategradients.
3. MULTIPLICATIVE SCHWARZ METHODS
We develop a preconditioning theory of product algorithms which establishes suffi-
cient conditions for producing SPD preconditioners. This theory is used to establish
sufficientSPDconditionsformultiplicativeDDandMGmethods.
SCHWARZMETHODS:TOSYMMETRIZEORNOTTOSYMMETRIZE 7
3.1. Aproductoperator. Consideraproductoperatoroftheform:
¯
E = I −BA = (I −B A)(I −B A)(I −B A), (3.1)
1 0 1
¯
where B ,B and B are linear operators on H, and where A is, as before, an SPD
1 0 1
¯
operator on H. We are interested in conditions for B ,B and B , which guarantee that
1 0 1
the implicitly defined operator B is self-adjoint and positive definite and, hence, can be
acceleratedbyusingtheconjugategradientmethod.
Lemma 3.1. Sufficient conditions for symmetry and positivity of operator B, implicitly
definedby(3.1),are:
(1) B¯ = BT ;
1 1
(2) B = BT ;
0 0
(3) (cid:107)I −B A(cid:107) < 1;
1 A
(4) B non-negativeonH .
0
Proof. By Lemma 2.5, in order to prove symmetry of B, it is sufficient to prove that E
isA-self-adjoint. Byusing(2.1),weget
E∗ = A−1ETA
= A−1(I −ABT)(I −ABT)(I −AB¯T)A
1 0 1
= (I −BTA)(I −BTA)(I −B¯TA)
1 0 1
¯
= (I −B A)(I −B A)(I −B A) = E,
1 0 1
whichfollowsfromconditions1and2.
Next, we prove that (Bu,u) > 0, ∀u ∈ H, u (cid:54)= 0. Since A is non-singular, this is
equivalenttoprovingthat(BAu,Au) > 0. Usingcondition1,wehavethat
(BAu,Au) = ((I −E)u,Au)
= (u,Au)−((I −BTA)(I −B A)(I −B A)u,Au)
1 0 1
= (u,Au)−((I −B A)(I −B A)u,A(I −B A)u)
0 1 1
= (u,Au)−((I −B A)u,A(I −B A)u)+(B w,w),
1 1 0
where w = A(I − B A)u. By condition 4, we have that (B w,w) ≥ 0. Condition 3
1 0
impliesthat((I −B A)u,A(I −B A)u) < (u,Au)foru (cid:54)= 0. Thus,thefirsttwoterms
1 1
in the sum above are together positive, while the third one is non-negative, so that B is
positive. (cid:3)
Corollary3.2. IfB = BT,thencondition3inLemma3.1isequivalenttoρ(I−B A) <
1 1 1
1.
Proof. ThisfollowsdirectlyfromLemma2.1andLemma2.5. (cid:3)
3.2. Multiplicativedomaindecomposition. Giventhefinite-dimensionalHilbertspace
H, consider J spaces H , k = 1,...,J, together with linear operators I ∈ L(H ,H),
k k k
null(I ) = {0}, such that I H ⊆ H = (cid:80)J I H . We also assume the existence
k k k k=1 k k
of another space H , an associated operator I such that I H ⊆ H, and some linear
0 0 0 0
operatorsIk ∈ L(H,H ),k = 0,...,J. Fornotationalconvenience,weshalldenotethe
k
inner-products on H by (·,·) (without explicit reference to the particular space). Note
k
thattheinner-productsondifferentspacesneednotberelated.
In a domain decomposition context, the spaces H , k = 1,...,J, are typically as-
k
sociated with local subdomains of the original domain on which the partial differential
equation is defined. The space H is then a space associated with some global coarse
0
mesh. The operators I ,k = 1,...,J, are usually inclusion operators, while I is an
k 0
8 M.HOLSTANDS.VANDEWALLE
interpolation or prolongation operator (as in a two-level MG method). The operators
Ik,k = 1,...,J, are usually orthogonal projection operators, while I0 is a restriction
operator(again,asinatwo-levelMGmethod).
The error propagator of a multiplicative DD method on the space H employing the
subspacesI H hasthegeneralform[9]:
k k
E = I −BA = (I −I R¯ IJA)···(I −I R I0A)···(I −I R IJA), (3.2)
J J 0 0 J J
¯
where R and R , k = 1,...,J, are linear operators on H , and R is a linear operator
k k k 0
on H . Usually the operators R¯ and R are constructed so that R¯ ≈ A−1 and R ≈
0 k k k k k
A−1, where A is the operator defining the subdomain problem in H . Similarly, R is
k k k 0
constructed so that R ≈ A−1. Actually, quite often R is a “direct solve”, i.e., R =
0 0 0 0
A−1. The subdomain problem operator A is related to the restriction of A to H . We
0 k k
say that A satisfies the Galerkin conditions or, in a finite element setting, that it is
k
variationallydefinedwhen
A = IkAI , Ik = IT. (3.3)
k k k
Recallthatthesuperscript“T”istobeinterpretedastheadjointinthesenseof(2.2),i.e.,
withrespecttotheinner-productsinH andH .
k
In the case of finite element, finite volume, or finite difference discretization of an
elliptic problem, conditions (3.3) can be shown to hold naturally for both the matrices
and the abstract weak form operators for all subdomains k = 1,...,J. For the coarse
spaceH ,often(3.3)mustbeimposedalgebraically.
0
Propagator(3.2)canbethoughtofastheproductoperator(3.1),bychoosing
1 J
(cid:89) (cid:89)
I −B¯ A = (I −I R¯ IkA), B = I R I0 , I −B A = (I −I R IkA),
1 k k 0 0 0 1 k k
k=J k=1
¯
where B and B are known only implicitly. (Note that we take the convention that the
1 1
first term in the product appears on the left.) This identification allows for the use of
¯
Lemma 3.1 to establish sufficient conditions on the subdomain operators R , R and R
k k 0
toguaranteethatmultiplicativedomaindecompositionyieldsanSPDoperatorB.
Theorem 3.3. Sufficient conditions for symmetry and positivity of the multiplicative do-
maindecompositionoperatorB,implicitlydefinedby(3.2),are:
(1) Ik = c IT , c > 0, k = 0,··· ,J ;
k k k
(2) R¯ = RT , k = 1,··· ,J ;
k k
(3) R = RT ;
(cid:13)0 0 (cid:13)
(4) (cid:13)(cid:81)J (I −I R IkA)(cid:13) < 1;
(cid:13) k=1 k k (cid:13)
A
(5) R non-negativeonH .
0 0
Proof. WeshowthatthesufficientconditionsofLemma3.1aresatisfied. First,weprove
that B¯ = BT, which, by Lemma 2.5, is equivalent to proving that (I − B A)∗ =
1 1 1
¯
(I −B A). Byusing(2.1),wehave
1
(cid:32) (cid:33)∗ (cid:32) (cid:33)T
J J 1
(cid:89) (cid:89) (cid:89)
(I −I R IkA) = A−1 (I −I R IkA) A = (I−(Ik)TRT(I )TA),
k k k k k k
k=1 k=1 k=J
¯
which equals (I −B A) under conditions 1 and 2 of the theorem. The symmetry of B
1 0
followsimmediatelyfromconditions1and3;indeed,
BT = (I R I0)T = (I0)TRT(I )T = (c I )R (c−1I0) = I R I0 = B .
0 0 0 0 0 0 0 0 0 0 0 0
SCHWARZMETHODS:TOSYMMETRIZEORNOTTOSYMMETRIZE 9
Bycondition4ofthetheorem,condition3ofLemma3.1holdstrivially. Thetheorem
followsbyrealizingthatcondition4ofLemma3.1isalsosatisfied,since,
(B u,u) = (I R I0u,u) = (R I0u,ITu) = c−1(R I0u,I0u) ≥ 0, ∀u ∈ H .
0 0 0 0 0 0 0
(cid:3)
Remark 3.4. Note that one sweep through the subdomains, followed by a coarse prob-
lem solve, followed by another sweep through the subdomains in reversed order, gives
rise an error propagator of the form (3.2). Also, note that no conditions are imposed on
the nature of the operators A associated with each subdomain. In particular, the theo-
k
rem does not require that the variational conditions are satisfied. While it is natural for
condition(3.3)toholdbetweenthefinespaceandthespacesassociatedwitheachsubdo-
main, these conditions are often difficult to enforce for the coarse problem. Violation of
variational conditions can occur, for example, when complex coefficient discontinuities
donotliealongelementboundariesonthecoarsemesh(wepresentnumericalresultsfor
such a problem in §6). The theorem also does not require that the overall multiplicative
DDmethodbeconvergent.
Remark 3.5. The results of the theorem apply for abstract operators on general finite-
dimensional Hilbert spaces with arbitrary inner-products. They hold in particular for
matrix operators on Rn, equipped with the Euclidean inner-product, or the discrete L2
inner-product. In the former case, the superscript “T” corresponds to the standard ma-
trix transpose. In the latter case, the matrix representation of the adjoint is a scalar
multipleofthematrixtranspose;thescalarmaybedifferentfromunitywhentheadjoint
involves two different spaces, and in the case of prolongation and restriction. This pos-
sible constant in the case of the discrete L2 inner-product is absorbed in the factor c in
k
condition 1. This allows for an easy verification of the conditions of the theorem in an
actual implementation, where the operators are represented as matrices, and where the
inner-productsdonotexplicitlyappearinthealgorithm.
Remark 3.6. Condition 1 of the theorem (with c = 1) for k = 1,...,J is usually
k
satisfied trivially for domain decomposition methods. For k = 0, it may have to be
imposedexplicitly. Condition2ofthetheoremallowsforseveralalternativeswhichgive
rise to an SPD preconditioner, namely: (1) use of exact subdomain solvers (if A is a
k
symmetric operator); (2) use of identical symmetric subdomain solvers in the forward
and backward sweeps; (3) use of the adjoint of the subdomain solver on the second
sweep. Condition 3 is satisfied when the coarse problem is symmetric and the solve
is an exact one, which is usually the case. If not, the coarse problem solve has to be
symmetric. Condition 4 in Theorem 3.3 is clearly a non-trivial one; it is essentially
the assumption that the multiplicative DD method without a coarse space is convergent.
ConvergencetheoriesforDDmethodscanbequitetechnicalanddependonsuchthings
as the discretization, the subdomain number, shape, and size, and the regularity of the
solution [5, 9, 27]. However, since variational conditions hold naturally between the
fine space and each subdomain space for nearly any formulation of a DD method, very
general convergence theorems can be derived, if one is not concerned about the actual
rate of convergence. Using the Schwarz theory framework in any of [5, 9, 27], it can
be shown that Condition 4 in Theorem 3.3 (convergence of multiplicative DD without a
coarsespace)holdsifthevariationalconditions(3.3)hold,andifthesubdomainsolvers
R are SPD. A proof of this result may be found for example in [14]. Condition 5 is
k
satisfiedforexamplewhenthecoarseproblemisSPDandthesolveisexact.
10 M.HOLSTANDS.VANDEWALLE
Consider now the case when the subspaces together do not span the entire space, ex-
cept when the coarse space is included. The above theorem can be applied with R = 0,
0
and by viewing the coarse space as simply one of the spaces H , k (cid:54)= 0. In this case, the
k
errorpropagationoperatorE takestheform:
I−BA = (I−I R¯ IJA)···(I−I R¯ I1A)(I−I R I1A)···(I−I R IJA). (3.4)
J J 1 1 1 1 J J
Thisleadstothefollowingcorollary.
Corollary 3.7. Sufficient conditions for symmetry and positivity of the multiplicative
domaindecompositionoperatorB,implicitlydefinedby(3.4),are:
(1) Ik = c IT , c > 0, k = 1,··· ,J ;
k k k
(2) R¯ = RT , k = 1,··· ,J ;
(cid:13)k k (cid:13)
(3) (cid:13)(cid:81)J (I −I R IkA)(cid:13) < 1.
(cid:13) k=1 k k (cid:13)
A
Remark 3.8. Condition 3 is equivalent to requiring convergence of the overall multi-
plicativeSchwarzmethod. Thisfollowsfromtherelationship
(cid:107)E(cid:107) = (cid:107)E¯∗E¯(cid:107) = (cid:107)E¯(cid:107)2 < 1,
A A A
whereE¯ = (cid:81)J (I −I R IkA).
k=1 k k
Remark3.9. If,inadditiontoconditionsofthecorollary,itholdsthatR = (I1AI )−1,
1 1
i.e., it corresponds to an exact solve with a variationally defined subspace problem op-
eratorinthesenseof(3.3),then
(I −I R¯ I1A)(I −I R I1A) = I −I R I1A,
1 1 1 1 1 1
since I −I (I1AI )−1I A is a projector. Therefore, space H (for example, the coarse
1 1 1 1
space)needstobevisitedonlyonceintheapplicationof(3.4).
3.3. Multiplicative multigrid. Consider the Hilbert space H, J spaces H together
k
with linear operators I ∈ L(H ,H), null(I ) = 0, such that the spaces I H are nested
k k k k k
and satisfy I H ⊆ I H ⊆ ··· ⊆ I H ⊆ H ≡ H. As before we denote the H -
1 1 2 2 J−1 J−1 J k
inner-products by (·,·), since it will be clear from the arguments which inner-product is
intended. Again, the inner-products are not necessarily related in any way. We assume
alsotheexistenceofoperatorsIk ∈ L(H,H ).
k
In a multigrid context, the spaces H are typically associated with a nested hierarchy
k
of successively refined meshes, with H being the coarsest mesh, and H being the fine
1 J
mesh on which the PDE solution is desired. The linear operators I are prolongation
k
operators, constructed from given interpolation or prolongation operators that operate
between subspaces, i.e., Ik ∈ L(H ,H ). The operator I is then constructed (only
k−1 k−1 k k
asatheoreticaltool)asacompositeoperator
I = IJ IJ−1···Ik+2Ik+1, k = 1,...,J −1. (3.5)
k J−1 J−2 k+1 k
ThecompositerestrictionoperatorsIk,k = 1,...,J −1,areconstructedsimilarlyfrom
somegivenrestrictionoperatorsIk−1 ∈ L(H ,H ).
k k k−1
The coarse problem operators A are related to the restriction of A to H . As in the
k k
case of DD methods, we say that A is variationally defined or satisfies the Galerkin
k
conditions when conditions (3.3) hold. It is not difficult to see that conditions (3.3) are
equivalenttothefollowingrecursivelydefinedvariationalconditions:
A = Ik A Ik+1, Ik = (Ik+1)T, (3.6)
k k+1 k+1 k k+1 k
whenthecompositeoperatorsI appearingin(3.3)aredefinedasin(3.5).
k