Introduction

The finite element (FE) method has proven to be one the most versatile numerical methods for solving computational problems in physics. Although it was initially applied to problems in mechanics, it has since its inception in the 1950s been used in virtually every field of engineering, including biomechanics and biophysics. Yet, at the time of FEBio’s initial development (ca 2005), no freely available open-source finite element software existed that specifically targeted the biomechanics and biophysics communities. FEBio was developed to fill this emerging need and it would accomplish this by focusing on features that are relevant in the field. The source code would be freely available and designed such that it would be easy for researchers to implement new algorithms (e.g. new constitutive models). A strong emphasis would also be placed on thorough documentation, support, and outreach to the community, which would make it much easier than before for researchers to develop new ideas and share them with others.

Initial Development (2005 – 2009)

In 1995 Dr. Jeffrey Weiss was involved in the development of a deformable image registration algorithm termed Hyperelastic Warping. This method uses a finite element method for solving the image registration problem. The result of this approach is a deformation that is physically realistic and can be used to calculate the strains and, assuming the constitutive properties of the materials are known, the stresses as well. This algorithm was originally implemented in NIKE3D, a nonlinear implicit finite element code for solving problems in solid mechanics, developed by the Lawrence Livermore National Laboratories (LLNL). Unfortunately, NIKE3D was not free and, as a government-owned program, the dissemination of the additions was complicated. Consequently, the algorithm was not widely used outside Dr. Weiss’s lab despite the fact that it had many valuable applications. To overcome this problem, the idea started to develop a new, custom, open-source, FE kernel that would serve as the engine that drives the future developments of the warping algorithm.

The development of this FE kernel commenced in 2005 at the same time that Dr. Steve Maas joined Dr. Weiss’s lab. However, it was quickly realized that the kernel could be used for solving traditional finite element problems as well, and in fact could replace the dependency on NIKE3D. This kernel code was greatly expanded in 2005 and 2006 and became the predecessor to FEBio. When other researchers started to express interest in this new code, the idea of FEBio was born: a new finite element code specifically designed for the biomechanics community. Dr. Gerard Ateshian became one of the earliest researchers to express interest in the FEBio project. With his expertise in mixture theory, Dr. Ateshian joined forces with Dr. Weiss and Dr. Maas and started working on the incorporation of biphasic modeling capabilities in FEBio, a capability that was of significant interest to the biomechanics community.
The first version of FEBio was released in December 2007. To support the creation of FEBio input files and the visualization of FEBio results, PreView and PostView were developed. Also as part of the first FEBio release, a forum was created where users could submit questions, report bugs or propose new features.

First Funding Period (2008 – 2012)

In 2008 the FEBio project received its first federally funded grant through the National Institute for General Medical Sciences in the U.S. National Institutes of Health. The main goal for this funding period was to expand FEBio’s capabilities for modeling problems in biomechanics. In particular, the initial implementation of biphasic (or poroelastic) materials was extended to mixtures with multiple constituents, using the general framework of mixture theory. A special emphasis was placed on solute-solid matrix interactions, including solute exclusion from a fraction of the matrix pore space (solubility) and frictional momentum exchange that produces solute hindrance and pumping under certain dynamic loading conditions. The finite element formulation implemented full coupling of mechanical and chemical effects, providing a framework where material properties and response functions can depend on solid matrix strain as well as solute concentration. Furthermore, to accommodate contact between solids and mixtures having a solid matrix, new contact algorithms were implemented that properly accounted for conservation of mass and momentum between contacting mixtures. These contact algorithms allowed both the solvent and solutes to flow across the contact interface.

Second Funding Period (2012 – 2016)

During the second funding period, the multiphasic framework was expanded to include modeling of any number of neutral or charged solutes in a porous deformable solid matrix that may carry electrical charge. Chemical reactions were also incorporated into the code to allow mass exchanges between solutes and solid-bound molecular species. Related to the expansion of the multiphasic framework, the capability to include body forces in biphasic materials was implemented, which served as a template for the implementation of active solute transport in the form of a solute momentum supply. The same reactive constrained mixture framework also allowed the modelling of damage mechanics in tissues, by modeling the permanent bond-breaking process as a reaction. The FEBio framework was designed so that users could combine any of these models as needed.

In this funding period a plugin framework was developed for FEBio that allowed users to create customized additions without the need for modifying or recompiling the FEBio source code. This plugin feature allows users to extend almost every part of FEBio, such as adding a new material model, customizing or extending the plot file output, coupling FEBio to other codes, and much more.

Initially, FEBio did not take advantage of parallel architectures on desktop computers, at least not directly. Support for parallelization was a result of using third-party parallel linear solver libraries. During this funding period, the serial parts of the code were parallelized using OpenMP. The parts of FEBio that were most time consuming were identified and, where possible, were rewritten with support for OpenMP. This initial effort resulted in some encouraging benchmarks and significant speedups for certain types of analysis.

Several quadratic tetrahedral formulations were also implemented during this grant period. These elements are an attractive alternative to linear tetrahedral elements since they maintain the advantages of tetrahedral mesh generation. In addition, they perform similarly and sometimes even better, both in terms of computational cost and in terms of accuracy, than the “gold standard” linear hexahedral element. This makes quadratic tetrahedral elements attractive for modeling the often complex geometries encountered in computational biomechanics problems.

Third Funding Period (2016 – 2000)

