C R Acad Sci Paris Ser I Analyse numérique Numerical Analysis

Yvon Maday

Description

Niveau: Supérieur
C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938 Analyse numérique/Numerical Analysis Coupling between scalar and vector potentials by the mortar element method Yvon Maday a, Francesca Rapetti b, Barbara I. Wohlmuth c a Laboratoire J.-L. Lions, CNRS & Université Pierre-et-Marie-Curie, boîte 187, 75252 Paris cedex 05, France b Laboratoire des mathématiques J.-A. Dieudonné, CNRS & Université de Nice et Sophia-Antipolis, Parc Valrose, 06108 Nice cedex 02, France c Universität Stuttgart, Mathematisches Institut A, Lehrstuhl 7, Pfaffenwaldring 57, 70569 Stuttgart, Germany Received 20 February 2002; accepted 4 March 2002 Note presented by Olivier Pironneau. Abstract The T – formulation of the magnetic field is widely used in magnetodynamics. It allows the use of a scalar function in the computational domain and a vector quantity only in the conducting parts. Here we propose to approximate these two quantities on different meshes and to couple them by means of the mortar element method. To cite this article: Y. Maday et al., C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938. ? 2002 Académie des sciences/Éditions scientifiques et médicales Elsevier SAS Couplage entre potentiels scalaire et vecteur par la méthode des éléments joints Résumé La formulation T? pour le champ magnétique est largement utilisée en magnétodyna- mique.

• ∂vc


• vc ?


• opérateur de passage d'informationh du maillage de ∂vc


• maillage de vc


• opérateur de relèvément harmonique


• condition ? ·

Sujets

• Elsevier
• Université de Stuttgart
• Universidad de Niza Sophia Antipolis

Informations

• Publié par: profil-nechor-2012
• Nombre de lectures: 58
• Langue: Français

Extrait

C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938

Analyse numérique/Numerical Analysis

Coupling between scalar and vector potentials by the mortar element method a bc Yvon Maday, Francesca Rapetti, Barbara I. Wohlmuth a Laboratoire J.L. Lions, CNRS & Université PierreetMarieCurie, boîte 187, 75252 Paris cedex 05, France b Laboratoire des mathématiques J.A. Dieudonné, CNRS & Université de Nice et SophiaAntipolis, Parc Valrose, 06108 Nice cedex 02, France c Universität Stuttgart, Mathematisches Institut A, Lehrstuhl 7, Pfaffenwaldring 57, 70569 Stuttgart, Germany Received 20 February 2002; accepted 4 March 2002 Note presented by Olivier Pironneau.

Abstract TheT–formulation of the magnetic field is widely used in magnetodynamics. It allows the use of a scalar function in the computational domain and a vector quantity only in the conducting parts. Here we propose to approximate these two quantities on different meshes and to couple them by means of the mortar element method.To cite this article: Y. Maday et al., C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938.
2002 Académie des sciences/Éditions scientifiques et médicales Elsevier SAS

Couplage entre potentiels scalaire et vecteur par la méthode des éléments joints Résumé La formulationT−pour le champ magnétique est largement utilisée en magnétodyna mique. Elle permet d’utiliser une fonction scalaire dans tout le domaine de calcul et une quantité vectorielle seulement dans les conducteurs. On propose ici d’approcher ces deux quantités sur des maillages différents et de les coupler par la méthode des éléments avec joints.: Y.Pour citer cet articleMaday et al., C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938.
2002 Académie des sciences/Éditions scientifiques et médicales Elsevier SAS

