INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION AT ELEMENT LEVEL FOR INCOMPRESSIBLE MEDIUM CHANG-CHUN WU*, LEI YUAN AND TOMONARI FURUKAWA Department of Modern Mechanics, ;niversity of Science and ¹echnology of China, Hefei 230026, China Department of Quantum Engineering and Systems and Science, ;niversity of ¹okyo, 7-3-1 Hongo, Bunkyo-ku, ¹okyo 113, Japan SUMMARY A deviatoric hybrid element approach, in which the deviatoric stress p, the pressure p and the displacement u are independently dealt with as the element variables, is suggested. The present approach is naturally universal for compressible and fully incompressible mediums. Moreover, it can be extended to the simulation of Stokes #ow directly. The resulting hybrid model is able to meet the zero volumetric strain constraint in terms of the incompatible displacement mode only. Therefore an incompressible elimination can be carried out within an individual element, and the complex system elimination for nodal displacements is then avoided. The present 3-"eld hybrid model maintains the important features of current hybrid stress elements*"nally resulting in a set of displacement-type discrete equations which can be easily solved, while not a set of u-p mixed-type equations resulted. Regarding the numerical stability of the element, an e!ective strategy is o!ered to suppress all the zero energy modes hidden in the model. Copyright 1999 John Wiley & Sons, Ltd. KEY WORDS: deviatoric hybrid element; incompressible medium; elimination at element level; zero energy mode 1. INTRODUCTION These years a signi"cant amount of research has been devoted to "nite element methods in incompressible problems (including solids and #uids) for both theoretical implications and engineering background underlying. Since Herrman's initial work was published [1], various of variational principles [2}4] and discrete models have been proposed for incompressible calculations as surveyed in [5]. Displacement-based "nite elements, especially isoparametric elements, have been historically used in incompressible analysis for their simple formulations and standard solution procedures. Displacement locking and instability of the stress solution may, however, arise in nearly * Correspondence to: Chang-Chun Wu, Department of Modern Mechanics, University of Science and Technology of China, Hefei 230026, China. E-mail: ccwu@ustc.edu.ch Contract/grant sponsor: National Natural Science Foundation of China CCC 0029-5981/99/290729}17$17.50 Copyright 1999 John Wiley & Sons, Ltd. Received 4 May 1998 Revised 17 September 1998 730 C.-C. WU, L. YUAN AND T. FURUKAWA incompressible situations. Approaches to overcome these problems are the penalty method with selective integration [6, 7] and the incompatible element method [8, 9], both of which can be applied in nearly incompressible analysis with satisfactory results, although not e!ective for fully incompressible analysis. The fully incompressible medium was "rst dealt with by Needleman and Shih [10], wherein the dependent displacement parameters stemming from the zero volumetric constraint (e"0) were eliminated from the system equations. In this approach, the band character of the global sti!ness matrix is however destroyed, therefore a complicated algorithm is required. Based on the topological analysis of the element, another displacement elimination was proposed by Zhong and Li [11, 12] to incorporate the incompressible condition beforehand. Both the elimination approaches nevertheless run into trouble in the estimation of stress. Meanwhile, the mixed element approximations with pressure and displacement (u}p) as independent "elds (see [1, 3, 7, 13, 14]) are suited to both the compressible and incompressible analyses. The approaches, however, always lead to a non-positive-de"nite element matrix with a zero diagonal block which a!ects the e$ciency for obtaining a solution signi"cantly. Furthermore, in the context of mixed element formulations based on the modi"ed potential energy principle, zero energy modes for stress (ZEM(p)) can be induced by the displacement continuity between elements ([15] for detail), and the stability di$culty mentioned above attracted much attention in the past few years [16]. An e!ective 3-"eld (u}p}e) mixed model was thus given by Simo et al. [17] and has been extended to the "nite deformation plasticity #ow with near incompressibility. In line with these investigations, it has been found that assumed stress hybrid elements do not usually lock and can o!er reasonable stress solutions at a nearly incompressible limit. In these models, stress parameters are eliminated at the element level such that the system of displacement sti!ness equations, which can be solved more easily, can be obtained. In spite of this fully incompressible materials with k"0)5, the #exible matrix H becomes singular, which results in failure for hybrid element models [18, 19]. For the numerical simulation of incompressible Stokes #ow, Tong [20] and Bratianu and Atluri [21] proposed an alternative approach where the stress "eld is split into the deviatoric stress "eld and the hydrostatic pressure "eld both of which are assumed independently in element approximation. Although this gives rise to the avoidance of the singularity of the #exible matrix, pressure variables cannot be eliminated at the element level, and thus the resultant equation system remains mixed. In this paper, a 3-"eld hybrid model, in which the incompatible displacement is introduced reasonably such that the pressure variable can be eliminated at the element level, is formulated. This model is based on an unconditional variational principle for both the compressible and the incompressible mediums, in which the complementary energy of the element is split into deviatoric and volumetric parts to be approximated respectively. The formulation therefore uni"es the compressible and incompressible (including Stokes #ow) calculations, and the band and sparse characters remain in the system equations. The outline of the paper is as follows. In Section 2, a deviatoric variational principle is deduced, based on which a hybrid element is proposed in Section 3. In Section 4, the internal variable elimination at the element level is discussed. As a concrete example, a plane hybrid element model is constructed in Section 5 whereas a technique for suppressing zero energy modes is introduced in Section 6. Numerical examples are presented in Section 7. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 731 2. DEVIATORIC VARIATIONAL PRINCIPLE The Hellinger}Reissner functional, in the case of small deformations and linear elastic materials, is given by 4 [pGH eGH!B(pGH)] dv!l.t. n (p , u )" 0 GH G (1) where e " (u #u ), B(p ) and l.t. are the strain, the complementary energy density, and the GH G H H G GH loading term, respectively. Meanwhile, the deviatoric stress and the deviatoric strain can be each expressed as p "p !pd (2) GH GH GH (3) e "e ! ed GH GH GH where p" p is the hydrostatic pressure, and e"u is the volume strain. Substitution of GG GG equations (2) and (3) into equation (1) yields a 3-"eld functional 4 [pGH eGH#pe!BQ (pGH)!BT (p)] dv!l.t. n(p , p, u )" GH G (4) where B (p ) and B (p) are the deviatoric complementary energy function and the volumetric one, Q GH T respectively, in the form of 1 1 B (p )" p p , G"shear modulus Q GH 2 2G GH GH 11 E B (p)" p, k" "bulk modulus T 2k 3(1!2k) (5) (6) It is easily certi"ed that if the prescribed displacement boundary condition is satis"ed, the stationary condition dn(p , p, u )"0 will lead to the equilibrium equation GH G p #p "0 in < (7) GHH G and the constitutive relationships: p e" in < k (8) 1 e " p in < GH 2G GH (9) In addition, a traction boundary condition can also be obtained. As for general compressible materials, the Hellinger}Reissner principle in relation to functional (1) is an unconditional variational problem [22]. The principle for incompressible materials is however converted into a conditional variational problem with the zero volumetric constraint e"u "0 (10) GG which causes an enormous amount of di$culty in unifying the compressible and the incompressible calculations. On the other hand, as an Euler equation, the e}p relation (8) derived from Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 732 C.-C. WU, L. YUAN AND T. FURUKAWA functional n(p , p, u ) (4) may naturally degenerate into zero volume constraint (10) when the GH G medium becomes incompressible. Being independent of compressibility, functional (4) is, hence, always an unconditional variational, which makes it possible to develop a universal numerical approach for both the compressible and incompressible analyses. The above deviatoric variational principle can also be used to solve the incompressible Stokes #ow, wherein u , p and p mean the velocity "eld, the deviatoric stress and the hydrostatic G GH pressure, respectively, while constant G in equation (5) represents the #uid viscosity. Stokes #ow allows the volume complementary energy to vanish, and energy functional (4) is therefore simpli"ed as 4 [pGH eGH#pe!B (pGH )] dv!l.t. n(p , p, u )" GH G (4b) 3. 3-FIELD HYBRID ELEMENT FORMULATION Similar to conventional hybrid element formulations, the element stress "elds for p and p are GH assumed as p"+p ," b (11) GH @ p" a (12) ? where a and b are parameters for the pressure and the deviatoric stress, respectively. The displacement "eld is assumed as q u"+u ,"uq#u "[Nq N ] G H H k (13) in which u is the conventional compatible displacement interpolation in terms of element nodal O displacement q uq"Nq q (14) On the other hand, the incompatible displacements are de"ned by parameter k within the element u "N k (15) H H and introduction of u may improve the performance of the element e!ectively [23, 24]. As we will H describe later, u plays a key role in the incompressible calculations. Usually the employed H incompatible displacement should satisfy the constant stress Patch Test Condition (PTC) [25}28]. In accordance with equation (13), the deviatoric and the volumetric strains can be, respectively, written as e"+e ,[Bq B ] GH H e"[( ( ] O H Copyright 1999 John Wiley & Sons, Ltd. q k q "e #e O H k (16) (17) Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 733 Substituting these assumed "elds into functional (4), we have q q 1 1 #a2[G G ] ! b2H b! a2H a!qT Q n(a, b, q, k)"b2[G G ] ?O ?H k @@ ?C @O @H k 2 2 (18) where T 2@ [BO BH ] dv [G G ]" 2 [( ( ] dv ?O ?H T ? O H 1 2 dv H " @@ 2G @ @ T 1 H " 2 dv ?? k @ @ T [G G ]" @O @H C C C C In equation (18), Q is the equivalent nodal loads. The stationary condition *n/*b"0 gives rise to q b"H\[G G ] @@ @O @H k (19) the substitution of equation (19) into equation (18) then yields 1 q 2 Kqq KTq 1 q q H #a2[G G ] ! a2H a n(a, q, k)" ?O ?H k ?? Kq K 2 2 k k H HH where the distortion sti!ness matrix is given by Kqq KTqj "[G G ]2 H\ [G G ] @O @H @@ @O @H Kqj K HH By using the stationary condition of n in equation (20) with respect to +q k,, we have Kqq KTqj GTaq Kqj Kjj GTaj Gaq Gaj !Haa q Q k " 0 a 0 (20) (21) (22) For the sake of convenience, let M" Kjj GTaj G aj !Haa (23) For a given compressible or nearly incompressible medium, H is positive de"nite, so that ?? "M"O0. Parameters j and a henceforth turn out that Kqj k "!M\ q a Gaq Copyright 1999 John Wiley & Sons, Ltd. (24) Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 734 C.-C. WU, L. YUAN AND T. FURUKAWA and can be eliminated from equation (22). Conclusively, equation (22) is transformed into the conventional sti!ness equation only with respect to q, i.e. K*qq q"Q (25) where K*qq"Kqq! Kqj 2 Kqj M\ Gaj Gaj (26) In accordance with the formulations, a computational procedure can be suggested as follows: (1) Assemble equation (25), and obtain nodal displacement q. (2) Compute element parameters +k a, by equation (24), and obtain the element pressure by equation (12). (3) Compute parameter b by equation (19), and obtain the element deviatoric stress by equation (11). 4. INCOMPRESSIBLE ELIMINATION AT ELEMENT LEVEL Although the hybrid element formulations in the previous section can be implemented to calculate compressible and nearly incompressible materials, it is still hard to use the formulations for solving fully incompressible problems, in which k"R and H "0, the matrix in equation ?? (23) being of the form K GT HH ?H G 0 ?H As will be shown later, matrix (27) is usually rank-de"cient. M" (27) ¹heorem. If incompatible displacement u passes P¹C, then it must satisfy the element incomH pressible condition in the sense of integration, speci,ed as T eH dv"0 (28) C where e "div u is the incompatible volume strain. H H Proof. Observing the fact that Du "0N div u "0 H H (29) DuH dv"0 N div uH dv"0 (30a) eH dv"0 (i.e. uH passes PTC)N eH dv"0 (30b) accordingly we have or equivalently Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 735 and the theorem is proved. In addition, incompressible condition (28), in view of e "W k in H H equation (17), can also be written as T (H dv"0 (31) C Without loss of generality, the assumed pressure trial function can be expressed in terms of the element's natural co-ordinates as a p" a"[1 m g f . . . ] a ? : (32) In accordance with equations (18) and (31), we have T 2? (H dv" G " ?H C i - 1- - e g m g j g f TC g f g k : h 0 -------------m i e (31) g gg ( dv ***P H j f (H dv f TC g g k :h (33) which shows that the row related to the constant pressure vanishes, the matrix M in equation (27) thereby becoming singular. Consequently, equation (24) will no longer hold, and the elimination of variables +k a, cannot be carried out. An e!ective strategy to avoid this problem is to introduce a small disturbance to the incompatible mode u , as shown in the next section, such that condition (28) is no longer met, and H the parameters +k a, can be eliminated within an element. 5. 4-NODE PLANE HYBRID MODEL Let us consider a 4-node plane hybrid model, the element geometry of which is the same as that of the well-known bilinear isoparametric element Q4. The compatible displacement is then naturally assumed as a standard bilinear interpolation. In order to pass PTC in this case, incompatible displacement u can be generated by using the general incompatible formulation [29], by which H an alternative incompatible model is initially assumed as j j [u v ]"[m!* "g#*" m#g!(m#g)] j j H H j j (34) where 2 J J m! g *" 3 J J Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 736 C.-C. WU, L. YUAN AND T. FURUKAWA Figure 1. Four-node plane hybrid element and the coe$cients J , J and J are de"ned by the Jacobian determinant as follows: "J""J #J m#J g (35) It can be veri"ed that the incompatible displacement u (34) satis"es incompressible condition H (28), and there must be de"ciency in the rank of matrix M in (27). To solve the singular problem, we introduce a small parameter d into the high-order component of equation (34), which gives j j [u v ]"[m!*"g#*"m#g!(m#g) (1#d)] j j H H j j (36) Incompatible function (36) does not meet equation (28) any longer and keeps the rank of matrix G su$cient. Moreover, in order to keep "M"O0, the following condition for matching ?H parameters should also be satis"ed: dim(k)*dim(a) (37) In accordance with inequality (37), a bilinear pressure trial function is adopted here: a : (38) a For the 3-"eld hybrid model based on functional (4), besides condition (37), there exists another parameter matching condition [31] that should be satis"ed to avoid the appearance of the Zero Energy Mode (ZEM), i.e. p"[1 m g mg] dim(b)#dim(a)*dim(q*)#dim(k) (39) in which q* is the vector of element nodal displacements excluding d.o.f. of the rigid body. This condition must be obeyed by the rationally assumed deviatoric stress "eld. However, many stress trial functions, such as the bilinear one, which satisfy condition (39), may still induce some zero energy displacement mode ZEM(u). Indeed, equation (39) is just a necessary condition for avoiding ZEM(u), and this ends up with the problem of how to select a rational deviatoric stress mode. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 737 6. SUPPRESSION OF ZERO ENERGY MODE General discussions for ZEM were given by Pian et al. [30, 15, 31], some of which are still applicable for present hybrid models. An alternatively more e!ective method suggested here to prevent all the possible ZEM(u) stems from the bilinear term C p e dv in functional (4). T GH GH The exclusion of the element rigid body d.o.f. from nodal displacement q yields q*, and deviatoric strain (16) then becomes e*"[B*q B ] H q* k (40) For the purpose of preventing ZEM(u), we de"ne the element deviatoric stresses in the manner of p" b"[B*q B ]b @ H (41) and this automatically satis"es dim(b)"dim(q*)#dim(k), and constrains (37) and (39). Note that the rigid body displacement does not give any contribution to energy functional (4), and in terms of (40) and (41), the term C p e dv in the functional can then be formulated as T GH GH T C T p2e dv" C p2e* dv"b2G* q* k (42) where T [B* BH ]2 [B* BH ] dv G*" q (43) q C Observing the positive-de"nite attribute of strain energy, we have T C e*2e* dv" q* 2 q* G* *0 k k (44) where matrix G* is positive de"nite, and it turns out that the resulting hybrid element must be free from ZEM(u). It is to be noted that e may be replaced with e due to the identity p2 e,p2e in equation (42). With the understanding of ZEM, take a rectangular element as an illustrating example since ZEM(u) is more easily induced by a regular element, for which *"0 in equation (36). Additionally, the small parameter d equation (36) can be ignored since it is independent of our element rank analysis. Thus the element displacements become u"+u v,"uq#u H a b : "[1 m g mg m g m#g!(m#g)] : a b Copyright 1999 John Wiley & Sons, Ltd. (45) Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 738 C.-C. WU, L. YUAN AND T. FURUKAWA where a and b are general parameters. For convenience, let the element natural co-ordinates G G coincide with the physical co-ordinates, such that the element strain can be formulated as i *u e g g *m g g 0 g *v g e"j f" 0 *g g g 0 g *u *v g g # g k *g *m h 1 0 g 2m 0 1!3m 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 m 0 2g 1!3g 0 1 m 0 2g 1!3g 0 1 0 g 2m 0 1!3m i g g g j g g g k a e : gg a g f b g : g g b h (46) Subtracting two zero columns and one of duplicate columns from the matrix of order (3;14), and imposing a linear transform to the matrix, the desired deviatoric stress mode can "nally be obtained: 1 0 0 g m 0 m 0 0 0 0 p" b" 0 1 0 0 0 0 0 m 0 g g @ 0 0 1 m 0 g g g m 0 m b : b (47) So far, all the element trial functions are determined as speci"ed by equations (13), (36), (47) and (38), respectively. The resultant deviatoric hybrid element is able to unify the compressible and the fully incompressible calculations. 7. NUMERICAL EXAMPLES In numerical calculations, the following 4-node plane elements are employed: (1) bilinear isoparametric element with 2;2 Gauss quadrature, denoted as Q4(2;2); (2) modi"ed bilinear element Q4(2;2 and 1;1), in which the deviatoric strain energy and the volumetric one takes 2;2 and 1;1 Gauss quadrature, respectively [7]; (3) enhanced strain incompatible element QUAD suggested by Simo and Rifai [32]; (4) the present 4-node deviatoric hybrid element with 3;3 Gauss quadrature. 7.1. Constant stress patch test Since disturbing parameter d was introduced in incompatible formulation (36), so that the elimination at the element level can be carried out, let us inspect the e!ect of the parameter on the constant stress patch test. A plane stress/strain panel with size 10;10 and unit thickness is shown in Figure 2. The panel was modelled by four irregular elements, for which the pure stretch patch test and the pure shear one were respectively considered. Some of the resultant stress solutions are listed in Table I. It is found that with the range of 10\)d)10\, regardless of the compressibility of the material, the present hybrid element is able to reproduce constant stress states for both p and q with high precision. In the following V VW computations, d was thus taken to be 10\. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 739 Figure 2. Constant stress patch test for plane stress/strain square panel (10;10): (a) pure stretch patch test; (b) pure shear patch test Table I. Patch test for the plane strain panel (Figure 2) by the present hybrid model Stretch stress p V Shear stress q VW d k"0)3 k"0)5 k"0)3 k"0)5 10\ 10\ 10\ 9)9958 9)9996 9)99996 9)9958 9)9996 9)9999 9)9915 9)9992 9)9999 9)9914 9)9992 9)9998 Exact 10 10 10 10 7.2. Numerical stability tests Figure 3 shows a plane strain circular disc with a pair of radial concentrated loads and its "nite element meshes. Three sets of stress solutions, given by Q4(2;2), Q4(2;2 and 1;1) and the proposed deviatoric hybrid element, are all listed in Table II together for comparison. Obviously, the isoparametric element Q4(2;2), provides us with very bad stress solutions for a nearly incompressible material (k"0)49), and moreover no solutions are provided for the fully incompressible one (k"0)5). The selective reduced integration element Q4(2;2 and 1;1) can provide rational stress solutions in the case of k"0)49, but no incompressible solutions either. From Table II, we "nd that only the proposed element can o!er rational stress solutions for compressible, nearly and fully incompressible materials, and exhibit a "ne numerical stability. Another illustrating example is the pure bending of a plane strain cantilever modelled by "ve elements as shown in Figure 4. The results in Table III indicate that some numerical problems, such as the displacement locking phenomenon, the poor stress solutions and no incompressible solutions, etc., cannot be avoided by Q4(2;2). The element Q4(2;2 and 1;1) provides very poor stress solutions too. The proposed element, however, uni"ed the compressible and incompressible calculations and produced excellent bending solutions for both the displacements and the bending stresses. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 740 C.-C. WU, L. YUAN AND T. FURUKAWA Figure 3. A plane strain circular disc with radial concentrated loads Table II. Stress solutions of a circular disc (Figure 3) Q4(2;2) Stress u"0)3 u"0)49 u"0)5 Exact p V p W p V Q4(2;2 and 1;1) p W 5)41 !26)1 !0)88 !12)6 !28)7 !59)2 8)51 !4)04 No solutions 6)37 !26)1 !1)28 !12)9 p V p W p V Present p W p V p W p V p W 6)10 6)10 !25)3 !1)89 !14)0 !25)3 !1)89 !13)5 No solutions 5)98 5)78 5)76 !26)1 !1)58 !13)0 !26)3 !1)88 !13)3 !26)3 !1)99 !13)4 6)37 !26)1 !1)28 !12)9 6)37 !26)1 !1)28 !12)9 Figure 4. Pure bending plane strain beam of "ve elements 7.3. Plastic analysis For plane strain problems, it is often hard for conventional isoparametric elements to "nd a rational plastic limit solution due to the assumption of plastic incompressibility. Let us consider Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 741 DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION Table III. Pure bending solutions of the plane strain beam (Figure 4) u"0)25 Element Q4(2;2) Q4(2;2 and 1;1) Present Exact u"0)499 < p /p V V! < p /p V V! 42)05 73)15 90)54 93)75 !1721/1710 !509/281 !2987/2989 !3000/3000 3)77 87)18 73)60 75)10 2327/4228 !526/268 !2993/2996 !3000/3000 u"0)5 < No No 73)54 75)00 p /p V V! solutions solutions !2994/2997 !3000/3000 Figure 5. Elastic}plastic thick cylinder under inner pressure P: (a) "nite element meshes; (b) loading}displacement curve by the Q4 element; (c) loading}displacement curves by the present element and QUAD [32] an elastic}plastic thick cylinder under inner pressure p. The analytical solution for a plastic limit load is 4802 N/cm. The discrete meshes for a quarter of the cross-sectional part of the cylinder are shown in Figure 5(a), and the usual incremental iterative formulae are adopted. As the plastic solutions, the loading-displacement curves by Q4(2;2) element are drawn in Figure 5(b). It is Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 742 C.-C. WU, L. YUAN AND T. FURUKAWA Table IV. Elastic and plastic limit loads of a thick-walled cylinder (Figure 5) Element Elastic limit load Plastic limit load 3000 2900 2900 2964 No solution 4750 4900 4802 Q4(2;2) QUAD[32] Present Exact Figure 6. Incompressible Couette #ow between in"nite plates: (a) computational model for Couette #ow (mesh 16;4); (b) velocity distributions on the middle section A}A (for the pressure radiant P"!8/!5/0/5/8) Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 743 Figure 7. Lubrication in bearing: (a) #ow in the wedge between the slide block and the guide surface (10;4 meshed); (b) pressure distribution over the slide block; (c) velocity pro"les in bearing shown that there exists no plastic limit load, and the real structural sti!ness is overestimated by the isoparametric element. On the other hand, as shown in Figure 5(c), the present hybrid element and the enhanced strain incompatible element QUAD o!er rational plastic solutions, and the obtained plastic limit loads are consistent with the analytical one. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) 744 C.-C. WU, L. YUAN AND T. FURUKAWA 7.4. Stokes yow analysis It should be pointed out that except the present deviatoric hybrid model, the other elements mentioned above, i.e. Q4(2;2), Q4(2;2 and 1;1) and QUAD, cannot be directly used in Stokes #ow calculation due to the strict incompressible constraint. For Stokes #ow analysis, two examples were o!ered here to verify the e$ciency of the proposed element. First, we considered the stationary Couette #ow between in"nite plates, which was approximated by a "nite #ow "eld (length 10 and high 1) with 16;4 "nite element meshes. As shown in Figure 6(a), the upper surface pulls the #uid #ow with velocity ;"1, and the change of the pressure radiant P on the left-hand side will result in the variation of velocity. Figure 6(b) shows "ve sets of velocity distributions along the middle section A}A, each of which corresponds to a certain given pressure P. As seen with the results, all of the numerical solutions are consistent with the analytical ones. Another Stokes #ow problem was the lubrication in a bearing. The essential features of this type of motion can be understood with the example of a slide block moving on a plane guide surface, (Figure 7(a)). In order to obtain a steady-state problem, let us assume that the slide block is at rest and that the plane guide is forced to move with a constant velocity ; with respect to it. Pressure P is then enforced at two sides of the bearing. A set of 10;4 discrete meshes were employed to simulate the #ow in the wedge between the slide block and the guide surface. The obtained pressure distribution over the slide block is drawn in Figure 7(b), and the related velocity pro"les are shown in Figure 7(c). All these numerical results demonstrate the e$ciency of the suggested deviatoric hybrid element approach for incompressible #uid analysis. 8. CONCLUSIONS The hybrid model based on the deviatoric variational formulation n(p , p, u ), which is universal GH G for compressible and incompressible mediums has been proposed. The greatest advantage of the model is that it is able to provide quite reliable fully incompressible solutions for both the stresses and the displacements. A dramatic result in the paper would be that the zero-volumetric strain constraints are enforced at element level by means of introducing the incompatible displacement function with a small disturbing parameter d. It follows that the pressure variables will not join the system discrete equations, and pure displacement formulations are "nally available. A practical strategy is suggested to suppress the zero energy displacement modes, ZEM(u), stemming from the bilinear term of p e dv in (4). This strategy, entails the deviatoric stress and T GH GH the deviatoric strain to take the same shape functions in element formulations. Another advantage of the present model is that it is able to unify the computations for both solids and #uids. The proposed hybrid element can be directly used in Stokes #ow calculations, and those troublesome numerical stability problems [16] will no longer appear. ACKNOWLEDGEMENTS The "nancial assistance from National Natural Science Foundation of China is gratefully acknowledged. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999) DEVIATORIC HYBRID MODEL AND MULTIVARIABLE ELIMINATION 745 REFERENCES 1. Herrmann LR. Elasticity equations for incompressible and nearly incompressible materials by variational theorem. AIAA Journal 1965; 3:1896}1900. 2. Taylor RL, Pister KS, Herrmann LR. On a variational theorem for incompressible and nearly incompressible elasticity. International Journal of Solids and Structures 1968; 4:875}883. 3. Key SW. A variational principle for incompressible and nearly incompressible anisotropic elasticity. International Journal of Solids and Structures 1969; 5:455}461. 4. Murakawa H, Atluri SN. Finite elasticity solutions using hybrid "nite elements based on a complementary energy principle, Part 2*Incompressible materials. Journal of Applied Mechanics 1979; 46:71}77. 5. Gadala MS. Numerical solutions of nonlinear problems of continua*II. Survey of incompressibility constraints and software aspects. Computers and Structures 1986; 22:841}855. 6. Malkus DS, Hughes TJR. Mixed "nite element methods*reduced and selective integration techniques: a uni"cation of concepts. Computer Methods in Applied Mechanics and Engineering 1978; 15:63}81. 7. Zienkiewicz OC, Taylor RL. ¹he Finite Element Method (4th edn), vol. 1. McGraw-Hill: UK, 1989. 8. Cook RD. Concepts and Applications of Finite Element Analysis (2nd edn). Wiley: New York, 1981. 9. Wu CC, Liu XY, Pian THH. Incompressible-incompatible deformation modes and plane strain plastic element. Computers and Structures 1991; 41:449}453. 10. Needleman A, Shih CF. Finite element method for plane strain deformations of incompressible solids. Computer Methods in Applied Mechanics and Engineering 1978; 15:223}240. 11. Zhang W, Li X. On the generalized displacement method in the "nite element analysis of incompressible materials. Journal of Dalan Institute of ¹echnology 1981; 20(4):1}10. in Chinese. 12. Li XK. Computerized generalized displacement method for "nite element analysis of incompressible material. Engineering and Computers 1984; 11:335}342. 13. Oden JT. Finite Element of Nonlinear Continua. McGraw-Hill: New York, 1972. 14. Thompson EC. Average and complete incompressibility in the "nite element method. International Journal for Numerical Methods in Engineering 1975; 9:925}932. 15. Wu CC. Dual zero energy modes in mixed/hybrid elements*de"nition, analysis and control. Computer Methods in Applied Mechanics and Engineering 1990; 81:39}56. 16. Girault V, Raviart PA. Finite Element Methods for Navier}Stokes Equations. Springer: Berlin, 1986. 17. Simo JC, Taylor RL, Pister KS. Variational and projection methods for the volume constraint in "nite deformation plasticity. Computer Methods in Applied Mechanics and Engineering 1985; 51:177}208. 18. Pian THH, Lee SW. Notes on "nite element for nearly incompressible materials. AIAA Journal 1976; 4:824}826. 19. Spilker RL. Improved hybrid stress axisymmetric elements including behavior for nearly incompressible materials. International Journal for Numerical Methods in Engineering, 1981; 17:483}501. 20. Tong P. An assumed stress hybrid "nite element method for an incompressible and nearly incompressible materials. International Journal of Solids and Structures 1969; 5:455}461. 21. Bratianu C, Atluri SN. A hybrid "nite element method for Stokes #ow: Part I*formulation and numerical studies. Computer Methods in Applied Mechanics and Engineering 1983; 36:23}37. 22. Hu HC. <ariational Principles in Elastic Mechanics and Applications. Science Press: Beijing, 1981. 23. Pian THH, Wu CC. A rational approach for choosing stress term for hybrid "nite element formulations. International Journal for Numerical Methods in Engineering 1988; 26:2331}2343. 24. Wu CC, Bu#er H. Multivariable "nite element: consistency and optimization. Science in China (A) 1990; (9):946}95 in Chinese, or 1991; 34(3):284}299 in English. 25. Strang G, Fix GJ. An Analysis of the Finite Element Method. Prentice-Hall: Englewood Cli!s, NJ, 1973. 26. Taylor RL, Beresford PJ, Wilson EL. A non-conforming element for stress analysis. International Journal for Numerical Methods in Engineering 1976; 10:1211}1219. 27. Wu CC, Cheung YK. The patch test condition in curvilinear coordinates*formulation and application. Science in China (A) 1992; (8):849}859 in Chinese, or 1993; 36(1):62}73 in English. 28. Wu CC, Pian THH. Incompatible Numerical Analysis and Hybrid Element Method. Science Press: Beijing, 1997. 29. Wu CC, Huang MG, Pian THH. Consistency condition and convergence criteria of incompatible elements, general formulation of incompatible functions and its application. Computers and Structures 1987; 27:639}644. 30. Pian THH, Chen DP. On the suppression of zero-energy deformation modes. International Journal for Numerical Methods in Engineering 1983; 19:1741}1752. 31. Cheung YK, Wu CC. A study on the stability of 3-"eld "nite element by theory of zero energy modes. International Journal of Solids and Structures 1992; 29:215}229. 32. Simo JC, Rifai MS. A class of mixed assumed strain methods and the method of incompatible modes. International Journal for Numerical Methods in Engineering 1990; 29:1595}1638. 33. Schlichting H. Boundary-¸ayer ¹heory (6th edn). McGraw-Hill: New York, 1968. Copyright 1999 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Engng. 46, 729}745 (1999)

1/--страниц