SOFTWARE

Understanding and visualizing drug distribution and transport effects in the human body as a whole is important to both medical device and drug manufacturers. This is particularly true in the case of innovative delivery methods such as stent-based or diffusive-patch delivery.

In 2004, Teague and Feezor described the problems inherent in creating complex finite element models of biological tissue, even in terms of creating a usable geometry file.¹ They described some available methods for coping with this problem, but focused mostly on difficulties in solving mechanical problems in a single tissue.

In 2001, Hwang et al. developed a numerical model to accompany laboratory testing of bovine aorta.² They were able to capture anisotropic diffusive properties in the arterial wall, as well as variation of transport properties from one drug to another, in their model. Their study also focused on effects in one tissue.

Although extensive research has focused on single-tissue or single-organ scale computational modeling, such work is limited in scope.^2-6 More precisely targeted therapies, such as transdermal patch (TDP)–delivered muscle pain drugs, for example, may require a more macro view of the human body to ensure that the desired effect is only realized at the desired site.⁷

Ansys finite element software and its proprietary finite-volume code, CFX, have been used to model hemodynamics and pulmonary delivery problems and some biomechanics or tissue mechanics problems. However, large-scale applications have not been sufficiently pursued. The work presented here explores and demonstrates how finite element methodology can be used for large-scale convective-diffusive transport simulations and includes an explanation of element types and necessary special settings. Several examples of small-scale problems are presented to clarify the mathematical adaptations necessary to build the systemic circuit and body model.

Case Studies

This article incorporates examples of several test cases. Three of these are simple examples, which were performed in order to compare results of this method with intuitive predictions and previous models. The fourth is a novel integration of a systemic-circuit model into a full-body soft tissue model, which can later be modified to include other organs. The test cases are as follows:

• A whole-body heat conduction problem.
• A simple arterial wall diffusion problem, as if by stent delivery.
• Integration of systemic-circuit and bulk soft tissue models, omitting other structures for the time being.
• A two-part stent problem: deployment (stent and wall modeled) and diffusion of the drug coating on the stent through the arterial tissues.

Finite Element Model Formulation

The Solid 70 element is a heat-transfer element that includes mass transport (diffusion) and nonlinear fluid flow (convection) as optional settings. These are built into the Solid 70 thermal elements purely out of convenience because the equations describing these effects are remarkably similar. However, as the bulk of transport in human tissue is attributed to convection-diffusion effects, the inclusion of these three groups of equations in the same element formulation is also convenient for solving transport problems in the human body.

Special element settings (KEYOPT (8) > 0) include mass transport effects (i.e., Fickian diffusion). The software even checks Peclet (Pe) number (a ratio of convection effects to diffusion effects) automatically to ensure that convective forces are relatively small and that the solution will, therefore, be physically valid. This is useful in the solution of the majority of flow-through-tissue problems.

Normally, Solid 70 thermal elements solve the Laplace equation for heat transfer problems:

with T, temperature gradient, described by:

where k represents thermal conductivity, A is cross-sectional area of the element, and q represents heat flux. However, when KEYOPT (9)>0, Solid 70 solves the Laplace equation instead with C, the concentration gradient for Fickian diffusion, as shown below:

To produce a dimensionally correct flow rate, the normal specific heat and density material property inputs are instead set to arbitrary values of unity, allowing diffusion coefficient (D) to replace thermal conductivity as the fundamental material property. Diffusion equation formulations are usually presented in dimensionless form, allowing the variables to be adapted more easily to specific problems. So, in this case, dimensionless concentration units are exchanged for temperature in the applied boundary conditions, and diffusivity in area per unit time replaces thermal conductivity.

The software also calculates a Pe number and ensures that it is less than one, as Pe > 1 would represent relatively weak diffusion compared with convection. This ensures that the solution will be physically valid, because convective forces would otherwise overwhelm the effects of diffusion.

The same thermal element also allows activation of a nonlinear fluid flow option, which is analogous to convective transport. Activation of KEYOPT (7) allows solution of pressure-driven, convection-dominated, high-Pe transport. This is done by a similar variable analogy:

Temperature gradient and its accompanying proportional tensors from the Fourier heat equation can thus be replaced by an appropriately grouped set of tensors and pressure gradients here (k is permeability, A is cross-sectional area, and µ is viscosity), as done above by concentration gradient and diffusivity.