The advancements to FEBio during the second funding period put us in a unique position to extend the FEBio framework to simulation of compressible and incompressible computational fluid dynamics (CFD), and to couple this CFD framework to existing mixture capabilities in FEBio to enable analysis of fluid-solid interaction (FSI) problems.  The FSI framework combined the existing modeling capabilities for solid mechanics and rigid body dynamics with our CFD solver. Our novel FSI solver is based on mixture theory, modeling the fluid domain as composed of both a massless solid constituent (which drives mesh deformation) and a viscous fluid constituent.

During this funding period, we also focused on developing new algorithmic, analysis and numerical capabilities by implementing iterative solvers and preconditioner, new nonlinear solution strategies, and adaptive meshing.  Since CFD analysis produces nonsymmetric stiffness matrices, we implemented Broyden’s method to complement the existing BFGS solver that is used for symmetric matrices.  We implemented a number of different iterative linear equation solvers to accommodate the different types of physics in FEBio, including conjugate gradient, FGMRES with ILU(x) preconditioners, adaptive multigrid solvers, and physics-based solvers that rely on a Schur-decomposition of the global stiffness matrix for CFD/FSI problems.  In addition, we implemented the Jacobian-Free-Newton-Krylov method as a nonlinear solver.  For adaptive remeshing, the algorithms for remeshing and the criteria that decide where the mesh needs to be refined are implemented as separate modules, giving users maximum flexibility in defining remeshing strategies. New remeshing algorithms and criteria can be easily implemented as plugins. We use a template based refinement method for hexahedral meshes, while for tetrahedral meshes we have both a template based approach and a whole-mesh refinement strategy. The latter uses the MMG library for remeshing the domain (www.mmgtools.org), while preserving its overall geometry.

Finally, we completely overhauled our pre- and postprocessing software.  The user interfaces were redesigned using the Qt widget toolkit. PreView’s modeling capabilities were extended to include all the new analysis types, boundary and loading conditions, and contact interfaces defined by the new CFD and FSI modules. A boundary-layer meshing tool was added to create biased meshes, which are necessary for CFD and FSI modeling.  Several new visualization features were added to PostView to aid in analyzing and visualizing fluid flow simulations, including streamline and particle flow plots.

Fourth Funding Period (2020 – 2024)

As of this writing, we are beginning the fourth funding period for the FEBio R01 grant.  In Aim 1, will formulate and implement a computationally efficient damage and fatigue failure framework for fibrous tissues. Damage mechanics plays a major role in research related to experimental and computational biomechanics, including traumatic brain injury, progression of osteoarthritis in cartilage, progression of bone damage, intervertebral disc degeneration, tendon rupture, rotator cuff damage, and thromboembolisms. The study of soft tissue damage under various activities, including single traumatic events as well as repetitive loading conditions, represents an increasingly active area of research in the biomechanics community.

In Aim 2, we will expand the multiphysics capabilities of FEBio to include solute transport, reactions involving solutes, and tissue growth and remodeling in fluid domains and at their interfaces.  This is highly significant from the perspective of modeling mechanobiology because it places a particular emphasis on interfacial transport and reactive phenomena between fluid and deformable hydrated solid tissue domains. The modeling of such phenomena, which represent essential physiological mechanisms in living organisms, has been a perpetual challenge in biomechanics and biophysics.

In Aim 3, we will integrate the use image data through the entire simulation pipeline, from model setup to model validation.  Application of image data in the biomechanics simulation pipeline has predominantly been limited to image segmentation to extract geometry for subsequent mesh generation.  Other uses of image data in biomechanics simulation, such as strain tracking, deformable image registration, and specification of model inputs (e.g., material properties and fiber orientation distributions), are at least as important for research but much less commonly used across the field. This limited use is primarily because most software packages for computational mechanics have little or no capability to work with image data beyond the segmentation stage.  We will develop and integrate new tools for image analysis, visualization, and strain tracking within the FEBio analysis pipeline.  This framework will also enable future extensions to other image based methods in computational mechanics such as deformable image registration.

Impact of FEBio

The impact of the FEBio project has surpassed all of our expectations.  To gather statistics from our user-base for monitoring the impact of our project, we require users to register at febio.org in order to download the software.  FEBio 1.0 was released in September 2007.  By March 2008 (at the time of our first NIH R01 submission), we had 230 unique registrants and 450 downloads of the FEBio software.  By October 2011 (second R01 submission), we had 1202 registered users and had recorded 4492 downloads of the FEBio executables, and 794 downloads of the source code.  By October 2015 (third R01 submission), we had 5153 registered users, 24884 downloads of the executables, and 3064 downloads of the source code.  As of February 2019, there were 9,993 registered users, 53,066 downloads of the installation packages, and 4,720 downloads of the source code.  The fact that downloads track closely with the number of registrants suggests that all registrants continue to show sustained interest with every new release.

In 2012, we published a peer-reviewed article to describe the development of FEBio and the features that were available as of that date.  We ask users to cite this paper when they publish a study that uses FEBio.  This article has been cited 411 times between 2012 and February 2019 (Google Scholar), averaging ~60 citations per year. We also track publications that use FEBio (https://febio.org/publications/), including peer-reviewed journal articles (407 to date), book chapters (40) and dissertations and theses (53).  We used to track conference abstracts that report using FEBio, but when the numbers rose into the low thousands the process became unwieldy.

These statistics show that FEBio has gained widespread use within the biomedical sciences community.  In fact, several publications have appeared, authored by independent groups that have compared FEBio to commercial finite element modeling software, focusing on specific simulation capabilities.  These types of independent verifications demonstrate that the biomechanics community has embraced FEBio.  The FEBio user-base is growing at a sustained rate of approximately 1,200 registered users per year.