As we have discussed in some of our previous blogs about Element Selection in Abaqus, choosing the appropriate finite element type is critical for obtaining accurate and reliable predictions. However, with the vast array of available element formulations in Abaqus, this can sometimes be challenging, particularly for new users of the software. Because many novice analysts are unaware of the limitations of certain element types, they sometimes encounter disagreement between traditional hand calculations and finite element analysis predictions. But how can that happen? Both are mathematical representations of the same problem, right? While true, it is important to remember that FEA is based on discretized calculations – meaning that equations are computed at specific points and values in between must be interpolated according to the element’s shape function – and this can lead to potential errors (particularly if a given element is not suitable for capturing the physics being simulated). One such example, and the topic of today’s article, is the use of hexahedral elements in bending-dominated problems.
Hexahedral Elements
One of the primary advantages of hexahedral elements is the combination of accuracy (when elements are well-shaped) and runtime efficiency; this often makes hexahedral elements more economical than other solid element formulations. However, this does not come without a price – there are three primary concerns one must consider when using hexahedral elements to simulate bending-dominant modes of deformation:
Hourglassing
Shear Locking
Volumetric Locking
Before we continue, let us first recall that elements can be either linear or quadratic (first- or second-order) and can use either reduced or full integration, which is discussed further in this article about Abaqus Elements. This is an important distinction because hourglassing, shear locking, and volumetric locking are only problematic under certain conditions.
Hourglassing
Hourglassing can occur in first-order, reduced-integration elements (C3D8R). Because there is only one integration point, neither of the visualization lines in the schematic below have changed in length, and the angle between them is also unchanged by the imparted deformation. This can result in zero-strain deformation modes that produce highly spurious results.
This issue, however, can usually be mitigated by using second-order, reduced-integration elements (C3D20R).
Shear Locking
Shear locking can occur in first-order, fully-integrated elements (C3D8) that are subjected to bending. This occurs when artificial shear strain develops due to an inability of the element edges to bend. As a result of this nonphysical shear strain, these elements can be too stiff when used in bending-dominant problems. Fortunately, though, this issue can also be address by using second-order elements (C3D20).
Volumetric Locking
Volumetric locking can occur in first- and second-order, fully-integrated elements when material behavior is nearly incompressible (e.g., rubber). This causes unrealistic pressure stresses to develop at the integration points, which can make the element behave too stiffly. Fortunately, first-order, fully-integrated hexahedral elements use selectively-reduced-integration (reduced integration on the volumetric terms only), so these elements do not experience volumetric locking with nearly-incompressible materials. (As a side note, one should utilize the Hybrid element formulation when working with nearly-incompressible materials).
Example Problem
Now that we understand some of the problems that can arise when working with solid, hexahedral elements, let’s take a look at a simple example to help illustrate the importance of proper element selection.
Let us consider a bending problem in which a cantilever beam measuring 10 mm x 10 mm x 100 mm is fixed on one side and loaded with 1000 N on the opposing end. What should our value of maximum deflection be?
Hand Calculation
First, we must calculate the area moment of inertia for our beam with a rectangular cross-section. Since the cross-sectional width and height are the same, we know that Iy = Iz.
Iy = Iz = (Base * Height^3) / 12 = (10 * 10^3) / 12 = 833.33 mm^4
Now that we know the area moment of inertia, we can use the following equation to calculate the maximum deflection in the beam:
Deflection = (Load * Length^3) / (3 * Young’s Modulus * Moment of Inertia)
Deflection = (1000 * 100^3) / (3 * 200,000 * 833.33) = 2.000 mm
FEA Results
Next, let’s build a finite element representation of the problem above. In order to highlight some common errors often made when using hexahedral elements, sensitivity studies will be performed with first-order elements only:
Mesh Sizes: 0.50 mm, 1.0 mm, 2.0 mm, 5.0 mm (Aspect Ratio = 2.5)
Integration Strategy: Reduced-Integration & Fully-Integrated
When using reduced-integration elements (C3D8R), the deflection is overpredicted, sometimes severely (32% in the example below). This occurs because reduced-integration elements are prone to hourglassing, a numerical error which tends to make elements too flexible. This phenomenon becomes more pronounced as element size increases and each element is responsible for representing a larger amount of bending. Based on the results shown below, a margin of error of < 5% is achievable with five elements through the thickness (although the Abaqus manual suggests sufficient accuracy is achievable with as few as four).
When using fully-integrated elements (C3D8), on the other hand, the deflection is underpredicted by just as large a margin (34.5%, actually). This time, though, it is caused by shear locking, a phenomenon which adds artificial stiffness to the element (due to parasitic shear strains). Once again, we see the phenomenon becoming more pronounced with increasing element size: more bending per element means larger angular deviations between integration points and more apparent shear locking.
Lastly, let’s take a look at Abaqus’ incompatible mode element (C3D8I), which is intended to provide improved accuracy in bending-dominant problems. As promised, the element seems to live up to the hype (based on the results of this study): our coarsest mesh, which has only two elements through the thickness, was able to match the hand calculation within 5% (compared to a 32% and 34.5% deviation when using reduced-integration and fully-integrated elements, respectively).
Final Thoughts
If using hexahedral elements, use C3D8I! Yes, there is an added computational expense; but, if you’re concerned about accurately capturing bending behavior (without modeling several layers of elements), you can’t beat the improved accuracy of the incompatible mode elements.
More importantly, as we’ve hopefully shown above, the discretization strategy and element type chosen for a given finite element analysis can have a significant impact on the prediction, particularly when using first-order hexahedral elements in bending problems. More specifically, our example above shows a potential deflection range of ~ 1.3 to ~ 2.6 mm between reduced-integration elements and fully-integrated elements, respectively – this is literally double! What does that mean if two different analysts are working on the same project, but using different modeling strategies? If you are making design decisions based on your CAE predictions (and you probably should be), this makes a big difference! In fact, it is the difference between robust, reliable data that can increase your profitability and a potential liability.
When accuracy matters, contact the experts at Fidelis to help you realize the true value of simulation!