We describe a preconditioner for the Kirchhoff plate equation 
for use of Bogner-Fox-Schmidt finite elements based on a 
hierarchical technique

Meyer, Arnd

English

Multimedia ONline ARchiv CHemnitz

Arnd MeyerHierarchical Preconditioners andAdaptivity for Kirchhoff–PlatesCSC/08-03Chemnitz Scientific ComputingPreprintsImpressum:Chemnitz Scientific Computing Preprints — ISSN 1864-0087(1995–2005: Preprintreihe des Chemnitzer SFB393)Herausgeber:Professuren fu¨rNumerische und Angewandte Mathematikan der Fakulta¨t fu¨r Mathematikder Technischen Universita¨t ChemnitzPostanschrift:TU Chemnitz, Fakulta¨t fu¨r Mathematik09107 ChemnitzSitz:Reichenhainer Str. 41, 09126 Chemnitzhttp://www.tu-chemnitz.de/mathematik/csc/Chemnitz Scientific ComputingPreprintsArnd MeyerHierarchical Preconditioners andAdaptivity for Kirchhoff–PlatesCSC/08-03CSC/08-03 ISSN 1864-0087 Sept. 2008Contents1 Introduction 12 The Kirchhoff–Plate Equation 13 Quick Solvers for the Stiffness Matrix 43.1 The Multiply with Q . . . . . . . . . . . . . . . . . . . . . . . . . 53.2 The Scaling Diagonal . . . . . . . . . . . . . . . . . . . . . . . . . 74 Numerical Example 85 Conclusion 8Author’s address:Arnd MeyerTU ChemnitzFakulta¨t fu¨r MathematikD-09107 Chemnitzhttp://www.tu-chemnitz.de/mathematik/1 IntroductionOur aim is the development of efficient numerical tools for calculating deforma-tions of thin structures as plates or shells considering finite (large) strain equa-tions. For the beginning we start with quick solving the finite element equationsof the well-known linear Kirchhoff-plate equation for the small strain case. Here,the main challenge is the fourth order biharmonic equation, that requires bothspecial finite elements (with C(1)-functions) and a special preconditioner for theresulting matrix equation.The paper is organized as follows. We start with a short sketch of deriving theKirchhoff plate equation and its conformal finite element discretization in Chapter2. In Chapter 3 we describe the hierarchical preconditioner for the PCGM solverof the discretized problem.2 The Kirchhoff–Plate EquationWe consider the deformation of a 3D-domain Ω3D = Ω× [−d2, d2] with thickness din z-direction and (x, y) ∈ Ω ⊂ R2 a plane. Here, d << diamΩ, so the well-knownKirchhoff hypothesis about an approximation to the 3D-displacement vector readsas:U(x, y, z) = −z wx−z wyw (1)with the unknown scalar function w(x, y) ∈ H2(Ω). The general 3D-solution of alinear elastic deformation of Ω3D reads asa(U ,V ) = f(V ) ∀V ∈ (H1D(Ω3D))3 (2)witha(U ,V ) =∫Ω3Dσ(U) : ǫ(V ) dΩdz (3)f(V ) =∫Ω3D~f · V dΩdz +∫Γ3DN~g · V dS (4)for given volume forces ~f and boundary tractions ~g on parts Γ3DN of ∂Ω3D. Inlinear elasticity we defineǫ(U ) =12(GradU +GradUT)(5)1andσ(U) = 2µ ǫ(U) + λ(tr ǫ)I (6)with the Lame´ constants µ and λ. Using Cartesian coordinates this meansǫij(U) =12(∂iU j + ∂jU i) (7)soa3D(U ,V ) =∫Ω3D2µ ǫij(U) · ǫij(V ) + λ divU divV dΩdz (8)The Sobolev space H1D(Ω3D) contains functions u ∈ H1(Ω3D) with zero values atthe Dirichlet boundary Γ3DD . We assume thatΓ3DD = ∂Ω ×[−d2,d2]andΓ3DN = Ω×{±d2}(hence: dS = dΩ).Under the restriction ~f = f(x, y)e3 and ~g = g(x, y)e3 with scalar functions fand g we can project the problem (2) into the subspace due to (1) and integrationover z ∈[−d2, d2]leads to the well-known Kirchhoff plate equationa(w, v) = f(v) ∀v ∈ H20(Ω) (9)witha(w, v) =d312∫Ω2µ(∇∇T w) : (∇∇T v) + λ(∆w)(∆v) dΩ,f(v) = d∫Ωf v dΩ +∫Ωg v dΩ.For the finite element equation it is convenient to defineǫ(w) = wxxwyy2wxy = D(∇)w (10)with the differential operatorD(∇) = (∂2∂x2,∂2∂y2, 2∂2∂x∂y)Tand the material matrixC = 2µ+ λ λ 0λ 2µ+ λ 00 0 µ ,2such thata(w, v) =d312∫Ω(Dv)TC(Dw) dΩ (11)From this abbreviation the structure of a finite element stiffness matrix is simplydeduced. Let (ϕ1, ..., ϕn) be the basis of the n finite element ansatz functions,thenD̂ = (Dϕ1... .......Dϕn)defines the stiffness matrixK =d312∫ΩD̂TCD̂ dΩ.Analogously we can use (ϕ1, ..., ϕn) as the n form functions of one element (witha much smaller n), thenKel =d312∫elD̂TCD̂ dΩ (12)is the element stiffness matrix.The appropriate finite element discretization of the equation (11) requires H2-functions for the global ansatz functions ϕi. This can be achieved with the so–called Bogner-Fox-Schmidt elements. Here, we assume Ω consisting of a unionof axis-parallel rectangles, then the 16 cubic form functions ϕabkl (x, y) (a, b, k, l ∈{0, 1}) can fulfill these requirements, if they are defined with the help of the1D–functions:pˆ00(s) =14(1− s)2(2 + s) pˆ01(s) = pˆ00(−s)pˆ10(s) =14(1− s)2(1 + s) pˆ11(s) = −pˆ10(−s)for s ∈ [−1, 1](13)After transformation of [−1, 1] onto the x-interval of the element we haveP = (p00(x), p01(x), p10(x), p11(x)) and after transformation of [−1, 1] onto the y-interval we have Q = (q00(y), q01(y), q10(y), q11(y)). Thenϕabkl (x, y) = pak(x) qbl (y) (14)or the basis Φ = (ϕ0000 , . . . , ϕ1111) can be written as the Kronecker–productQ(y)⊗P (x) of the one dimensional functions. This will be used later for derivingthe refinement formula of the preconditioner. The calculation of the elementmatrix is easily done from (12) usingD̂ = Q ⊗ ∂xx P∂yy Q ⊗ P∂y Q ⊗ ∂x P (15)3At each node of the finite element mesh, we have 4 ansatz functions with supportof the 4 adjacent elements. The dimension of the linear system is four timesthe number of nodes. The entries of the solution vector can be partitioned into4-dimensional sub-vectors belonging to each nodewnode = (w00node, w10node, w01node, w11node),which have the following meaning (due to the properties of the ansatz functions)wabnode =∂a∂xa∂b∂ybw|node (16)3 Quick Solvers for the Stiffness MatrixFor efficient solving the discretized finite element linear equations, we use the pre-conditioned conjugate gradient method (PCGM). Then the number of arithmeticoperations is proportional to the number of unknowns for each step of the itera-tion. It remains to construct a preconditioner that guarantees a small number ofiterations (until a given accuracy) that does not grow after mesh refinement. In[2, 4, 5] these preconditioners are derived as Additive Schwarz Method similar tothe BPX-preconditioner[1]. Here, we consider an easier approach based on thehierarchical basis. Even with a very slow growth of the number of iterations (aslogN), this would be an efficient solver from two reasons.First, the implementation is very cheap and second we will embed it into theadaptive mesh refinement, where the so-called cascadic effect dampens the re-quired iteration numbers additionally.The idea behind is the definition of the hierarchical basis of the finite elementspace VN of the ansatz functions. If Φ = (ϕ1, . . . , ϕN) is the original (nodal)basis in VN and Ψ = (ψ1, . . . , ψN ) the hierarchical one, then it exists an(N ×N)-matrix Q that transforms one into another,Ψ = Φ Q.From the fact that the stiffness matrix KH belonging to Ψ is much better condi-tioned than our usual stiffness matrix K belonging to Φ, we arrive at a precondi-tioner.If λi(KH) denote the eigenvalues of KH , thenλi(KH) = λi(QTKQ) = λi(QQTK),so C withC−1 = QQT4is a “good” preconditioner for K. This was first discovered by Yserentant for theLaplace equation on linear triangular finite elements [7]. The multiply with the(theoretically dense) matrices Q and QT is done in O(N) operations only fromthe hierarchical definition of the mesh. Analogously to this 2nd order equation,we can derive a hierarchical refinement formula for our cubic functions as well,which is used for the preconditioning action.3.1 The Multiply with QLetΦcoarsk = (ϕ00k , ϕ10k , ϕ01k , ϕ11k )be the 4 ansatz functions belonging to the mesh point k of a fixed mesh withmesh spacings h+x , h−x in x-direction and h+y , h−y in y-direction around the k-thnode. After mesh refinement as subdivision of all elements around this node, wehaveΦF inek = (f00k , f10k , f01k , f11k )the new (fine) functions around k-th node with half support in both directionsand 8 new quadruples of functions belonging to the 8 new nodes around thisnode. These are ΦF inej = (f00j , ... f11j ) with 8 new indices j ∈ Nk (the set of allnew neighbors of k). The refinement formulaΦcoarsk = ΦF inek +∑j∈NkΦF inej Aj (17)requires the (4x4)-matrices Aj , which can be deduced from the 1-dimensionalcase, following [4].We need to consider the refinement formula of the 2 functionsPCa (x) = (p0(x), p1(x)),defined on an interval [a, b] belonging to the left end point a. (So:PCa (a) = (1, 0), PCa (b) = (0, 0),∂xPCa (a) = (0, 1), ∂xPCa (b) = (0, 0),see Fig. 1)The subdivision of this interval is done by introducing c = a+b2as new grid point.Now we have the fine grid functionP Fa (x) = (q0(x), q1(x))5 0 0.2 0.4 0.6 0.8 1 a  c  bp0p1Figure 1: PCa (x)with support in [a, c] and a new fine function F (x) = (f 0(x), f 1(x)) with supportin [a, b] andF (c) = (1, 0), F (a) = F (b) = (0, 0)∂xF (c) = (0, 1), ∂xF (a) = ∂xF (b) = (0, 0)acording to Fig.2 and 3. 0 0.2 0.4 0.6 0.8 1 a  c  bq0q1 0 0.2 0.4 0.6 0.8 1 a  c  bf0f1Figure 2: P Fa (x) Figure 3: F (x)Then we easily calculatePCa (x) = PFa (x) + F (x) · A(+h)forA(h) =(1/2 h/8−3/(2h) −1/4)and h = b− a. (18)If the coarse mesh function PCb belongs to the right point b (instead of a) thesame formula with −h is correct:PCb (x) = PFb (x) + F (x) ·A(−h). (19)6Now, we can derive the 2D-refinement formula from the Kronecker-product defi-nition of Φcoarsk :Φcoarsk = QCk ⊗ PCk= (QFk + FF1 ·A(−hy) + FF2 · A(+hy))⊗ (P Fk + FF3 · A(−hx) + FF4 · A(+hx))= QFk ⊗ PFk + (QFk ⊗ FF3 )(I ⊗ A(−hx)) + ...(20)which is exactly the sum (17) with the matrices Aj asI ⊗ A(±hx), A(±hy)⊗ I, or A(±hy)⊗A(±hx)depending on the place of the 8 nodes j relative to the center node k. Theimplementation is very simple, if the data structure “Edge”, as the tree of alledges (from coarse to fine) is defined with 3 informations on each edge: “left”node number (l), “right” node number (r), number of mid-node (m). From thedefinitions above the following algorithm for the multiply w := Qy arises. Notethatw = (wT1 , ... , wTn )Twith wi ∈ R4 for n nodes of the actual mesh (and N = 4n).Algorithm: Multiply w := Qy1. w := y2. for each edge with nodes (l, r,m) from coarse to finest dohx = (x-coord. of node r)− (x-coord. of node l)hy = (y-coord. of node r)− (y-coord. of node l)if (|hy| > |hx|) then B+ := A(+hy)⊗ IB− := A(−hy)⊗ Ielse B+ := I ⊗A(+hx)B− := I ⊗ A(−hx)wm := wm +B+wr +B−wl(21)Note that hx and hy have a sign depending on the position of the endpoints ofthis edge and one of them is zero due to the axis-parallel mesh.3.2 The Scaling DiagonalIn the case of 2D calculations on 2nd order p.d.e.’s, the size of the elements of thestiffness matrix is fixed independent of the mesh size. From this reason we do not7need a scaling diagonal D in the definition of the preconditioner C−1 = QDQTin the original paper of Yserentant.Here, this is no longer true, because we consider a 4th order p.d.e. and the 4different degrees of freedom per node are of different kind (the value and deriva-tives, which scale differently w.r.t. h).The best way out is the definition of D inthe MDS (multiple diagonal scaling [6]) style, which means that D is the diagonalof the actual hierarchical stiffness matrix KH .This information can be derived and stored during the mesh refinement steps:We mark a node as “used” if its diagonal entry of a previous stiffness matrix hasbeen stored in D. Then after mesh refinement, we calculate the diagonal D˜ ofthe actual stiffness matrix and copy all entries of D˜ belonging to new nodes (not“used”) into D. With this simple procedure we have always the correct scalingfactor in D.4 Numerical ExampleWe demonstrate the power of this preconditioner at some simple rectangular plateexamples. We consider Ω = (0, 1.5) × (0, 1) with different boundary conditionsand constant force everywhere. Each of the 4 boundary edges are either fullclamped (“c”, means w = 0 and ∇w =−→0 ) or only softly supported (“s”, meansw = 0 and τ · ∇w = 0 with the tangential direction τ).The example “cccs” has soft support on one of the longer edges, “ccss” at twoadjacent edges and “cscs” at two opposite (longer) edges. The coarse mesh con-sists of one rectangle only, we start the calculation after 2 refinements with thegiven number of elements in the Table 1. The PCGM iteration starts on each newmesh with interpolated values from the next coarser mesh. The relative decreaseof the test value (r, w) under 10−8 has used as stopping criterion. The followingTable 1 contains the number of iterations required until this stopping criterion isfulfilled for these four examples.5 ConclusionThe most simple hierarchical technique for preconditioning the Kirchhoff–plateequation leads to a very efficient solver. It is necessary to consider the correctrefinement formula for the Bogner–Fox–Schmidt element functions as well asthe correct scaling diagonal as the main ingreedients for this preconditioner. In8# elem. # d.o.f.’s. #it. #it. #it. #it.“cccc” “cccs” “cscs” “ccss”16 100 8 12 12 2764 324 18 26 24 29256 1 156 25 27 26 261024 4 356 22 23 24 224096 16 900 18 19 21 1816 T 66 564 16 16 16 1665 T 264 196 14 15 15 14Table 1: Number of PCGM iterations on 4 examplescombination with the cascadic approach we obtain a small and near fixed numberof PCG–iterations required for some given accuracy.References[1] I. H. Bramble, J. E. Pasciak, J. Xu, Parallel Multilevel Preconditioners,Math. Comp. 55 (1990) 191, 1-22.[2] H.Matthes, Die nichtu¨berlappende Gebietszerlegungsmethode zur Paral-lelisierung und Vorkonditionierung iterativer Verfahren zur Lo¨sung vonPlatten– und Schalenproblemen, Dissertation, TU Chemnitz (1996).[3] P.Oswald, Multilevel Finite Element Approximation: Theory and Applica-tions, Teubner Skripten zur Numerik, B.G.Teubner Stuttgart 1994.[4] M.Thess, Parallel Multilevel Preconditioners for Thin Smooth Shell FiniteElement Analysis, Num.Lin.Alg.with Appl. 5 (1999) 5 ,401-440(1999)[5] M.Thess, Parallel Multilevel Preconditioners for Thin Smooth Shell Prob-lems, Dissertation, TU Chemnitz (1998).[6] J. Xu, Iterative Methods by Space Decomposition and Subspace Correction,SIAM Rev. 34 (1992) 4, 163-184.[7] H. Yserentant, Two Preconditioners Based on the Multilevel Splitting ofFinite Element Spaces, Numer. Math. 58 (1990) 163-184.9Some titles of the former CSC preprint series:06-01 T. Eibner, J. M. Melenk. p-FEM quadrature error analysis on tetrahedra. Oc-tober 2006.06-02 P. Benner, H. Faßbender. On the solution of the rational matrix equation X =Q+ LX−1L⊤. September 2006.06-03 P. Benner, H. Mena, J. Saak. On the Parameter Selection Problem in theNewton-ADI Iteration for Large Scale Riccati Equations. October 2006.06-04 J. M. Bad´ıa, P. Benner, R. Mayo, E. S. Quintana-Ort´ı, G. Quintana-Ort´ı, A.Remo´n. Balanced Truncation Model Reduction of Large and Sparse GeneralizedLinear Systems. November 2006.07-01 U. Baur, P. Benner. Gramian-Based Model Reduction for Data-Sparse Systems.February 2007.07-02 A. Meyer. Grundgleichungen und adaptive Finite-Elemente-Simulation bei GroßenDeformationen. Februar 2007.07-03 P. Steinhorst. Rotationssymmetrie fu¨r piezoelektrische Probleme. Februar 2007.07-04 S. Beuchler, T. Eibner, U. Langer. Primal and Dual Interface ConcentratedIterative Substructuring Methods. April 2007.07-05 T. Hein, M. Meyer. Simultane Identifikation voneinander unabha¨ngiger Materi-alparameter - numerische Studien. Juni 2007.07-06 A. Bucher, U.-J. Go¨rke, P. Steinhorst, R. Kreißig, A. Meyer. Ein Beitrag zuradaptiven gemischten Finite-Elemente-Formulierung der nahezu inkompressiblenElastizita¨t bei großen Verzerrungen. September 2007.07-07 U.-J. Go¨rke, A. Bucher, R. Kreißig Zur Numerik der inversen Aufgabe fu¨r gemis-chte (u/p) Formulierungen am Beispiel der nahezu inkompressiblen Elastizita¨tbei großen Verzerrungen. October 2007.07-08 A. Meyer, P. Steinhorst. Betrachtungen zur Spektrala¨quivalenz fu¨r das Schurkom-plement im Bramble-Pasciak-CG bei piezoelektrischen Problemen. Oktober 2007.07-09 T. Hein, M. Meyer. Identification of material parameters in linear elasticity -some numerical results. November 2007.07-10 T. Hein. On solving implicitly defined inverse problems by SQP-approaches.December 2007.The complete list of SFB393 preprints is available viahttp://www.tu-chemnitz.de/sfb393/preprints.html.Chemnitz Scientific Computing Preprints – ISSN 1864-0087

Hierarchical Preconditioners and Adaptivity for Kirchhoff-Plates

HSSS - Hochschulschriftenserver der SLUB

Qucosa

https://monarch.qucosa.de/api/qucosa%3A18972/attachment/ATT-0/

Hierarchical Preconditioners and Adaptivity for Kirchhoff-Plates

Abstract

Similar works

Full text

Available Versions

Multimedia ONline ARchiv CHemnitz

HSSS - Hochschulschriftenserver der SLUB

Qucosa