One-Point Quadrature Shell Element
In this section, a one-point quadrature shell element formulation will be developed from the general formulation described in the previous section. It is based on the Physical Stabilization method which explicitly computes the stabilization terms in making some simplifications.
The following formulations will be written in the local coordinate system [ ] (the circumflex in the co-rotational system notation has been omitted for convenience).
Kinematic Approximation
The projection consists in eliminating the nodal drilling rotations in order to reinforce Mindlin's kinematic condition at the nodes. It has been pointed out 1 that without this projection, the element is too flexible and cannot pass the Twisted Beam test.
This projection has a drawback of changing the invariance property to rigid body motion, that is, a warped element being invariant to rigid body rotation will now strains under to rigid body rotation if the drill projection is applied. To overcome this problem, a full projection proposed by Belytschko-Leviathan 2 which free either drilling rotation or rigid motion should be used. This full projection is only used for QEPH element.
In-plane Strain-rate Construction
Constant part
with:
is the area of element.
The derivation of the shape functions is given by:
The advantage of this shape function form is that a linear field expressed with Cartesian coordinates and a bilinear field expressed with Natural coordinates is decomposed so that the constant part is directly formulated with the Cartesian coordinates, and the non-constant part is to be approached separately.
with:
with
The parameter is a measure of the warping of the element.
The first term is neglectable. You have verified that the order of is times the second term of with which vanishes when the element is rectangular. Thus, this term is not used in the program.
The constant part of the in-plane rate-of-deformation formulation without the term is consistent with the result of Belytschko's family shell element 4, 6though this part has been obtained in a very different manner. Letellier has given the same result in his thesis 5, and studies were also made of the quadratic terms with respect to .
Non-constant Part
The main simplification for the non-constant part formulation, in order to overcome the difficulties described above, is:
The element is considered to be flat.
and:
This will lead to 'membrane locking' (the membrane strain will not vanish under a constant bending loading). According to the general formulation, the coupling is presented in the bending terms not in the membrane terms, yet the normal translation components in do not vanish for a warped element due to the tangent vectors which differ from (0,0).
Out-plane Strain-rate Construction
The straightforward form of [] is obtained using one additional simplification:
which is true for a parallelogram element. Although this simplification is not necessary, it is justified by the fact that the transverse shear terms serve mainly as a penalty function.
Explicit Integration of the Nodal Internal Force Vector
Elastic Case
The elementary nodal force vector is computed by:
- Constant in-plane fields along with the non-constant ones
- Membrane and the bending
- In-plane fields and the out-plane fields
- Decomposed non-constant out-plane fields
Resulting in:
With the constant part being computed with one-point quadrature, and
These can be evaluated explicitly.
The rate-of-deformations will be written explicitly.
The rate form of the constitutive relation is written as (stress plane for in-plane terms):
With the assumption: the spin is constant within the element, the objectivity principle will be satisfied. The incremental computation is performed with the hourglass generalized rate-of-deformation :
Noting that if is considered as constant over a time step, it is equivalent to the incremental stress computation.
Physical Nonlinear Case
Now consider an elastoplastic problem.
The elementary nodal internal force vector is now computed by:
The constitutive relation is written by either a tangent form: , or a projection form:
Where, is the history-dependent tangent tensor; is the trial stress, and the function consists of projecting the trial stress on the updated yield surface.
The decomposed form of the last equation is written:
The constant part is computed by integration at each integration point through the thickness.
- Keep the same orthogonalities as in the elastic case, and
- Use the unidimensional tangent modulus to evaluate the non-constant rate-stress, that
is,
(19) Where,- Young's modulus
- Matrix of elastic moduli
Thus, the elastic case easily extends to the nonlinear case.
Where, is obtained by the constant stress incremental computation along the thickness and in the elastic zone, and is the average value of and .
For the QPH shell, = 1
The key orthogonalities has been maintained without any significant deterioration in performance, although the first two orthogonalities might have been slightly violated. In fact, it is simply due to these orthogonalities that a one-point quadrature element dramatically reduces the computation cost; otherwise you return to the full integration scheme.
Most of the physical stabilization elements have incorporated the following assumption.
The material response is constant within the element.
- Take the elastic matrix directly:
which means that the plastic rate-of-deformation is constant within the element; or
- Take the tangent matrix [] ( constant).
Since the components of [] are generally functions of the updated stress (precisely, the stress deviator for an elastoplastic problem with an associative flow rule, which means that the stress is constant within the element.
Neither of these alternatives is theoretically perfect. Note that the option results in a contradiction with the stress computation (which yields different results for the constant part and the non-constant part); it is more expensive and the tangent form is not generally used for constant stress computations within an explicit program. Hence, the approximation based on the above assumption is not necessary.
The choice of the moduli for the nonlinear case has not been studied for Belytschko-Leviathan's element 6, and it has been shown that this choice has little effect on the result of the "Cylindrical panel test". In 7, the elastic tangent matrix has been used for the evaluation of the stabilizing forces.
QPH,QPPS have shown often stiffer behaviors and sometimes have certain numerical problems in crash simulations.