This means that in certain tissues such as the vitreous humor, where convection plays a much greater role (up to 40% of total mass transport in this particular case), a submodel can be used to simulate the sum of convective and diffusive forces by solving in two steps.⁸

The problem of modeling the circulatory system can be addressed as well. The software includes an element type called Fluid 116, one-dimensional (beamlike) coupled physics elements that solve pressure and temperature equations simultaneously, representing simple pipelike flows that may also be conducting heat. The temperature solution can be adapted to chemical transport problems by the same analogy discussed above. The pressure solution allows appropriate, pulsed-pressure boundary conditions for mass transfer from tissue to blood and vice versa.

Fluid 116 elements solve the following:

N is the number of flow channels, C is specific heat, Kt is thermal conductivity, and KP is pressure conductivity. In converting to dimensionless concentration units as above, specific heat must be eliminated from the calculations by similar means, and Qg must be set to zero (internal heat generation). By these means, Fluid 116 elements can complement the ability of Solid 70 elements to form soft-tissue models.

The system may be modeled such that the pressure drop across capillary beds is fully taken into account. This is done by modifying N, the number of flow channels, and by modifying the programmed radius of the channels. By reducing the diameter significantly and increasing the number of flow channels by one to two orders of magnitude, different capillary layers may be modeled. This allows a clear distinction between injection on the high-pressure side and the low-pressure side of the systemic circuit.

This is realistic because capillaries, the smallest vessels in the circulatory system, can be in the range of one red-blood-cell diameter themselves, whereas the aorta is the largest individual vessel. These capillaries are often ~2.5 cm in diameter. More than realistic, however, it is also beneficial in forcing the pressure drop necessary to distinguish effects in the low-pressure side (as a point of drug entry, for example) from those in the high-pressure side.

The pulsing pressure load applied by the heart can also be modeled. By setting time steps to roughly 30 seconds, a loop can be used to vary pressure between 120 and 80 mmHg with each time step by a simple on-off algorithm that returns 0 for an even number time step or 1 for an odd number step. This results in a heart rate of about 60 beats per minute, within the range of normal resting heart rates.

Heat Transfer and Diffusion Tests

The first and simplest verification problem was a test of Solid 70 thermal elements over a transient load case, because the nonangular surfaces (which necessitate a tetrahedral mesh) often cause significant errors in transient thermal analyses. These errors are often manifested by results that fall outside the boundaries defined in the problem.

Because an error was not observed in this case study, the body mesh was deemed viable. In fact, the inclusion of all three sets of equations as options for the Solid 70 element type was especially useful here, because the body model did not have to be remeshed after checking results. Only option settings for the existing mesh needed to be changed to the mass transport option for the circulatory system simulation to follow (see Figures 1 and 2).

Another simple verification case was a symmetric model showing three-fourths of an artery wall, with drug sources applied as if by stent coating. Results of this simple case indicated that Solid 70 elements can effectively capture the heterogeneity of final drug concentration due to anisotropic transport properties (see Figure 3) as observed in smooth muscle tissue, such as artery walls. Concentration contours begin with a nonuniform distribution on the inner surface of the artery but then merge and become a uniform radial gradient, owing to the larger diffusion coefficients DO and DZ as compared with DR.

Body with Circulatory System Model

In this case study, the Fluid 116 elements were added to the bulk soft tissue model. Because these are line elements, this step required the addition of small volumes within the existing solids, in order for the solids forming the body tissue to see the line elements as part of their boundaries. In this way, the lines to be meshed were part of the boundary of the solid that makes up the soft tissue model, ensuring temperature communication between the two element types. A further advantage of this approach is that it does not require complex operations like multiple steps, contact elements, or convection surface elements.

As seen in the isosurface contour plots in Figures 4 and 5, the model with the circulatory system has numerous advantages over the extremely coarse, single-material model used in the initial heat-conduction problem. Concentration gradients can now be seen moving to and around the heart. As the model increases in complexity, capillary beds will be included, further taking advantage of the flexibility of this element type and providing an accurate picture of redistribution by the body's active transport mechanism.

Two-Part Stent Case Study

One other example that illustrates the flexibility of these tools with regard to biomedical and biotechnology problems is the following: A simplified stent case study, where the software initially is used to solve the structural problem of the expansion of the stent into the artery wall by angioplasty balloon, and then it is used to solve the subsequent drug transport problem. This stent case represented a reality check, demonstrating that the element formulation can solve a more-traditional problem but will require refinements for large-scale studies.