Version française abrégée La simulation numérique de phénomènes electromagnétiques basse fréquence fait souvent appel au modèle des courants de Foucault. On s’interesse ici à la formulation en champs magnétique dans un domaine borné deRnotéVcomposé de deux parties : un conducteurV c et un isolantV\V c . Les équations sont alors (1) dans le conducteur et (2) dans l’isolant. Afin de ramener le calcul à des quantités scalaires,  quand c’est possible, on a recours à l’introduction d’un potentiel scalairedéfini dansVet un potentiel   vecteurTlimité àV c tel queH=T− ∇. Partant des équations variationelles définissantH(3) et (4), 1   on obtient la formulation (7) pour déterminerTet. L’espace H 0 (curl;V c )est défini en (5) et H(V )est 0 l’espace classique de Sobolev. Ce problème ne possèdant pas de solution unique, on choisit une condition 1  (de type jauge) qui impose à=+φ,φ∈H(V c ), d’être harmonique dansV c . On aboutit alors à la 0  1 formulation variationelle (10) oùb c (W, v):= −W∇Hv, pourv∈H(V )etW∈H 0 (curl;V c )etHest 0 V c l’opérateur de relèvément harmonique défini en (9). On montre d’abord que (11) a lieu pour une constante 0< γ 1 <1, puis le lemme suivant sur la formea g ((W, w), (V , v)):=a c (W, V )+b c (W, v)+b c (V , w)+ 1 a(w, v),pourV , W∈H 0 (curl;V c )etv, w∈H(V ): 0

Email addresses:[e-mail] (Y. Maday); [e-mail] (F. Rapetti); [e-mail]stuttgart.de (B.I. Wohlmuth).
Tous droits réservés2002 Académie des sciences/Éditions scientifiques et médicales Elsevier SAS. S 1 6 3 1  0 7 3 X ( 0 2 ) 0 2 3 5 3  1 /FLA

933

Y. Maday et al. / C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938

1 EMME L 0.1.–La forme bilinéairea g (∙,∙)est continue et elliptique surH 0 (curl;V c )×H(V ). 0 On en déduit que le problème (10) possède une solution et une seule. Ce problème faisant intervenir une quantité scalaire surVet une quantité vectorielle surV c , on peut penser à utiliser, dans le cadre d’une discrétisation en éléments finis, deux maillages différents : l’un surVpour calculeret indexé parHet l’autre surV c pour calculerTet indexé parh. On peut même envisager que les deux maillages soient contruit indépendamment (ce sera d’ailleurs naturellement le cas si le conducteur bouge dansVau cours du temps). La discrétisation de ce problème peut alors se faire en utilisant des éléments finis de type Nédélec pourT(espace notéX 0;h (V c )) et des éléments finis standards pour(espace notéS 0;H (V )). La formeb c faisant intervenir un relèvement harmonique, doit être discrétisée. On doit tout d’abord définir un tel relèvement discretH h défini en (12) mais aussi un opérateur de passage d’information h du maillage de∂ V c hérité du maillage deVà celui de∂ V c hérité du maillage deV c . Nous proposons de faire cette opération de manière optimale comme indiqué par (13), inspiré par la méthode des éléments avec joint, soit directement en suivant [2], soit à partir de [4] qui conduit à une définition de h par matrice diagonale. On obtient, dès lors que le rapporth/Hest assez petit, le fait que 1 3 le problème (14) possède une solution unique et l’estimation d’erreur (15), valable dès queT∈(H(V c ))  1 3β avec∇ ×T∈(H(V c ))et∈H(V ), pour un réelβ2, oùb h (W, v):= −W∇H h  h v, pour tout V c v∈S 0;H (V )etW∈X 0;h (V c ). La forme matricielle du problème (14) est donnée en (16) oùQest la matrice associée à h ,Sest associée au relèvement harmonique etBcelle associée à la formeb h . La méthode itérative de Gauss–Seidel (17) converge sans paramètre de relaxation.

