C. R. Acad. Sci. Paris, t. 328, Skrie I, p. 73-80, Analyse numtkiquelNumerica/ Analysis
1999
Domain
methods
decomposition
Jacques-Louis
LIONS
a CollBge
de France,
” Institut 4, place
universitaire Jussieu,
(Rqu
et accept6
‘,
3, rue
Olivier d’Ulm,
de France 75231 Paris 1~ 23 novembre
PIRONNEAU 75231
Paris
for CAD
b
redex
et laboratoire d’Analysr cedex 05, France
05,
France numtrique,
Universitk
Pierre-et-Marie-Curie,
1998)
Cette Note ainsi que les deux precedentes sont dedie’es a la memoire de Jean Leray
Abstract.
Constructive solid geometry (CSG) in CAD leads to Domain decompositxons which are based on primitive shapes, as briefly explained in the introduction below. A special role is played by “holes”, which can be viewed in several different ways. We combine this remark with the method of Virtual controls. In a previous note [7] where this method was introduced, the distributions resulting from the decomposition and Green’s formula where thought of as virtual controls - here the virtual controls are introduced a priori with support in the holes, or outside the domain, or in the intersections of the domains. They are then chosen so as to (approximately) satisfy the Boundary conditions, which is possible by virtue of approximate controlability re.sults (which have to be proven!). These ideas are presented here on an example (Section 1) which is certainly not the most general one but which is sufficient to show how everything extends to very many other situations. Algorithms, following the Approximate controllability lemma, are given in Section 2 and a numerical example is presented in Section 4. 0 Academic des Sciences/Elsevier, Paris
M&odes
RkisumC.
Note
prkentie
de d&composition
de domuine
pour
la CA0
En CA0 la description des solides pur Glome’trie constructive conduit a des methodes de decomposition de domaine base’es sur les <(formes primitives M, comme il est explique’ brievement darts I’introduction ci-dessous. Un role special est joue par les R trous B et ils peuvent &tre approches de multiples faEons. Nous combinons cette observation avec la methode des Controles virtuels, introduite dans une note precedente [7], ok la decomposition et la ,formule de Green s’interpretaient comme des contrciles virtuels -ici les contriiles virtuels sont introduits a priori sur des supports dans les trous, ou en dehors. ON dans les intersections des domaines. 11s sont choisis pour que les conditions aux limites soient satisfaites module les erreurs d’approximations, ce qui est possible grace ti UIZ resultat de controlabilite approchee (ce qu’il j&t done demontrer !). Ces id&es sont presentees ici sur un exemple (paragraphe 1) qui n’est certuinement pas lc plus gPnt+al mais qui est .s@sant pour comprendre que tout se g&eralise a beaucoup d’autres situations. Ensuite le lemme de Controlabilitt; upprochee et les algorithmes sont donnks dans le paragraphe 2 ainsi qu’un exemple nume’rique au paragraphe 4. 0 Academic des Sciences/Elsevier, Paris
par
0764~4442/99/03280073
Olivier
PIRONNEAU.
0 Acadkmie
des Sciences/Elsevier,
Paris
73
J.-L. lions,
0. Pironneau
Introduction This study is motivated by two applications: 1. In Image synthesis and Virtual reality (VR) (ser [l] for instance) scenes are entered by Constructive solid geometry, i.e. each object of the scene is described by set operations on primitive shapeslike cubes, cylinders, spheres and cones. For instance a cubic room with a table iinside can be described in VRML (the language of VR) as a cube with four cylinder and a brick inside it, the table’s plateau. These elementary objects are positioned at their proper places but never intersected. Rendering on a computer screen is done by triangulation of the surface of each elementary shapeand the painter’s algorithm with a Z-buffer. Scientific computing (like the computation of the temperature in the room) in such domains requires a translation from VR data structure to CAD data structure, but this operation is difficult. We propose an alternative by domain decomposition. 2. The Chimera method [9] in CFD uses several non-intersecting structured meshesto avoid the problem of mesh generation on difficult geometries. A typical example is a wing and its engine each having its own mesh. Our algorithm proves the method. In both casesthe domain 12is a union and substraction of simple domains. We assumethat each domain can be meshedseparately and that if a domain is subtracted it is possible to surround it by an admissible mesh. We will compute the solution of a PDE on 12by an iterative algorithm involving computation of the PDE on each domain only. Such algorithms already exist: Schwarz [6] if 12is a union of sets; the fictitious domain method [5], [2], [8] if 12 has holes only. But for the general case the problem seems1.0be open. The method is presented in general and it is tested in 2D. There are of course many possibilities and we have been guided by computational efficiency: 1. We assume that we have an efficient general interpolator [3], [4] to compute a function irrespectively of the mesh on which it is originally defined. 2. For a problem discretized on a mesh 5” we avoid computing integrals on curves which are not union of edges of T.
1. Statement
of the problem
Let 12 be a bounded open set of R”, d = 2 or 3, such that 0 = f21 u 122, lli = open sel, 111n
it2 # ti.
Moreover, we assumethat 121has a hole Cl, the position of 0i and C1 being representedon Figure 1.
Figure
I
Remark 1. - The set Cl does not intersect !I*. We can extend what follows to the case where C1 n Q # 8, as we will show elsewhere. 74
Domain
decomposition
methods
for CAD
Remark 2. - All what follows readily extend to the cases where there is also a hole Cz in 02, which does not intersect RI fl 02. Remark 3. - We shall show elsewhere how all what follows extends to the case of a collection of domains R = dII U . . U f1,, with arbitrarily disposed holes. The boundary of 0 consists of I’1 U r2 U X1, with the notations of Figure 1 (rI = that part of i)R1 which is outside flZ and the other way around for r~). Let 4~ *Ij.) = C ,12 (J,)i3. Tc cdn: 7:. i..; which are zero on i?lt). We assume that E L’“(!J),
be defined on the Sobolev space Hh(12) (functions in H’(12)
not necessarikysymmetric
(1)
the n,,‘s are .smoothenough such that the unique continuation theorem holds true.
(2)
(J;j
c
~&mJ
o > 0.
2 (1cc,‘,
(L;j
Let f be given in L’(12) (to fix ideas). We consider the unique .solution U, in HA(G) a(~, ii) = (f, ,C, Vii E H,!,(R),
(3)
We want to give a completely general method of approximation of u based on boundary value problems in domains without holes.
2. Virtual
controls
and approximate
controllability
We are going to use (cfi Figure 2) RI
U
C1 (domain without the hole C,),
‘Dr = neighborhood of ??I such that ??I c V,,
VI n (0, n 0,) = 0.
and we denote by S1 = that part of OR1 which is inside & (and the other way around for S,).
Figure
2
We define VI = {VI ) 711E H1(R1 U C1),vl
= 0
V, = (~2 1 v2 E H1(&,):712 = 0
on
on T’,}, I’,},
75
J.-L. lions,
0. Pironneau
where uij is extended by continuity into ??I in such a way that the ellipticity condition (1) and the uniquenessproperty (2) are valid in Cr. We define nz(ua, &a) in the same way (no extension is needed here since Rz has no hole) for us, 6,~ E V2. In order to avoid any confusion, we define for ‘u11:61 E H1 (Dr ),
The algorithms presented here use only the bilinear forms al ~aIDI, a2 and their adjoints. We now introduce the virtual controls. Let 0 be an open set contained in RI n $12,and let 81 be an open set contained in VI, G1n Cl # (d.
Figure
3
Figure
4
(0 and Gr are shaded on the (magnified) Figures 3 and 4). We set 1s = characteristic function of a set X. We introduce the virtual controls Xl,A2
E L2(0),
x1 +x2 = 0,
pl,ul
E L2(G1),
pl + v1 = 0
on
(4)
G1n pl\cl).
Remark 4. - We could certainly save the index “1” in 1’1YDr , G’r: ,UI! ~1 but we keep it here so as to make the extension to the case where there is also a hole in 622completely straightforward. We still denote by f any extension of f inside the hole Cl, such that f E L2(f2 u cl) (we can extend it by 0 inside Cr) and we decompose ,f in f
=
Sl
+
flDI
+
f2r
fl
E
,.mL2(%
u
Cl),
sm,
E L2(Q):
f2
E
L”(622)
(5)
(where the functions are extended by 0 where they have to be extended). We now define
711 Evi, WI E H'(D,),
u2 E V,
(6)
as the solutions of:
For a given choice of the virtual controls A;, ,~i: v,, the equations (7) (which can be solved in parallel) uniquely define the “state” ~1: lol, ~2.
76
Domain
decomposition
methods
for CAD
Let us admit for a moment the APPROXIMATE CONTROLLABILITY LEMMA. - When the virtual controls &, pi: vi spun L2(0), L2(G1) subject to conditions (4), the functions (~1 + wl)I~c,, ulll~~, UI Is1, u2 Is2 span a dense s&set of L2(dCr) x L2(3Dl) x L2(&) x L2(&).
Then, one can find X1, X2, /rr, ~1 such that: u1 + ‘wl N 0 on w1 N 0
aC, (Y 0 meansof L2(3C1 ) norm arbitrarily small), u120
on dDr,
on
Sr.
1~2~0
on
(8)
S2.
But by virtue of (7) one has 5 = 0 on Sr, % = 0 on aDr, % = 0 on S2 (where & = co-normal derivative attached to A = -c &(Uij&)) so that if we extend ~1 (resp. ‘112,resp. wr) by 0 outside 01 (resp. Ra, resp. Dr) we have “approximately”, outside Si and i3Dr that u1 E H2(R), 1~2E H2(R), w1 E H2(R). According to (7), we have: AvI = jl + X110 + /-~rl~, ,
Awl
=
flz?,
+
hlg,
Au2 =
1
f2
+
X210.
Adding up and using (4), we obtain A(ul + UJ~+ 1~2) = f, ul + IQ + 7/.2 N 0 on i3C1., and ~1 + WI + 2~2= 0 on rl U 172,so that u1 + w1 + ULL:! is an approximation of the solution u. Before showing in Section 3 below how all this can be applied to find an iterative algorithm to compute the virtual controls leading to (8), let us give a hint of the proof of the approximate controllability lemma. By translation, we can assume that fr, flD,f2 are zero. Let (gr, 1~r,II! Z2) be in L2(dC,) x L2(L)Dr) x L2(Sr) x L2(S2) such that (9) (where the measure on the surfaces are not written). We define the functions JYQ!ql~p2 pl E VI, q1 E Hl(Dr), 1)~ E I& as the solutions of 4hth)
G&1&)
=
Ji,,,
.rl&
+
6,
~11%
=
Ji,,,
i/l&
+
i,,
hi1
6
E
&I
v,>
E
4(pshz)
=
L2
hh
V@2
E
Vi,
(1o)
HIP’l):
where a: denotes the adjoint of a,. By taking $1 = 1~1,& = WI, $2 = us in (IO), and using (7) (where the ,f’s are zero), and also using (4), we obtain : Yl = 0
in
61 nG,
Pl = 41
in
G1n
ql=O
(Dl\cl)!
in
p1 =p2
GlnCr: in
0.
(11)
But A*pl = 0 inside Cr and using the uniquenessproperty (we use (2) here), it follows that ~1 = 0 in
Cr, and (analogously) q1 = 0 in
Cr.
(12)
By using the same uniquenesstheorem (in various domains) we obtain that Pl = (I1
in
%\CI,
111=p2
in
El f31Z2.
(13)
But p1 is regular enough such that, by (12), nlac,
= 0
(14)
If we then set 7r = PI in 01, p2 in 02 (they adjust on RI f~ 1~12 by (13)), one has ,4*~ = 0 in 11, r = 0 on rr U r2, r = 0 on dCr , by (14), so that r = 0 in 62. Therefore pr = 0 in Qr u Cr, 77
J.-L. lions,
0. Pironneau
p2 = 0 in f/2 and therefore g1 = 0,/l = 0, /2 = 0. But by (13) q1 = p1 in DD1\C1so that Q = 0 in Vl\Cl and in Cl, hence /I,~ = 0. Remark 5. - Properties analogous to (2), needed here to establish (12), for finit.e element approximations, are known to be an open question in general.
The CSG algorithms
(Constructive
solid geometry
algorithms)
We penalize conditions (8). It leads to the introduction of the cost function
(15)
and we introduce a gradient algorithm for the approximation of
(subject to (4)). We introduce jol: 91,p2 by:
Then, for p > 0 sufficiently small, one defines:
(17)
One computes u;“, ,~o~+l, ZL;~+~(in parallel) by (7) where X1 = X;+l (in parallel) $+‘, qE+l, $+l by (15), and one proceeds.
etc, and then one computes
R.&ark 6. - Everything applies with trivial changesto the caseswhere other Boundary conditions are used on I1 U P2 and on X1. We shall come back to that. Remark 7. - All what has been presented extends to higher order boundary value problems and to systems of equations. Remark 8. - Problems are (of course) much more difficult in the nonlinear cases.
78
Domain
I
,
60+!\
“;*
decomposition
methods
for CAD
/..
t -i Ftgure 5. - The spanner is the union of a circle 92 (upper left) with a rectangle i>t (upper right) minus two wholes Cr. 6’2 which, for computational purpose, are surrounded by triangulated artifical domains (lower left triangulation for the mouth and lower right for the spanner hole). The solution of the PDE after 11 iterations is displaid in the: center: no jump discontinuities are visible. A convergence history of the control cost function .I and of the L’ error with the standard FEM solution is given on the lower left graph. On the lower right graph a one dimensional plot of the solution is shown on a line joining the uppet corner of the spanner mouth with the center of the left hole; iterations I. 2 and the standard FEM solution are shown. the difference with the solution at iteration 11 is not visible.
We cunjecrure that the approximate controllability lemma is still (in general) valid for nonl.inear problems, but this has not been proven. Remurk 9. - For control problems (i.e. problems where there is a “real” control, say ,u) one can still apply the above remarks. We deal then with virtual and real controls. The cost function (15)
79
J.-L. lions,
0. Pironneau
becomes now a function of the virtual controls and of the real ones and one can apply conjugate gradients algorithms. With different techniques (based on Domain decomposition methods) this idea has been used in the Note [7] of the authors. Remark 10. - Of course the conditions (8) could be treated by other methods than penalty, such as Lagrange multipliers, augmented Lagrangians etc...
4. A numerical
example
We wish to compute the temperature U(X) in a spanner. The geometry R is entered by set operations on a rectangle Ui. a cercle &, a hole Cr and the mouth Cz (see Figure 5):
The partial differential equation is: -Au = 0 in 52, z/,Ic~~1Ua~2 = 100, u[~s~-~c,~~c~ =: 0, where Ci is the mouth of the spanner and Cz the hole. In order to account for the holes we construct two arbitrary domains 0s and 624around Ci and C;lz and included respectively in Ri and Ra. So we use the Virtual control algorithm with four subdomains. Discretization is done with the Finite element method of order one on triangles; we have used freefem+ (see [4]) to handle the interpolations on different meshesin an efficient way. Figure 5 shows the 4 subdomains, the solution of the Laplace equation obtained, a convergence history and a one-dimensional plot of X H u(x(X), v(X)) w here X w {z(X), :c/(X)} is the line joining the upper right corner of the spanner mouth with the center of the hole.
Conclusion This example showsthat the method is numerically efficient. However more work needsto be done, one to justify convergence at the discrete level (an inf-sup condition is likely to be necessary for compatibility between meshes)and two to treat extrusions like the mouth of the spannerhere, because in practice extrusion of a set D from R is done by R\(R fl D) while here we require that .D c R.
References [I] Burden G., Coiffe Ph., Virtual Reality Technology, Wiley, New York, 1994. [2] Glowinski R., Pan T., Periaux J., A fictitious domain method for external incompressible flows modelled by Navier-Stokes equations, Univ. Houston internal report, 1992. [3] Hecht F., Lions J.-L., Pironneau O., Algorithmes de points fixes et de controle pour la CA0 par CSG, (to appear). [4] Hecht F., Pironneau O., Handling multiple triangulations in freefem+, INRIA report, (to appear). [S] Kuznetsov Y., Trufanov D., Two-stage fictitious components method for solving the wave Helmoltz equation, Sov. J. Numer. Anal. Math. Mode. III (5) (1988) 371-391. ]6] Lions P.-L., On the Schwarz alternating method I, in: First International Symposium on Domain Decomposition Methods for Partial Differential Equations (Paris, 1987). SIAM, Philadelphia, PA, 1988, pp. 142. 17) Lions J.-L., Pironneau O., Algorithmes paralleles pour la solution de problemes aux limites. C. R. Acad. Sci. Paris 327 Serie I ( 1998) 947-952. ]8] Pironneau 0.. Fictitious Domain Methods versus Boundary Fitted Meshes. Proc Regional Conference in Applied Math and Engineering, E. OAate (Ed.). La Corufia, Spain, 1993. [9] Steger J.L., The Chimera method of flow simulation. Workshop on applied CFD, Univ. of Tennessee SpaceInstitute. August I99 I.
80