In this model, a stent was expanded into the artery wall by a uniform fixed displacement, which is reasonably accurate given that angioplasty balloons are inflated with saline water solution and have fairly rigid walls. The displacement forces the stent into the artery walls and stretches the artery walls somewhat. At this point, the structural problem has been solved.

Many finite element analysis tools also have the capability to display additional information, such as exact penetration depth at points of contact and contact pressure. As an example of the utility of this information, significant contact pressure may warrant modifications to the problem, such as creating a localized increased-convection condition due to the pressure on the drug-coating layer from the stent strut, or modifying transport properties of the artery wall in regions of high distention.

In this case study, the tissue modeling method was applied. The structural elements used to model the contact problem of stent expansion were converted to thermal elements, with the appropriate settings discussed above. Contact elements were also reset from structural-contact to thermal-contact equations (see Figures 6–9).

Discussion

Although the stent study is effective proof that the problem may become manageable, the model requires a great deal of refinement. It requires the addition of layers of tissue and modeling of components through which diffusion is significantly reduced (bones, cartilage, and selective or semipermeable membranes). It also requires the addition of capillary bed models, and other, more-complicated tissues such as electrophysiological pathways.

Beyond mere geometric refinements, laboratory testing to determine diffusivity of particular tissues will be necessary for accurate modeling of soft tissue transport. Of course, diffusivity is not a material property but a relationship—a given tissue will yield different diffusivities for different compounds. However, broad guidelines could be developed to use as a range. Diffusion coefficients in a given medium often vary in an almost linear fashion with molecular weight—the smaller the molecule, the more easily it is conducted.

That said, the software's postprocessing capabilities can make results very easy to interpret and discuss, especially if one were to complicate the problem with time dependency, such as other drug inputs delivered by a balloon catheter that is removed after expansion, for example. As research continues, the spectrum of problems and accuracy of this method will improve.

Conclusion

A novel method for applying the Ansys finite element analysis code to large-scale studies of drug distribution and chemical transport in the human body has been presented.

Positive initial results suggest that in the future, incorporation of diffusivity relationships from simple laboratory testing could allow characterization of new drugs, delivered by new techniques. This will benefit new device, drug, and delivery designs in the same way that advanced simulations have done in many other industries. More sophisticated and advanced features will be incorporated into the current model in the near future.

References

1. Chris Teague and Chris Feezor, “Understanding Soft Tissue and Stent Design Behavior,” Medical Device & Diagnostic Industry 26, no. 1 (2004): 112–117.
2. Chao-Wei Hwang, David Wu, and Elazer R Edelman, “Physiological Transport Forces Govern Drug Distribution for Stent-Based Delivery,” Circulation 104, no. 5 (2001): 600–605.
3. Yu-Hua Jiang, Michael G Klein, and Martin F Schneider, “Numerical Simulation of Ca2+ ‘Sparks' in Skeletal Muscle,” Biophysical Journal 77, no. 5 (1999): 2333–2357.
4. CR Ethier, “Computational Modeling of Mass Transfer and Links to Atherosclerosis,” Annals of Biomedical Engineering 30, no. 4 (2002): 461–471.
5. G Nucci and C Cobelli, “Models of Subcutaneous Insulin Kinetics: A Critical Review,” Computer Methods and Programs in Biomedicine 62, no. 3 (2000): 249–257.
6. J Hrabe, S Hrabetova, and K Segeth, “A Model of Effective Diffusion and Tortuosity in the Extracellular Space of the Brain,” Biophysics Journal 87, no. 3. (2004): 1607–1617.
7. Tsung-Min (Richard) Hsu et al., “Transdermal Delivery of Hydrophilic Drugs: The Current Status and Potential for Future Delivery,” Drug Delivery Technology 4, no. 2 (2004): 58–60.
8. J Xu et al., “Permeability and Diffusion in the Vitreous Humor: Implications for Drug Delivery,” Pharmaceutical Research 17, no. 6 (2000): 664–669.

Doruk Seren is a student at Pennsylvania State University. Contact Seren via e-mail at [email protected]. Metin Ozen, PhD, is principal of Ozen Engineering Inc. (Moffett Field, CA).

Copyright ©2006 Medical Device & Diagnostic Industry