2.3.5 Performance of continuum and shell elements for linear analysis of bending problems

Products: ABAQUS/Standard  ABAQUS/Explicit  

It is well known that fully integrated linear isoparametric continuum elements, in both two and three dimensions, are too stiff in modeling the simple flexural deformation of a beam. Similarly, fully integrated standard displacement formulations for 4-node shell elements are too stiff in modeling bending about an axis perpendicular to the plane of the shell; i.e., in-plane bending. Full integration refers to the Gauss integration order required for exact integration of the polynomial of the order being integrated when the element is rectangular. Although full-order integration elements can represent rigid-body and constant-strain displacement fields exactly, they tend to “lock” in bending problems because a disproportionately large shear-related strain energy arises, which greatly increases the flexural rigidity of the model.

Problem description

Geometry and models

Three examples are considered here to illustrate the behavior of the various elements in modeling bending behavior.

The same number of first- and second-order elements is used for each mesh layout. Obviously, for an equal number of elements, the mesh for the second-order elements will contain more nodes and, thus, have a larger number of degrees of freedom. However, the objective is to show that, even with fewer degrees of freedom, the incompatible mode elements and—to a certain extent—the reduced-integration, linear elements give comparably accurate results with respect to the second-order elements. Hence, the meshes included here are coarse and should not be taken as good modeling practice.

Cantilever beam with end shear load

The first example is the linear static analysis of a cantilever with an end shear load. The geometry is shown in Figure 2.3.5–1. The Young's modulus is 689.5 GPa (10 × 107 lb/in2), and Poisson's ratio is 0. In ABAQUS/Standard and ABAQUS/Explicit the following four meshes are used: 1, 2, 4, or 8 elements through the depth of the beam, combined with 4 or 16 elements along the beam. The 8 × 16 mesh is included to examine the “converged” solution. The resulting end deflections are normalized with respect to the Euler-Bernoulli beam theory prediction of 2.74 mm (0.108 in) and tabulated in Table 2.3.5–1. Although this example is a beam problem, the thickness is chosen to be small so that the problem can be reasonably modeled with both continuum and shell elements. Due to the small thickness, warning messages may be issued for poor aspect ratios for analyses involving coarser meshes of three-dimensional continuum elements. The analytical solution depends only on the thickness through the combination Young's modulus times thickness (). Hence, any thickness can be used so long as the product remains fixed; the solution remains the same with the exception of the tetrahedral elements, where the thickness influences the tetrahedral quality measure and, hence, the displacement solution.

In the ABAQUS/Explicit analyses the load is applied to the structure using a SMOOTH STEP amplitude curve to minimize the dynamic effects.

Skew sensitivity analysis of cantilever beam

The second example examines the skew sensitivity of the elements with respect to two shapes: parallelogram and trapezoidal. The cantilever beam and loading are the same as in the first example. A 1 × 8 mesh is used for both shapes when the elements are quadrilateral or hexahedral, as shown in Figure 2.3.5–2 and Figure 2.3.5–3. When the problem is modeled using triangular or tetrahedral elements, the basic parallelogram or trapezoidal shape meshed with one quadrilateral or hexahedral element is filled with either two triangular elements in the two-dimensional case or five tetrahedral elements in the three-dimensional case. Therefore, the triangular and tetrahedral elements do not have a true parallelogram or trapezoidal shape. The following skew angles are tested for the ABAQUS/Standard and ABAQUS/Explicit runs: 0°, 15°, 30°, and 45°. The resulting end deflections are again normalized with respect to the Euler-Bernoulli beam theory prediction of 2.74 mm (0.108 in) and tabulated in Table 2.3.5–2 and Table 2.3.5–3. The use of a single layer of elements precludes the testing of the reduced-integration linear isoparametric elements (except when using enhanced hourglass control) since a minimum of four layers is required for acceptable results.

Buckling analysis of ring under external pressure

