Pure beam elements may be 2D or 3D elements, depending on whether or not they are employed for 2D or 3D analysis problems. Beam elements are defined by two nodes. In a 2D analysis all beam elements in the model lie in a single plane, assumed here to be the x-y plane. No bending or displacement in the out of plane direction is allowed. Each node of the beam element admits three degrees of freedom: translation in the x direction, translation in the y direction, and rotation about the z axis. Figure 2 shows an example of a 2D beam element. Note the rotational degrees of freedom q1 and q2, in addition to the translational degrees of freedom u1, v1, u2, and v2. For the 2D beam element, three geometric parameters of the element must be specified: cross sectional area A, the area moment of inertia of the element’s cross section about the z axis Izz, and the section height h.

In contrast to the uniaxial element, beam elements admit no stretching or compression stiffness. Thus, a force applied along the axial direction of a beam element is not resisted by the element. Only the component of force (or translation) perpendicular to the element’s axial direction results in elastic deformation and hence strain and stress.

Three-dimensional pure beam elements also are defined by two nodes but have six degrees of freedom per node corresponding to translation in the three coordinate directions and rotations about each of the coordinate axis. The element has stiffness due to bending about the two perpendicular axes that lie in the plane of the element’s cross section, as well as stiffness due to torsion or twisting of the element about its axial axis. As a minimum, the geometric information necessary to define the element include its cross sectional area A, two area moments of inertia corresponding to bending about each of the axis that lie in the plane of the element’s cross section, the moment of inertia about the element’s axial axis, and the height and width of the element’s cross section.  

3) Combined uniaxial and pure bending line elements (frame or standard beam elements)

This element is a direct superposition of the uniaxial tension-compression element and the pure beam bending element. Consequently, the frame element has both uniaxial and bending stiffnesses. The 2D and 3D forms of this element have identical degrees of freedom as the 2D and 3D pure beam elements. Thus, there is very little additional programming complexity or computational expense to add uniaxial stiffness capabilities to a pure beam element to create this more general element. Consequently, for many commercial finite element codes this beam element may be the default beam element offered, and a pure bending-only beam element may not even be offered.

Surface or Area Elements

There are two types of surface or area elements: two-dimensional planar elements and 3D shell or plate elements.

1. Planar Elements

It is cost efficient to analyze some physical systems using two-dimensional planar elements. Since the real world is inherently three-dimensional, such a two-dimensional analysis represents a simplification or idealization of the analysis problem. The extent to which the system can be simplified to a two-dimensional problem depends on factors including the objectives of the analysis, the accuracy and resolution requirements, and how much the system deviates from behaving in a purely two-dimensional manner. The following conditions are required for 2D analysis.

• the geometry of the system must be uniform or nearly uniform in the out-of-plane direction. In the engineering CAD world, solid shapes are typically created by first creating the planar 2D geometry and then extruding or sweeping this geometry in the out-of-plane direction to create the 3D volume.
• loads must lie in or nearly lie in the 2D plane of the analysis.
• boundary conditions must be uniform or nearly uniform in the out-of-plane direction.
• material behavior (i.e. material properties) must be uniform or nearly uniform in the out-of-plane direction.

Planar elements permit only two degrees of freedom per node: displacement in the x direction and translation or displacement in the y direction.

In solid mechanics there are three different types of planar analysis problems: 2D plane strain , 2D plane stress, and 2D axisymmetric. In 2D plane strain analysis, there are no out-of-plane displacements (and hence no out-of-plane strains), and there are no shear strains or shear stresses acting on planes perpendicular to the 2D plane that is being analyzed. The lack of out-of-plane displacements implies that either the object's out-of-plane dimension is much, much larger than its in-plane dimensions (i.e. a very thick object), or the out-of-plane ends of the object are fixed against out-of- plane motion.

In 2D plane stress analysis, there is no out-of-plane stress, or shear strains or shear stresses acting on planes perpendicular to the 2D plane that is being analyzed. Physically, this implies the object's out-of-plane dimension is much, much smaller than its in-plane dimensions, i.e. a very thin object.

In 2D axisymmetric analysis the 2D region being analyzed represents the cross section of an object which, when revolved around an axis of revolution, results in the 3D geometry of the system. Thus, axisymmetric problems are applicable for bodies of revolution, such as cylinders, tubes, spheres, etc. In this case the out-of-plane direction is the circumferential direction and in-plane directions correspond to the radial and axial directions. There are no displacements in the circumferential direction, yet circumferential strains and stresses do in the form of hoop strains and stresses. The hoop strain at any point is directly proportional to the radial displacement of the point. Thus, the 2D finite element displacement solution is used to compute hoop strains (and then hoop stresses by using material property relationships).

