Belytschko T. - Introduction (779635), страница 92
Текст из файла (страница 92)
In most cases, the above is an equality, but it is possible, even with regularquadrature schemes and undistorted elements, to lose the equality, i.e. to have linearlydependent rows in the B matrix.The rank-sufficiency of QUAD4 will now be examined for various quadratureschemes. The element has 4 nodes with 2 degrees of freedom at each node, so dim (K e) =8. The number of rigid body modes is 3: translation in the x and y directions and rotationin the (x,y) plane. By Eq. (26a), the proper rank of K e = 5.The most widely used quadrature scheme is 2x2 Gauss quadrature. The number ofquadrature points n Q=4, the number of rows in each B(x α)=3, so the number of rows inoB =12 This exceeds the proper rank.
However, based on the linear completeness of thequadrilateral, it will be shown later that (see Section 1.4)ux x,y = α ox + α 1x x + α 2x y + α 3x h(1.2.33a)uy x,y = α oy + α 1y x + α 2y y + α 3y h(1.2.33b)h ≡ ξηThen27α 1x + α 3xoe = Du =α 2y + α 3yα 2x + α 1y + α 3x∂h∂x∂h(1.2.34)∂y∂h∂y+ α 3y∂h∂xExamination of the above shows that the strain-rate field contains 5 linearly independentfunctions: α 1x , α 2y , α 3x∂h/∂x, α 3y∂h/∂y, and α 2x + α 1y .
Note that the two constants inthe shear strain field must be considered as a single independent field and the functionα 3x∂h/∂y + α 3y∂h/∂x cannot be considered linearly independent because it is a combinationof the two functions that have already been included in the list. Thus dim Du = 5 and sinceorows in B = 12, it follows from (32) thatorank(B )=5oIt may be concluded that for any quadrature scheme the rank of B for QUAD4 cannotexceed 5.The rank of the element stiffness of QUAD4 for one-point quadrature can beoascertained similarly. In one-point quadrature, B consists of B evaluated at a single point,so its rank is 3.
Therefore rank K e is 3 by Eq. (30), and Eq. (26b) indicates that theelement has a rank deficiency of 2. This rank deficiency can cause serious difficultiesunless it is corrected. Such corrective procedures are described later.1 . 2 . 4 . Nodal Forces and B-Matrix for One-Point Quadrature Element. Priorto describing procedures for correcting the rank deficiency of QUAD4 with one pointquadrature, it is worthwhile to develop the one-point quadrature formulas in detail. Theseformulas will then provide the framework for the development of the rank correctionprocedures, which in QUAD4 are often called hourglass control.The internal nodal forces are given by (31b). When one-point quadrature is used, thequadrature point is selected to be the origin of the coordinate system in the reference plane:ξ = η = 0. Evaluating the Jacobian at this point yields (see Eq.
(9)) givesJ(0 ) = 1 x31 y42 + x24 y31 = A(1.2.35)84where A is the area of the element. The expression for the internal nodal forces nowbecomesf int = 4B t (0)s (0)J(0) = AB t (0)s (0)(1.2.36)Evaluating the B matrix from Eq. (10) at ξ = η = 0 is a simple algebraic processwhich gives28b tx00btybtybtxB(0) =(1.2.37)wherebtx ≡ bt1 = 1 y24 , y31 , y42 , y132A(1.2.38)ttby ≡ b2 = 1 x42 , x13 , x24 , x312ASince one-point quadrature is used, the nodal forces are simply the product of the area andthe integrand evaluated at ξ = η = 0, which using (37) and (39) gives b 0 by σ x f int = A x(1.2.39) σ y = AB t (0)s (0)0bbyx σ xyThe element stiffness matrix for the underintegrated element can be obtained by usingthe stress-strain laws = Ce(1.2.40)in conjunction with Eqs. (37), (39) and e=Bd.
This giveswheref int = K ed(1.2.41)K e = AB t (0)CB(0)(1.2.42)The element stiffness matrix could also be obtained from (28) by using one-pointquadrature and the values of B and C at the quadrature point.1.2.5 Spurious Singular Modes (Hourglass) The presence and shape of the spurioussingular modes of the one-point quadrature QUAD4 element will now be demonstrated.Any nodal displacement dH that is not a rigid body motion but results in no straining of theelement is a spurious singular mode.
From (43) it can be seen that such nodaldisplacements will not generate any nodal forces, i.e. they will not be resisted by theelement, since in the absence of strains, the stresses will also vanish, so f int = 0.Consider the nodal displacementsdHx = h0dHy = 0h(1.2.43)ht = +1, -1, +1, -1It can easily be verified thatbtx h = 0bty h = 0Therefore, it follows from Eqs. (37), (43) and (44) that(1.2.44)29B(0)dHx = 0(1.2.45a)B(0)dHy = 0(1.2.45b)Fig.
3 shows rectangular elements with the spurious singular modes of deformation dHxand dHy, and also both modes simultaneously. In the rectangle, it can be seen that thehourglass modes are associated with the bilinear term in the displacement field. Thedeformed configuration of the mesh in the spurious singular modes is shown in Fig. 4. Avertical pair of elements in this x-mode looks like an hourglass, an ancient device formeasuring time by the flow of sand from the top element to the bottom. For this reason,this spurious singular mode is often called hourglassing or the hourglass mode. Becauseeach element is in an hourglass mode, the entire mesh can deform as shown without anyresisting forces from the element.d = dHxd = dHyd= hhFigure 3.
Hourglass modes of deformationThe problem of hourglassing first appeared in finite difference hydrodynamicsprograms in which the derivatives were evaluated by transforming them to contour integralsby means of the divergence theorem; see for example Wilkins and Blum (1975). Thisprocedure tacitly assumed that the derivatives are constant in the domain enclosed by eachcontour. This assumption is equivalent to the constant strain (and stress) assumptionwhich is associated with one-point quadrature. The equivalence of these contour-integralfinite difference methods was demonstrated by Belytschko et al.
(1975); also seeBelytschko (1983). Many ad hoc procedures for hourglass control were developed byfinite difference workers. The procedures focused on controlling the relative rotations ofelement sides; no consideration was given to maintaining consistency.Figure 4. Mesh in hourglass mode of deformation30This phenomenon occurs in many other settings, so a variety of names have evolved.For example, they occur frequently in mixed or hybrid elements, where they are calledzero-energy modes or spurious zero-energy modes. Hourglass modes are zero-energymodes, since they don't result in any strain at the points in the element which are sampled.Therefore they do no work and (dHx) t fint = ( dHy) t fint = 0In structural analysis, spurious singular modes arise when there is insufficientredundancy, i.e.
the number of structural members is insufficient to preclude rigid bodymotion of part of the structure. Such modes often occur in three dimensional trussstructures. In structural analysis, they are called kinematic modes, and because of the closerelationship between the structural analysis and finite element communities, this name hasalso been applied to spurious singular modes. Other names which have been applied to thisphenomenon are: keystoning (Key et al. (1978)), chickenwiring, and mesh instability.For finite element discretizations of partial differential equations, spurious singularmodes appears to be the most accurate term for this phenomenon, so we shall use thatname. For example, the terms kinematic modes or zero-energy modes are not appropriatefor the Laplace equation.
In elements where the spurious singular mode has a distinctiveappearance, such as the hourglass pattern in QUAD4, we shall also use that name. Thecondition which leads to spurious singular modes is rank deficiency of the element stiffnessmatrix.When rank deficient elements are assembled, the system stiffness will often besingular or nearly singular.
Therefore, in matrix methods, the presence of spurioussingular modes can be detected by the presence of zero or very small pivots in the totalstiffness. If the pivots are zero, the stiffness will be singular and not invertible. If thepivots are very small, the total stiffness is near-singular, and the displacement solutionswill be oscillatory in space, i.e., they will exhibit the hourglass mode.Because a system stiffness matrix is never assembled in explicit methods, thesingularity cannot be readily detected. In iterative solvers, the presence of spurioussingular modes will often lead to divergence of the solution. With explicit integrators,singular modes are not readily detectable without plots of the deformed configuration.
Thisis also true for matrix dynamic methods, since the mass matrix then renders the systemmatrix nonsingular even when the stiffness matrix is singular.The evolution of an hourglass mode in a transient problem is shown in Fig. 5. In thisproblem, the beam was supported at a single node to facilitate the appearance of thehourglass mode. If all nodes at the left-hand end of the beam were fixed to simulate aclamped support condition, the hourglass mode would not appear.
Although rank-deficientelements may sometimes appear to work, they should not be used without an appropriatecorrection.31Figure 5. Four views of a simply supported beam showing the evolution of the hourglassmode (due to symmetry, only half the beam was modeled)1 .
3 Perturbation Hourglass StabilizationThe simplest way to control spurious singular modes without impairing convergenceis to augment the rank of the element stiffness without disturbing the linear completeness(consistency) of the isoparametric element. One approach to this task is to augment theB(0) matrix of the one-point quadrature element by two rows which are linearlyindependent of the other three. These additional two rows consist of a g vector that will bederived subsequently.