The final example is the buckling analysis of a ring in a plane under external pressure. The geometry and material are the same as those used in Buckling of a ring in a plane under external pressure, Section 1.2.2, and are shown in Figure 2.3.5–4. The ring buckling problem requires a rather fine mesh in the circumferential direction, presumably to model the strain gradients accurately in the buckling mode. Using symmetry boundary conditions, only a 45° sector of the ring is modeled, which is enough to reproduce the primary buckling mode. Table 2.3.5–5 gives the solutions obtained with various models for this case. The exact solution, based on Euler-Bernoulli beam theory, is a critical buckling pressure of 51.71 KPa (7.5 lb/in2).

Results and discussion

Additional modeling considerations

These analyses illustrate the use of continuum elements to simulate bending behavior of thin structures. Similarly, these analyses illustrate the use of shell element S4 for bending in the plane of the element. In general, the use of beam or shell elements (with out-of-plane bending) is recommended for thin structures; the difficulties discussed here would be encountered only if the analyst could not use the more appropriate bending elements for some reason. The results presented here do not reflect the usefulness or importance of continuum elements in other types of problems.

We have seen that the incompatible mode elements perform almost as well as the second-order elements in many situations if the elements have an approximately rectangular shape. The performance is considerably poorer if the elements have a parallelogram shape and quickly becomes unacceptable with trapezoidal element shapes. In addition, these elements offer no advantages when a side is degenerated or collapsed into a node. The degenerated elements can only represent a constant strain field, and the incompatible modes cannot improve on such fields.

Due to the internal degrees of freedom (4 for CPS4I; 5 for CPE4I, CAX4I, and CPEG4I; and 13 for C3D8I) the incompatible mode elements are somewhat more expensive than regular displacement elements but are more economical than the second-order elements. The additional degrees of freedom do not increase the wavefront size substantially, since they can be eliminated immediately. In addition, it is not necessary to use selectively reduced integration, which partially offsets the cost of the additional degrees of freedom.

The reduced-integration, linear elements also give satisfactory solutions for the set of problems attempted here when a minimum of four layers is used. However, there are cases where the elements may yield a nearly singular stiffness and, hence, physically unreasonable solutions. This is especially true when large-strain problems are analyzed. Thus, as with all modeling decisions, analysts must develop their discretization with careful testing of the effectiveness of the elements for a particular application.

Further examples

The use of the elements discussed above is illustrated further in the following example problems:

Input files

References

Tables

Table 2.3.5–1 Normalized tip deflection of a cantilever beam (). Bernoulli-Euler theory prediction:  2.743 mm (0.108 in).

Element Type Mesh Size (Depth × Length)
1 × 4 2 × 44 × 48 × 16
CPS30.0120.0120.0120.159
CPS4I0.9850.9850.985 1.000
S40.8990.9430.9370.966
S4R**0.9850.9850.9851.000
S4R***0.9850.9850.9851.000
CPS40.0340.0340.0340.363
CPS4R*1.1510.9441.008
CPS4R**0.9850.9850.9851.000
CPS4R***0.9850.9850.9851.000
C3D40.0010.0010.0020.065
C3D8I0.9850.9850.9851.000
C3D80.0350.0340.0340.364
C3D8R*1.3061.0501.016
C3D8R**0.9840.9850.9851.000
C3D8R***0.9850.9850.9851.000
CPS60.9860.9860.9861.000
CPS6M0.9400.9460.9470.999
CPS80.9870.9870.9871.000
CPS8R1.0011.0011.0011.001
C3D100.9850.9850.9851.000
C3D10M1.0210.9850.9691.000
C3D200.9870.9870.9881.000
C3D20R1.0011.0011.0011.001
* yields singular stiffness matrix
** ABAQUS/Explicit with enhanced hourglass control
*** ABAQUS/Standard with enhanced hourglass control

Table 2.3.5–2 Normalized tip deflection of a cantilever beam with parallelogram-shaped elements (). Bernoulli-Euler theory prediction:  2.743 mm (0.108 in).