2. Shell or Plate Elements

Shell or plate elements are used to model thin 3D structures, typically those that are acted upon by bending-type loads. The element uses a different stiffness formulation than a standard solid element that permits it to model bending deformations much more accurately and with fewer degrees of freedom than solid elements. Further, since the element is a surface element, only surfaces are meshed to create the finite element model. A thickness value (or even a thickness value for each node of the element) is associated with each shell element corresponding to the thickness of the structure in the direction normal to the element (i.e. in the thin direction).

Ideally, the surface that is meshed is what is called the “neutral surface or plane.” For most problems this surface is the mid-plane located halfway between the two surfaces defining the volume. Thus, for a skull model the mid-plane is the surface halfway between the outer and inner surfaces of the skull. Some CAD and FEA tools have automatic mid-plane extraction algorithms to create this surface from a volumetric model. However, we have yet to find a mid-plane extraction algorithm robust enough to handle the geometric complexity of thin structures, such as the bat skull. This is indeed unfortunate, as a mid-plane surface (once constructed) would be easy to automatically mesh with quadrilateral or triangular-shaped shell elements. The resulting finite element model would have significantly fewer total degrees of freedom than a volumetric tetrahedral element model, and yet would be more accurate for virtually any type of transverse loading scenario. Since many, if not most, biological systems exhibit a thinness property, it is our hope that next generation FEA tools will exploit this fact by including more robust mid-plane extraction algorithms.

Solid Elements

4-noded linear tetrahedral

The simplest solid element is the four-noded linear tetrahedral element, as shown in the above figure. The element is defined by four nodes with each node having three degrees of freedom, the displacements in the x, y, and z direction at the node. The element is a linear element because the displacement fields that the element admits are linear functions in terms of the coordinates x, y, and z. Thus, if we let u denote the displacement field in the x direction, v denote the displacement field in the y direction, and w denote the displacement field in the z direction, then u = f(x,y,z), v = g(x,y,z) and w = h(x,y,z) . The functions f, g, and h for each element are complete linear polynomial functions in terms of spatial coordinates x, y, and z. Note that it takes four terms to define a complete linear polynomial in a 3D space, i.e.

(1)

and the element has precisely four nodes. This means no additional ‘higher order’ terms are needed beyond the linear terms shown in the equation above.

In the finite element formulation, the displacement field through out the element is obtained by interpolation of nodal displacements using linear interpolating polynomial functions, which are also called shape functions in the FEA literature:

(2)

Knowing the nodal coordinates, the linear interpolating functions N i(x,y,z) are easy to construct, and they have the following properties:

(3)

where the subscript indicates the appropriate local node number of the element (i.e. (x 1,y 1,z 1) are the nodal coordinates for node 1 of the element). Since strains are given by the partial derivatives of the displacement fields and the displacement fields are linear functions with no additional high order terms, this means that the strain fields in the element are constant. Thus, the four-noded tetrahedral element is a constant strain element; each of the six components of the 3D strain tensor (the three normal strains and 3 shear strains ) is a single value over the element volume. Since thin objects under bending-type loads exhibit strains and stresses that vary linearly through the thickness direction, it takes at least a couple of tetrahedral elements meshed through the thickness of the structure to roughly approximate the type of stress distribution that occurs in thin objects subjected to bending type loads.

Thus, the four-noded linear tetrahedral element performs poorly for bending type load problems. For force loading applications, the FE model will be overly stiff and under-predict both deflection and bending stresses. However, with sufficient number of elements, and therefore nodal degrees of freedom, the theoretically exact solution admitted by the Theory of Elasticity can be approached to any desired level of accuracy.

Despite these drawbacks, four-noded tetrahedral element models are quite common in models of biological systems. The main reason for this is that automatic tetrahedral mesh generators are fairly efficient and robust. Almost any complex, but properly defined, 3D volume can be meshed with tetrahedral elements if the element size is made small enough.

2) 10-noded linear tetrahedral

A tetrahedron has six edges. If a node is added to the middle of each edge of a four-noded tetrahedral element, a ten-noded tetrahedral element is created. With 3 degrees of freedom per node, the element now has a total of 30 degrees of freedom. The six additional nodes enable 6 additional terms to be appended to the four terms that define a linear polynomial functions as shown in Equation (1). Serendipitously, a complete quadratic polynomial in terms of spatial coordinates x, y, and z with no additional cubic or higher-order terms has exactly 10 terms. Thus, the ten-noded tetrahedral element employs complete quadratic interpolating functions to interpolate the displacement field from nodal values. This means that the strain tensor for the element is fully linear (i.e. each of the various strain components can vary as a linear polynomial over the element volume), and the element performs much better under bending types of loads. It also turns out that the element can now take on a quadratic shape. Thus, for example, the edges of the element may be curved to fit a portion of a circular arc.