1. Introduction Low frequency electromagnetic devices are often modeled numerically on the basis of the eddy current formulation [1]. Two main families of formulations are widely used, the one based on magnetic and the one 3 based on electric fields. Here, we restrict ourselves to the magnetic field approach. The entire spaceR, 3 is decomposed in the conducting regionV c and the external regionR\V c . Denoting byH,B,JandE, respectively, the magnetic field, the magnetic flux density, the current density and the electric field, the quasistationary Maxwell equations restricted to the conducting regionV c read as follows: ∇ ×H=J,∇ ×E= −∂ t B,∇ ∙B=0.(1) The densities and the fields are linked by the constitutive properties, i.e.,J=σ E,B=µH, whereµis the magnetic permeability andσσ 0 >0 stands for the electric conductivity. Moreover, we assume that the material parameters are time independent and associated with a linear isotropic material. We suppose that 3 no external sourceJ s is situated within the conducting regions. So, in the external insulating regionR\V c , whereσis zero, we obtain the following field equations: ∇ ×H=J s ,∇ ∙B=0, B=µH.(2) The problem is well posed by imposing regularity conditions at infinity and suitable interface conditions on∂ V c . In particular,[H] c ×n c =0,[B] c ∙n c =0, wheren c is the outer normal to∂ V c and[v] c stands for the jump of the quantityvat∂ V c . Clearly, such a type of interface conditions has also to be verified at any surface whereσorµis discontinuous. Additionally to the boundary conditions, we have to impose suitable initial values for the vector fields at a given timet 0 . In particular, the initial condition onBhas to give∇ ∙B=0 and[B] ∙n=0 at any interface. We point out the fact that the vector fieldsJand∂ t B are automatically forced to be soleinodal by Eqs. (1). Of course, the condition∇ ∙B=0 is satisfied at any time provided that it is verified by the initial condition. By introducing artificial boundary conditions, we can work on a bounded domainV. For simplicity, we assume thatV c is a simply connected polyhedral subdomain ofVandV c ⊂V. In a weak form, equations∇ ∙B=0 andB=µHinVcan be written as follows:  1 µH∇v=0,∀v∈H(V ),(3) 0 V

934

To cite this article: Y. Maday et al., C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938

1 22 3 where H(V )= {v∈L(V )| ∇v∈(L, v(V )) |∂ V =0}. From Eqs. (1) inV c , we can write 0   ∇ ×H∇ ×W+σ ∂ t (µH )W=0,∀W∈H 0 (curl, V c ),(4) V c V c where     3 3 2 2 H 0 (curl, V c )=W∈L(V c )| ∇ ×W∈L(V c W) ,×n c =0.(5)  For the current densityJ, the condition∇ ∙J=0 suggests the introduction of a vector potentialT   such thatJ= ∇ ×T .Then, inV c , the difference between the vector potentialTand the magnetic    fieldHcan be written as the gradient of a scalar function, i.e.,H=T− ∇.A similar argument holds for the insulating region where we assume knowing a vector potentialT s such thatJ s = ∇ ×T s .   Combining external and conducting regions we writeHinVasH=T− ∇,inV c andH=T s − ∇ inV\V c . By eliminating the magnetic fieldHin (3) and (4), we obtain a coupled eddy current problem in terms   of the electric potentialTdefined only in the conducting regionV c and the scalar potentialdefined everywhere inV. This system is completed with appropriate interface conditions over∂ V c stating, e.g.,     thatis continuous. This is nevertheless not enough to defineandTuniquely. In fact∇ ∙Tis not  specified and thus there are many different gauge possibilities. One of them is to require thatThas  the same divergence asHinV c but this eliminatesonV c . We prefer another condition, stated in Section 2.

2. Variationalproblem

  In this section, we define a variational formulation based on the decomposition ofHintoTand∇. We rightaway consider the system obtained after discretization in time of (3) and (4). Only implicit time discretization schemes can satisfy the stability condition. Having discretized the time derivative by means of a finite difference scheme of time stepδt, then at each time step, we have to face a boundary value problem of the form   1   ∇∇v−T∇v=f∇v,∀v∈H(V ), 0 V V c V\V c   (6)   α∇ ×T∇ ×W+T W− ∇W=f c W,∀W∈H 0 (curl, V c ), V c V c V c