Element Type Skew Angle
15°30°45°
CPS3 0.042 0.032 0.022 0.017
CPS4I 0.996 0.898 0.791 0.742
S4 0.903 0.833 0.470 0.226
S4R**0.9960.8980.7910.742
S4R***0.9960.8980.7910.742
CPS40.1250.1100.0790.049
CPS4R** 0.996 0.898 0.791 0.742
CPS4R***0.9960.8980.7910.742
C3D4 0.001 0.001 0.002 0.002
C3D8I 0.997 0.898 0.791 0.742
C3D80.1320.1210.0930.061
C3D8R** 0.996 0.897 0.790 0.742
C3D8R***0.9960.8970.7910.742
CPS6 0.997 0.982 0.931 0.821
CPS6M 0.991 0.985 0.965 0.926
CPS8 0.998 0.998 0.996 0.988
CPS8R 1.001 1.001 1.000 0.997
C3D10 0.997 0.897 0.711 0.484
C3D10M 1.040 1.004 0.920 0.814
C3D20 0.998 0.998 0.983 0.961
C3D20R 1.001 0.988 0.980 0.896
SC8R***0.9960.8980.7910.742
** ABAQUS/Explicit with enhanced hourglass control
*** ABAQUS/Standard with enhanced hourglass control

Table 2.3.5–3 Normalized tip deflection of a cantilever beam with trapezoidal-shaped elements (). Bernoulli-Euler theory prediction:  2.743 mm (0.108 in).

Element Type Skew Angle
15°30°45°
CPS3 0.042 0.041 0.034 0.025
CPS4I 0.997 0.411 0.140 0.067
S4 0.903 0.469 0.169 0.102
S4R**0.9970.4110.1400.067
S4R***0.9960.4110.1400.067
CPS40.1250.1020.0600.035
CPS4R** 0.997 0.411 0.140 0.067
CPS4R***0.9960.4110.1400.067
C3D4 0.001 0.002 0.006 0.010
C3D8I 0.997 0.411 0.140 0.067
C3D80.1320.1080.0630.037
C3D8R** 0.996 0.410 0.140 0.067
C3D8R***0.9960.4110.1400.067
CPS6 0.997 0.997 0.994 0.986
CPS6M 0.991 0.990 0.990 0.990
CPS8 0.998 0.959 0.985 0.915
CPS8R 1.001 0.971 0.996 0.981
C3D10 0.997 0.995 0.986 0.963
C3D10M 1.040 1.038 1.035 1.024
C3D20 0.998 0.959 0.985 0.914
C3D20R 1.001 0.956 0.984 0.974
SC8R***0.9960.4110.1400.067
** ABAQUS/Explicit with enhanced hourglass control
*** ABAQUS/Standard with enhanced hourglass control

Table 2.3.5–4 Normalized tip deflection of a cantilever beam using continuum shell elements with skewing in the thickness direction of the element (). Bernoulli-Euler theory prediction:  2.743 mm (0.108 in).

Mesh Skew Angle
15°30°45°
Parallelogram-shaped elements1.00 1.001.001.00
Trapezoidal-shaped elements1.001.001.001.00

Table 2.3.5–5 Normalized eigenvalue pressure estimates ( Bernoulli-Euler theory prediction:  51.71 KPa (7.5 lb/in2).

Element Type Mesh Size (Depth × Length)
4 × 104 × 20
CPS393.8324.27
CPS4I1.0051.001
CPS432.008.72
CPS4R1.0590.968
CPS4R*1.0061.001
C3D4 64.27 16.82
C3D8I 1.005 1.001
C3D8 31.98 8.70
C3D8R 1.043 0.966
C3D8R*1.0061.001
CPS6 1.010 1.001
CPS6M 1.002 1.000
CPS8 1.027 1.002
CPS8R 1.000 1.000
C3D10 1.013 1.001
C3D10M 1.014 0.998
C3D20 1.027 1.002
C3D20R 1.000 1.000
* ABAQUS/Standard with enhanced hourglass control


Figures

Figure 2.3.5–1 Typical mesh and loading.

Figure 2.3.5–2 Typical parallelogram-shaped elements.

Figure 2.3.5–3 Typical trapezoidal-shaped elements.

Figure 2.3.5–4 Ring buckling problem.

Figure 2.3.5–5 Parallelogram-shaped continuum elements: tip displacement versus skew angle.

Figure 2.3.5–6 Trapezoidal-shaped continuum elements: tip displacement versus skew angle.