Virtually all FE tools will automatically generate a mesh of 10-noded quadratic tetrahedral elements from a mesh of 4-noded linear tetrahedral elements by inserting nodes at the mid-point on every element edge. The model is then resolved, and the results compared to the previous results to determine if the solution has changed significantly. Increasing the order of the displacement field admitted by the element by adding additional nodes without changing the number of elements is called p- refinement, as in the literature the variable p is often used to denote the order of the polynomial function describing the element displacement field. Increasing the number of elements without the changing the element order by decreasing the element size and remeshing is called h-refinement, as the variable h is often used in the literature to designate the element size. Research has shown that a combination of h- and p- refinement, called h-p refinement, will result in the fastest convergence rate, especially if the refinement is done locally where most needed, such as in a region of high strains.

While the 10-noded tetrahedral element performs pretty well, the downside is the additional degrees of freedom required. Further, it takes a significant number of tetrahedral elements to mesh a volume. Consider a brick-like object. A minimum of 6 tetrahedral elements are required to mesh it.

3) 8-noded linear hexahedral or brick element

Figure 5 illustrates an 8-noded linear hexahedral or brick element. Because the element has eight nodes, 4 additional higher order terms are included in the element’s displacement field, so that for example the displacement field in the x-direction is a polynomial function of the form

(4)

Thus, unlike the linear tetrahedral element, the linear hexahedral is not a constant strain element, as the addition of these higher order terms enables the strain field in the element to vary linearly as a function of spatial coordinates. The result is the element performs much better, especially under bending types of loads, than the linear tetrahedral element. Despite this fact, linear hexahedral elements still is not able to undertake the proper bending shape. When a structure bends, it exhibits curvature as it deforms. Yet the edges of the linear hexahedral elements remain straight lines, even when pure bending type loads are applied to the element. Essentially, the element absorbs some of the work exerted on it by pure bending loads in shear strain energy instead of all of the work being absorbed in pure bending strain energy. This is called the parasitic shear effect. To improve the element performance in bending applications, many commercial FEA codes modify or augment the element formulation to include nine “nodeless” internal degrees of freedoms by adding additional polynomial functions (i.e. shape functions) to the element formulation. The result is that the element is able to admit a displacement field with curvature. These internal “nodeless” degrees of freedom are hidden from the user, but the documentation for the element may provide the option for the user to include these “extra shape functions,” so that the element can perform better under bending-type loads. Of course, even without these extra shape functions, a mesh of hexahedral elements will converge to the theoretically exact solution admitted by the Theory of Elasticity with sufficient mesh refinement. Further, the unmodified linear hexahedral element significantly out performs the linear tetrahedral element in bending applications.

Given this fact and the fact that it takes significantly fewer hexahedral elements to mesh a volume than tetrahedral elements, the “Holy Grail” of automatic FEA mesh generation research and development has been a robust hexahedral automatic mesh generator. Fairly robust 2D quadrilateral mesh generators now exist. However, the extension of algorithms that enable a 2D quadrilateral mesh to be automatically generated on complex 2D geometry to automatic 3D hexahedral mesh generation has proven to be much more challenging. Many tools will claim to have a fully automatic 3D hexahedral mesh generator, but fail when faced with the geometric complexity of biological systems.

In general, the ideal hexahedral element is a perfect cube. The element can also perform well if interior angles are right angles (i.e. a rectangle with reasonable aspect ratio (typically less than 10 to 1). When the interior angles become too acute or obtuse, the accuracy of the element degrades and many tools will warn users that the element is overly distorted. Sometimes this is not an issue, especially for linear problems in which the location of the distorted elements is not in a region of interest.

4) 20-noded quadratic hexahedral or brick element

A hexahedron has 12 edges. If a node is added to each edge of the 8-noded hexahedral element, the result is a 20-noded quadratic hexahedral element. The element has 60 degrees of freedom. The element displacement field is a complete quadratic function and also includes some additional cubic and higher-order terms. The element is able to take on the correct deformation shape in bending applications and therefore performs very well. Further, the addition of the mid-edge nodes enables the element to have a curved shape, which might be useful to model curved geometry. Virtually all FE tools will automatically generate a mesh of 20-noded quadratic hexahedral elements from a mesh of 8-noded linear hexahedral elements by inserting nodes at the mid-point on every element edge.

The main downside of the element is the large number of degrees of freedom resulting in a 20-noded hexahedral element mesh. Depending on the type of solver the FEA program uses, the result can be a solution time at least four times longer than an identical mesh of linear hexahedral elements.

 

Copyright - biomesh.org, 2007