Finite element analysis (FEA) is one of the mostly widely used engineering analysis techniques in the world today. Engineers employ FEA to simulate how a physical system (usually an engineered product or manufacturing process) will respond to expected loading conditions. Practical applications of FEA include structural analysis of bridges and buildings, impact or crash analysis of automobiles, aerodynamic analysis of airplanes and airfoils, electromagnetic analysis of AC and DC motors, injection molding simulation of plastic parts, fluid flow analysis in channels and pipes and heat transfer through residential and commercial buildings. Finite element analysis is based on the fundamental physical principles that govern the behavior of these physical systems. Since biological systems must obey the same fundamental physical principles, the physical response of biological systems to known loading conditions can also be predicted using FEA. For simplicity this brief introduction to finite element modeling will be restricted to elasticity or structural analysis applications.
The following sections provide an Introduction to the basic principles of finite element analysis and a summary of the Methods used to build and analyze finite element models of biological structures.
Like all analysis models, a finite element model is an abstraction of a more complicated physical system. The physical world is invariably too complex to model at every level of detail. For example, to predict how a long slender object might deflect due to an applied transverse load, the object might be modeled as a “beam” which involves a number of simplifying assumptions according to the mechanics of beam theory as to how the object will deform under an applied load. In the case of beam theory, these assumptions enable analytical equations to be derived that give the beam response to loads, such as shown in Figure 1 below.
|Figure 1: Bending of a simple beam|
However, for more complicated problems analytical solutions cannot be obtained, or if obtained, lack sufficient accuracy to be of any value. Computer-based techniques must be employed to simulate how such systems will behave (i.e. physically respond) under known imposed loads. Consider the ladder shown in Figure 2. A finite element model using “beam” finite elements is shown next to the ladder. The beam finite elements, shown in blue and purple, are represented as straight lines defined by points, called nodes, shown in green in the figure. Each beam element can have distinct material properties (e.g. Young’s modulus of elasticity), a distinct cross sectional area, and distinct moments of inertia about its bending axes. Forces and/or moments may be applied to a set of nodes and a different set of nodes may be constrained against deflections and/or rotations. A finite element analysis (FEA) program can “solve” this model and provide a “solution” that predict the deflections, moments, and stresses throughout the structure.
|Figure 2: A ladder and a simple finite element model|
|Figure 3: Examples of some planar types of finite Elements|
Here are many different types of finite elements that can be used to enable model physical systems at various levels of detail. For example, the beam element model of the ladder above does not represent geometric details such as how the rungs of the ladder are attached to the ladder rails. There are many other types of finite elements, each of which is optimal for different situations. Shell (or plate) elements are useful for modeling thin, three-dimensional structures. Planar elements are used to model systems that are assumed to respond to in-plane loading but exhibit little variation in geometry, material properties, and constraints in the out-of-plane direction (Figure 3). Corresponding three dimensional solid element versions of these elements are tetrahedral (4-sided) elements, brick (8-sided) elements, and hexahedral (6 sided) elements.
|Figure 4: Object modeled using tetrahedral elements, showing differences in the density of element packing|
Figure 4 illustrates an object modeled with tetrahedral elements. The subdivision of the object into a contiguous set of finite elements connected together by nodes is called a finite element mesh. Note that size of the finite elements can vary spatially; smaller elements (higher element density) are often used to improve the accuracy of FEA solution in regions where the stress gradients are predicted to be high. Most commercial programs contain tools to generate FE meshes with non-uniform mesh density for certain types of elements. In particular, there are robust automatic mesh generators for meshing complex two dimensional geometry with triangular plate elements, and for meshing complex three dimensional geometry using tetrahedral elements. However, the price often paid for automatically meshing three dimensional objects with simple tetrahedral elements is loss of accuracy and/or large increases in the computational resources required to solve the problem.
Automatic meshing of thin, geometrically complex, objects is a particular problem. Imagine that the object in Figure 4 were hollow with a wall thickness less than or equal to the minimum distance between any pair of adjacent nodes in the region of high element ndensity. In this case, the tetrahedral elements in the lower element density regions would have very high aspect ratios. This can significantly reduce the accuracy of the solution. Further, thin structures exhibit bending behavior. To accurately capture bending behavior with four-noded solid tetrahedral elements, the object would need to be meshed such that the average size of each tetrahedral element is 1/4th the wall thickness or less. This means many more nodes and, typically, elements would be required, which would result in a dramatic increase in analysis time and may even exceed available hardware and software resources.
|Figure 5: I-beam meshed using shell elements|
One solution to this problem would be to mesh a thin structure with shell or plate elements, which are specifically designed to accurately model bending in thin three-dimensional structures. An example of a shell element mesh of an I-Beam is illustrated in Figure 5. Here, shell elements were placed along the mid-plane surface of the I-Beam (an imaginary surface lying half way between the outer and inner surfaces), and each element was assigned a thickness to define the wall thickness through out the I-Beam. Unfortunately, most biological structures are complex, thin, three-dimensional and digitally reconstructed from 2D slices. At best, it is highly problematic to automatically determine mid-plane surfaces in these kinds of models.
The complex and free-form organic shapes of biological structures pose a challenge in finite element modeling. It is extremely difficult, if not impossible, to manually construct complex organic shapes with FE or CAD software. These engineering tools are not designed to handle the highly irregular, complex geometric shapes of organics systems. The most efficient way to generate FE models of biological systems is to begin with serial images derived using computerized topography (CT), magnetic resonance imaging (MRI), confocal microscopy, or even histological preparations.
Figure 6: A CT slice through the skull, forearms, and fingers of the wrinkle-faced bat, Centurio senex
Recent innovations have revolutionized biological imaging. It is now possible to capture serial sections of virtually any structure and generate exquisitely detailed 3D reconstructions. We use 3D surface reconstructions created from CT scans as templates for 3D finite element models. Figure 6 illustrates a typical CT slice. We work on small animals and are fortunate to have access to a micro-CT scanner in Dr. J. W. Hagadorn’s lab at Amherst College. Other well-known resources for micro-CT scanning are the Harvard Center for Nanoscale Systems and the High-Resolution X-Ray Computed Tomography Facility at the University at Austin (UTCT)
Initial 3D surface reconstructions are typically quite rough and require significant editing before they can be imported into an FE tool and successfully meshed as a finite element model. 3D surface representations can be saved in a number of different file formats depending on the software that you use. In our experience, the Standard Tessellation Language (stl) file format is the most portable and thus the easiest means of taking 3D surface representations from a reconstruction program to the next phase of FE model building – editing the 3D Image.
At this point in time, editing 3D images is the most time intensive step in building FE models of biological structures. It is critical to think though the sequence of editing tools and file formats you plan to use before starting the process. File format incompatibilities and software interoperability issues are common (if not rampant). We experimented with many different combinations of image reconstruction, surface editing, and FE analysis packages before identifying a set of file formats and programs that can be used to produce a FE model in a series of relatively simple steps.
We have had good success with keeping files in .stl format and using Mimics for building 3D surface models from ct scans, Geomagic Studio for “cleaning” the 3D images, Mimics for generating finite element meshes, and Strand7 for running finite element analyses. There are undoubtedly other suites of tools and file formats that “play well” together; finding them and determining how to use them efficiently simply requires time and a good supply of patience.
The ultimate goal of 3D image processing is to generate a “water-tight” surface model that can be imported into and successfully manipulated in FE software. To envision a water-tight model, think of a 3D surface representation as a skin that covers all surfaces of a model. In the case of a skull, this would include not only the external surfaces but also the internal surfaces such as the nasal and middle ear cavities and any sinuses or canals. In reality, the spaces between the surfaces are filled with bone. The surfaces must be continuous in order to model the volume they define using finite elements. To complete the analogy, if you could immerse the surface model in water, the spaces between the surfaces should remain dry (i.e., - the model should be “water-tight”).
3D surface representations derived from modern serial imaging techniques can be spectacularly detailed – a morphologist’s dream. The unfortunate reality is that the geometric complexity of these images are invariably too much for current FE modeling and analysis software to handle. We have found that the most efficient way to deal with these problems is to simplify the 3D representations using a combination of 3D reconstruction and reverse engineering software (i.e., Mimics and Geomagic Studio).
The most important step in the simplification process involves smoothing and removing details in selected areas of the model (Figure 7). The decision of where and to what extent to simplify a model rests with the investigator. In general, any area that is not likely to be load carrying is a good candidate for simplification. For example, in developing models of skulls we routinely remove the nasal conchae, smooth the walls of the nasal cavity, and treat the semicircular canals and much of the middle ear cavity as solid structures. We are willing to assume that these alterations will not affect the results of our loading experiments. At the same time we retain complex structures that we suspect are structurally significant. These include structures such as the sinuses, nasal septum and the pterygoid plates.
|Figure 7 (click to enlarge): Progressive simplification and decimation of the skull of Centurio senex, the wrinkle-faced bat|
Another important goal of simplification is to significantly reduce the overall size of surface representations. This vastly decreases the computational resources required to manipulate the FE model and conduct loading experiments. 3D surface representations are composed of connected polygons (triangles) and are often referred to as ‘polygon models’. The more polygons a model contains, the greater it’s fidelity to the object it represents but the larger its size. The easiest way to reduce, or decimate, the number of polygons in a model is to use the polygon reduction/decimation tools available in either 3D reconstruction or reverse engineering software (Figure 7).
Determining the extent to which a surface model should be decimated requires balancing the capabilities of the FE software against the need for a geometrically accurate model. There is a strong correlation between the number of polygons in a surface model and the number of elements in, and therefore the size of, the subsequent FE model. Therefore, it’s important to know the limits of your FE software and computational resources. (Our favorite FEA tool, Strand7 easily handles FE models containing more than 1 million elements.)
Even if computational resources are not an issue, larger models will almost certainly contain larger numbers of errors that will need to be fixed during FE mesh construction. Moreover, it is important to keep in mind that complex models are not necessarily “better” than simpler models. Finite element models are analysis models, and like any model a finite element model constitutes an idealized (i.e. simplified) representation of the physical world. The goal is not to develop the most accurate model possible. Rather, the goal should be to develop the simplest model that still represents structural relationships with accuracy and/or resolution required to answer the question of interest. Again, the decision of how simple to make a model lies with the investigator. .
Image processing is the most labor-intensive aspect of conducting FE analyses of biological structures, so choosing the right software for smoothing and decimating surface models in preparation for FE model-building is absolutely critical. One of our goals has been to streamline this process in order to make FE analysis more readily available to comparative biologists. Our software recommendations are based on our experiences in working to attain that goal.
In our experience, Mimics is a good tool for simplifying areas of complex morphology when it is useful to refer to cross-sectional images to differentiate between scanning artifacts and morphological details. We find Geomagic Studio to be indispensable when it comes to filling holes in, smoothing surfaces on, and decimating surface models. We have taken the approach of starting with as rich a data set as possible and then simplifying it as needed in each subsequent step. Therefore, we decimate our polygon models as a final step before generating finite element meshes. We use Mimics to accomplish these tasks. Finally, we import finite element meshes and complete analyses in Strand7.
Finite Element Modeling
The modeling process is made up of 4 important steps:
- Importing 3D representations into FEA tools
- Assigning material properties
- Constraining the model
- Loading the model
The ease or difficulty of importing a 3D surface model into FEA software and generating an FE mesh depends upon the quality of your 3D surface model and your choice of FE software.
The first hurdle in importing a 3D surface model is file format compatibility. Most FEA tools can import water-tight surface meshes that are saved in stl (Standard Tesselation Language) format. STL models are typically recognized as plate models from which a solid mesh can be generated. Many FE tools also have the ability to import models proprietary Computer Aided Design (CAD) formats. This is useful if you are building FE models with CAD software rather than directly from surface representations derived from CT-scans. CAD models are generally recognized by FE software as “geometric models”, which means that their shapes are defined by mathematical formulae (parameterized). Biologists are often surprised to learn that models dervied from serial images do not have "geometry" because they are not parameterized.
Once an STL surface model is imported successfully into FE software, flaws in the surface model can cause problems. These flaws can include poorly formed plates (angles excessively acute or obtuse), overlapping plates/nodes, t-junctions between plates, duplicate plates/nodes, and small gaps between adjacent plates. In our experience these kinds of errors are common in surface models built from CT-scans. In the worst-case scenario, some FE software is very sensitive to errors and simply can’t open models containing errors. In other cases, stl files with errors can be opened but attempts to generate a solid mesh from them will fail. The key to solving these problems is to use FE software that can either fix or help to pinpoint the problems. Mimics, and to some extent Strand7, have built-in mesh cleaning tools that remove duplicate plates/nodes, zipping together adjacent plates, eliminate t-junctions, and control the aspect ratios of plates within the surface models. Once a solid mesh is generated and an FE analysis is underway, Strand7 (and most other FE tools) provides warnings to identify malformed element in the mesh. If they cause significant errors, they must be fixed by making changes to the STL model and generating a new mesh.
Material property information for the system being modeled must be known a priori. It is an input - not an output - to a finite element analysis.
Any FE model can include elements with different material properties. Therefore it is possible to model virtually any complex composite structure. Most commercial finite element tools support the specification of many different material behavior models. These include: linear elasticity, nonlinear elasticity, nonlinear inelastic, viscoelastic, isotropic, orthotropic, anisotropic, etc. Some tools (including Strand7) even permit the user to write macros or subroutines for the purposes of defining unique material models. Currently, no commercial FE software provides default material properties databases for tissues such as bone, ligament, cartilage or muscle. We have assembled a database of material properties for biological tissues in the "Materials Database" section of this website.
All finite element models must be kinematically constrained before they can be analyzed. In rigid body static analysis, we use Newton’s laws to deduce the set of forces necessary to keep a body in static equilibrium and often graphically represent this with a rigid “free body diagram.” The deformation of the body is not of interest in rigid body static analysis. In contrast, predicting deformation is the goal of finite element analysis and it is mathematically impossible to predict deformation without imposing kinematic constraints. FE models must be fixed in space in order to predict how they will deform under load (i.e. before one can compute stresses and strains). For FE analysis this means that there must be sufficient constraints on a model to prevent all possible modes of rigid body motion. A rigid body in space has 6 rigid body modes of motion - translation in each of three mutually orthogonal directions and rotation about the axis that define each of the three mutually orthogonal directions.
Choices about constraints can have a significant impact on the patterns of stress and strain predicted through FE analyses. Figure 8 compares two FE analyses of a skull that differ only in how constraints were applied. The skull on the left was constrained at 3 nodes: one at each of the centers of the left and right jaw joints and another on the upper left molar. The model on the right was constrained at the same node on the upper left molar and at one node on each of the occipital condyles. In both cases, the node at the tooth was constrained against all movement while the other two nodes were free to rotate but fixed against displacement.
|Figure 8: Two FE analyses run on the same model, differing only in how constraints were applied|
There are clear differences in both the magnitude and distribution of stress under the two loading regimes. If one is interested in stresses generated during feeding, we would argue that the constraints on the left (one at each jaw joint and one at the tooth) are “better” because they more closely model how the lower jaw contacts the skull.
Note that defining kinematic constraints, like model geometry, is another example of an analysis idealization. When a single point defined by a node is rigidly constrained against motion, stresses predicted at this point are not realistic. The constraint results in what is called a “stress singularity.” This means that if a series of analyses were conducted with increasing mesh density surrounding the constrained node, the predicted stresses at this point would keep increasing without bound. The singularity in stresses merely reflects the physical impossibility of constraining a real physical object at a point. Where to place constraints on an FE model is clearly an important decision that should be made carefully and with as much prior knowledge about how a system works as possible.
Loading the model
The final step required for FE analysis is to apply loads to the model. As the primary goal of FE analyses is to predict how systems will perform, it is important that loads be applied as realistically as possible. For example, engineers who want to use FEA to predict how an earthquake will affect a new building design would need to define the amplitude, frequency and direction of shock waves from a seismic event. Similarly, a biomechanist interested in the peak compressive stresses within the femur of a runner would need to define the magnitude and direction of forces acting on the bone. FE tools offer many different ways to apply loads. These include the ability to apply loads as forces acting on individual nodes (point-loads), forces distributed over groups of nodes or element surfaces, pressures, loads per unit volume such as due to gravity or acceleration, or initial stresses or strains that result in nodal loads. The problem in modeling biological systems isn’t so much what method to use, but knowing how the loads should be applied.
In contrast to many man-made structures, the forces acting on biological systems can be extremely complex. Take, for example, the relatively simple situation of a muscle applying a load to a bone. For an FE model, a researcher would want to know something about the areas over which the muscle is attached, the magnitude of the force produced by the muscle, the direction in which muscle force is applied, and loads imposed by ground reaction forces (in the case of locomotion) or contact with other bones.
A great deal of thought and, in some cases effort, need to be put into loading FE models of biological systems. In our lab we have addressed the problem of loads imposed by muscles that wrap around their surfaces by writing a program. The program requires basic input regarding muscle forces and directions and can apply a combination of tangential tractions and normal loads due to muscle attaching to and wrapping around irregularly shaped bone surfaces. You can download the program from the "Resources" page of this website.
Once an FEA solution has been obtained, it is critical to assess its accuracy. Inaccuracies due to incorrect modeling assumptions are difficult if not impossible to assess without valid empirical data or expert knowledge about the analysis.
|Figure 9: Bar affixed at the left face under a distributed axial load|
Other sources of errors can be effectively assessed and controlled. As a numerical method, finite element analysis results are approximations of theoretically exact (and usually intractable) solutions to the underlying equilibrium equations that govern the physical behavior of a system. In reality, the continuum of a system has an infinite number of degrees of freedom. The discretization of the continuum with a finite element mesh results in a finite number of nodal degrees of freedom that approximates the solution behavior.
Using separate analyes of a series of properly constructed meshes with increasingly more nodal degrees of freedom, one can come arbitrarily close to the theoretically exact solution. This process of controlling the inherent finite element discretization error is called a convergence study. Figure 9 illustrates a simple bar that is fixed at the left face and under a distributed axial load (i.e. load per unit length) Below, Figure 10 shows the uniaxial deflection and stress in the bar for the exact solution compared to a single finite element solution.
|Figure 10: Comparison of uniaxial deflection and stress for the exact solution versus the one-element solution|
Figure 11 compares the exact solutions to the finite element solutions when the bar is meshed with either two or four elements. Note that the finite element solution approaches the exact solution as the bar is meshed with more elements, and hence more nodal degrees of freedom. (Also note that the displacement solution is more accurate for a given mesh than the stress solution. This is typically the case.) Typically, one expects a quantity of interest, such as maximum deflection, principal strain, or von Mises stress, to converge asymptotically to the theoretically exact value. This means that to demonstrate convergence to an asymptotic value, at least three different analysys are required to generate a plot of the quantity of interest versus nodal degrees of freedom.
|Figure 11: Comparison of exact solutions versus finite element solutions for a bar meshed with two and four elements|
Given the large number of nodes and elements required to model some biological systems, it may not be possible to develop a series of models with increasing mesh refinement. While convergence studies and comparison to experimental data are always preferred methods of validating FE analyses, we can still look for evidence that confirms results derived from single analyses. For example, one can examine the results of a finite element model for lack of local equilibrium. On a free external surface (i.e. an external surface with no forces acting on it), the analysis should predict that the component of stress normal to the surface is essentially zero. Similarly, although each finite element admits different states of stress, certain stress components must be continuous across inter-element boundaries in a theoretically exact solution. Inspecting the lack of continuity in certain stress components across inter-element boundaries is another way to help assess how accurate the model is.
Copyright - biomesh.org, 2010