where the coefficientα >0 is constant. In the more general approach, it is uniformly positive definite and depends on the material parametersσ,µas well as on the time step (e.g.,α=δt /µσinV c ). Note that the   unknownsTanddenote the approximation at the current timestep,f c depends on the approximations   ofTandat the previous timestep, andfdenotes the scaled source term depending onT s . Additionally,   Tandhave to satisfy the interface conditions at each time step: we choose hereT∈H 0 (curl;V c )and 1    ∈H(V ). In our approach, the strong coupling betweenTandat the interface will be replaced 0 1   by a weaker one. We consider the following discrete problem: find(T , )∈H 0 (curl;V c )×H(V )such 0 that     1  ˆ  a ,v+b c T , v=f∇v,∀v∈H(V ), 0 V\V c  (7)     ˆ  a c T , W+b c W, =f c W,∀W∈H 0 (curl;V c ). V c

935

Y. Maday et al. / C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938

Here, the bilinear forms are defined by       a c T , W:= α∇ ×T∇ ×W+,T W∀T , W∈H 0 (curl;V c ), V c  1 ˆ H (8) b c (W, v):= −W∇v,∀W∈H 0 (curl;V c ),∀v∈ 0 (V ), V c    1   a ,v:= ∇∇v,∀, v∈H(V ). 0 V The bilinear forms are continuous with respect to the associated norms. 1    Now, it is easy to see that if(T , )is a solution of (7), then(T+ ∇φ, +φ),φ∈H(V c ), is a solution 0  as well. To achieve uniqueness, we chooseφsuch that=+φis harmonic onV c ; to this purpose, we 1 1 introduce the harmonic extension operatorH:H(V c )→H(V c )verifying  1 H:=v | ,∇Hv∇w=0,∀w∈H(V c ),(9) v | ∂ Vc∂ Vc 0 V c 1 and state the modified variational problem: find(T , )∈H 0 (curl;V c )×H(V )such that 0  1 a(, v)+b c (T , v)=f∇v,∀v∈H(V ), 0 V\V c  (10) a c (T , W )+b c (W, )=f c W,∀W∈H 0 (curl;V c ), V c  1 whereb c (W, v):= −W∇Hv, forv∈H(V )andW∈H 0 (curl;V c ).In the following, ∙  0;D stands 0 V c 2 for the Lnorm on the open bounded setD. L 2.1.–There exists a constant0< γ 1 <1depending only onV c andVsuch that EMMA 1 ∇Hv 0;V c γ 1 ∇v 0;V ,∀v∈H(V ).(11) 0 For homogeneous Dirichlet boundary conditions, Lemma 2.1 is sufficient to establish the ellipticity ofa g ((W, w), (V , v)):=a c (W, V )+b c (W, v)+b c (V , w)+a(w, v),whereV , W∈H 0 (curl;V c )and 1 v, w∈H(V ), which is the bilinear form associated with the variational problem (10). 0 1 )×H. L EMMA 2.2. –The bilinear forma g (∙,∙)is continuous and elliptic onH 0 (curl;V c 0 (V ) 1 Proof. –We start by consideringa g ((W, w), (W, w)), for(W, w)∈H 0 (curl;V c )×H(V ), in more detail. 0 With the Hilbert spacenorm H 0 (curl;V c ), we associate the energy| ∙ | V c defined by 2 2 2   W + c ). |W| := 0;V c α∇ ×W 0 ,∀W∈H 0 (curl;V ;V c V c 2 This norm is equivalent to the standard Hilbert space norm. The Lorthogonality ofHwyields  2 2  )= |W| +∇w−2W∇Hw a g ((W, w), (W, w) V c  0;V V c 2 2  w W∇ + 0 −2W 0;V ∇Hw. V c ;V c 0;V c Now, we can use (11) to bound∇Hwn terms of∇  0;V c iw 0;V and we get   2 2    a (W,w), (W, w) ∇w g W + 0;V −2γ 1 W ∇w 0;V V c V c   2 2   (1−γ 1 )W∇ +w. 0;V V c The last inequality shows that the ellipticity constant only depends on the constant in (11) and that to obtain a “good” ellipticity constant, it is sufficient thatγ 1 is not too close to one.✷ The unique solvability of the variational problem (10) is a consequence of Lemma 2.2 and of the continuity of the bilinear forms. The first equation in (10) and the definition of the harmonic extension

936

Pour citer cet article : Y. Maday et al., C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938

1 1 yieldsb c (T , v)=0,∀v∈H(V c ). Moreover, takingW= ∇v, for a functionv∈H(V c ), in the second 0 0 equation of (10), we get thatTis divergence free as soon asf c is divergence free. Hence,Tis implicitly gauged, andrestricted toV c is harmonic.

Remark1. – The constantγ 1 does not depend on the shape regularity ofV\V c andV c . In particular for time dependent problems whereV c may be moving, this might be important. More precisely,γ 1 does not depend on the distance between∂ V c and∂ V. In the limit case that this distance tends to zero,γ 1 is still bounded away from one. In this situation we work with a ficticious domainV g such thatV⊂V g andV g \V 1 is shape regular. Due to the homogeneous Dirichlet boundary conditions ofv∈H(V )we extendvby zero 0 1 to an element in H(V g ). Moreover, the ellipticity constant does not degenerate if the distance between∂ V c 0 and∂ Vtends to zero.

3. Discretization

For the discretization of the vector fieldT, we propose to use the lowest order Nédélec finite elements, also known as edge elements [4]. We denoteX 0;h (V c )the associated finite element space. For a given quasi uniform simplicial triangulationT h ofV c , the finite element functions are curlconforming and the degrees of freedom are vector circulations over the edgeseofT h . We note thatT∈X 0;h (V c )has a vanishing tangential component on∂ V c . Thus the degrees of freedom ofX 0;h (V c )are associated with the interior edges ofT h . The domainVis associated with a second quasiuniform simplicial triangulationT H . The discretization ofonVis based on standard conforming elements of the lowest order and we denote by S 0;H (V )the associated finite element space. The finite element functions are now piecewise linear and 1 H conformingand the degrees of freedom are associated with the interior vertices of the triangulationT H . We denoteS h (V c )the space of standard conforming elements of first order associated withT h onV c . We remark that no boundary conditions are imposed onS h (V c ). Moreover, we assume that∂ V c can be written as union of faces inT h . The trace space ofS h (V c )on∂ V c is calledW h (∂ V c ). To formulate the discrete version of the variational problem (10), we have to specify a discrete operator replacing the harmonic extension (9). A natural choice is to involve the discrete harmonic extensionH h defined as a mapH h :W h (∂ V c )→S h (V c )verifying  1 H h v:=v | ,∇H h v∇w=0,∀w∈S h (V c )∩H(V c ).(12) | ∂ Vc∂ Vc 0 V c However in general, the restriction ofv∈S 0;H (V )on∂ V c is not an element inW h (∂ V c ). Thus, we cannot apply directly the discrete harmonic extension to the restriction ofv∈S 0;H (V )on∂ V c . To overcome this difficulty, we introduce a projection operator on the boundary. This operator is well known in the mortar finite element context [2] and can be defined in terms of a Lagrange multiplier spaceM h (∂ V c ):   2  h :L(∂ V c )−→W h (∂ V c ),  h vν=vν,∀ν∈M h (∂ V c ).(13) ∂ V c ∂ V c To obtain a well defined operator h , the Lagrange multiplier spaceM h (∂ V c )has to be well chosen. There are many possibilities but, for simplicity reasons, we restrict ourselves to two choices. In the first case, M h (∂ V c ):=W h (∂ V c ), as in [2]. In contrast to mortar finite element methods with many subdomains, no modification ofW h (∂ V c )is required since∂ V c is a closed surface without crosspoints. Then, the 2 operator h is a Lprojection and, numerically, a mass matrix system has to be solved. In a second case, we replace the standard hat functions by piecewise linear dual basis functions [5] in the definition ofM h (∂ V c ). 2 Then, h projection having the same qualitative stability properties as before. The advantageis a quasi L s is that the mass matrix system is diagonal. These choices guarantee that the operator h is H stable for 1/2−1/2 0snorm and in the H1. Furthermore, it satisfies the approximation property in the Hnorm. In terms of the discrete harmonic extension and the operator h , we formulate the new discrete variational problem: find(T h ,  H )∈X 0;h (V c )×S 0;H (V )such that

937

Y. Maday et al. / C. R. Acad. Sci. Paris, Ser. I 334 (2002) 933–938  a( H , v)+b h (T h , v)=f∇v,∀v∈S 0;H (V ), V\V c  (14) a c (T h , W )+b h (W,  H )=f c W,∀W∈X 0;h (V c ), V c  whereb h (W, v):= −W∇H h  h v, forv∈S 0;H (V )andW∈X 0;h (V c ).Lemma 2.1 has an analogous V c discrete version that yields to prove an optimal error estimate (see[3] for further details). L EMMA 3.1. –Forh/Hsmall enough, problem(14)has a unique solution and there exists a constantC 1 31 3 independent of the meshsize and a realβ2such that, forT∈(H(V c ))with∇ ×T∈(H(V c ))and β ∈H(V ), we have    β−1   −T +−T ∇×T +H.(15) T h H 1;V C h 1;V c +  1;V c β;V V c 4. Algorithmics We conclude our analysis by presenting a numerical algorithm to solve the discrete problem (14). Let us denote byA c andAthe standard stiffness matrices associated with the bilinear formsa c (∙,∙)anda(∙,∙) respectively. LetQandSbe the matrices associated, the first, with the mortar projection h fromS 0;H (V ) toW h (∂ V c )and, the second, with the harmonic extensionH h fromW h (∂ V c )toS h (V c ). The algebraic form 1 of the discrete problem (14) reads:find two vectorsT h and H solution of the linear system t tt A c T h +BSQ H =F c , A H +Q S B T h =F,(16) t whereBis a rectangular stiffness matrix associated with the bilinear formb h (T h ,  H )and denotesthe transpose operator. As a linear iterative solver for (16), we suggest the use of a block Gauss–Seidel method, n n+1n+1 i.e., starting from, first we computeTand thenby means of H hH n+1n n+1 t t t n+1 A c T+BSQ=F c , A+Q S B T=F.(17) h H H h This algorithm converges thanks to the following lemma, whose proof relies on Lemma 2.1 and on the continuity and coerciveness of the bilinear formsa(∙,∙)anda c (∙,∙). n n EMMA tor .Th L 4.1.–Let us denote byethe error H −at the iterationnon the vec H e mapping H n n+1 e→eis a strict contraction overS H (V ):there exists a constant0<< θ1such that    n+1n+1n n a e, e< θa e ,e . n n t t The convergence ofto H yields the one ofTtoT h . In (16), the application ofB,Q,QandB H h t is standard. We do not assembleSandSbut solve two homogeneous Dirichlet boundary value problems. Thus, at each stepn, we have to solve two Dirichlet boundary value problems onS h (V c ), one curl problem onV c and one Dirichlet problem onS 0;H (V ).

1 In what follows, we use the same lettersTandfor the discrete functions and their vectors of unknowns in the h H appropriated respective basis.

References [1] R. Albanese, G. Rubinacci, Fomulation of the eddycurrent problem, IEE Proc. A 137 (1990) 16–22. [2] C. Bernardi, Y. Maday, A. Patera, A new nonconforming approach to domain decomposition: The mortar element method, in: H. Brezis, J. Lions (Eds.), Nonlinear Partial Differential Equations and Their Applications, Pitman, 1994, pp. 13–51. [3] Y. Maday, F. Rapetti, I.B. Wohlmuth, Mortar element coupling between global scalar and local vector potentials to solve eddy current problems, Enumath Conference Proceedings, Springer, 2002, to appear. 3 [4] J.C. Nédélec, Mixed finite elements inR, Numer. Math. 35 (1980) 315–341. [5] I.B. Wohlmuth, Discretization Methods and Iterative Solvers Based on Domain Decomposition, Lecture Notes in Comput. Sci. and Engrg., Vol. 17, Springer, 2001.

938