A finite element modeling of the human lumbar unit including the spinal cord ()
1. INTRODUCTION
Spinal cord injury (SCI) is usually a consequence of a traumatic movement resulting in a deformation of the neural and vascular structure of the spinal cord exceeding their structural and physiological limits. As a matter of fact, the neurological damage is due to different loading combination. The identification of the parameters, leading to neurological deficit during walking (Neutrogenic claudication) can help to understand and prevent all these phenomena. Based on clinical observations, the physiological changes of soft tissues throughout life is another factor to be taken into account. They could reduce the available space for the neural structures (spinal cord, medullar cone and roots) [1]. This effect is amplified by vertebrae slip due to the weakness of vertebral facets, ligamentous ossification, vertebral joint fusion and also osteoporosis. This is emphasises the need of the correlation between the vertebral elements pressure and the one taken on the spinal cord element.
In general, the increase of stenosis during walking and prolonged standing position leads to neutrogenic claudication and in common to sphincter disorder. In some cases, these degenerative lesions exist whereas the clinical phenomena did not appear during rest.
In medical revues, the spinal component displacement (mainly the intervertebral disc) have been often related to neurological deficits or complications after lumbar spine manipulations. These observations have been made in [2]. Also in [3], the authors showed that the compression of the spinal cord of a fresh corpse produces a transverse crack in the spinal cord. The authors demonstrated that the damaged surfaces of this crack deviate. It is then essential to connect biomechanical observations obtained through computations to these various physiological observations.
From biomechanical point of view, many authors provided some answer to how different components of the spinal column interact with each other under different loading situations but none of these works examined the spinal cord behavior in the lumbar unit. One example of these works is the one presented in [4] where a three dimensional nonlinear finite element analysis of the mechanical behavior of the L2-L3 disc-vertebrae unit is performed with a focus in the use of tissue engineered intervertebral discs under complex loads. Also, in [5], authors presented an original work on the formulation of wrapping elements sliding over solid body edges and used to study the loadbearing capacity of simplified beamrigid body thoracolumbar and lumbosacral) spines. In [6], the authors presented a two-dimensional axisymmetric finite element model of a spinal motion segment consisting of the first lumbar vertebral body and adjacent intervertebral disc to allow the inclusion of a centrally located tumour in the vertebral body. In [7], the author presented a three-dimensional finite element model to study the static and vibrationnal loads and to analyze the stress distribution on the intervertebral disc of a L4-L5 lumbar unit. In their work taken from [8], the authors presented recent experimental evidences concerning the distribution of forces and moments acting on the lumbar spine and neural arch. The authors mainly conducted experimental works on cadaveric motion segments to investigate the distribution of the stresses on the lumbar spinal unit. All these works detailed the behavior of different components of the lumbar unit without including the spinal cord. The major article related to the study of the spinal cord is the one presented in [9]. In this article, the author stated that the loss of normal pattern of spinal motion causes pain and/or neurological dysfunction and that the system of spine may be divided into three subsytems 1) the spinal column; 2) the spinal muscle; and 3) the neural control unit. The authors precised that the two first subsystems have been widely studied in the litterature, wheareas, the neural zone was found more sensitive and is less studied in the biomechanical litterature.
Finally, the biomechnical study of cervical flexion myelopathy throughout a finite element model presented in [10] is another important work found in the litterature dealing with the computations of stresses on the spinal cord. However, the presented results have been obtained for the upper part of the spine (the cervical spine) whereas our model deals with the lower part of the spine (the lumbar vertebrea).
A numerical biomechanical model will then be presented in this article to describe the mechanical behavior of a complete functional human spinal unit; this anatomical model consists in five lumbar vertebrae, four intervertebral discs, the physiological ligaments, capsular articular parts, and spinal cord (terminal cone).
The biomechanical study of this model has been established to study the influence of the stresses on all the elements of this model and to explain advanced state of neurological deficit signs.
2. MATERIAL & METHOD
2.1. Finite Element Method
Three dimensional finite element model is built up using the computed tomography (CT) of the L1-L5 lumbar unit including mechanical properties designed for the five vertebrae, the four intervertebral discs, the ligaments (anterior, posterior, flavum, interspinous and supraspinous), articular and capsular parts, and spinal cord (cauda equina).
The numerical modelling based on a finite element methods (FEM) simplifies the structure whether it is anatomical or not by reflecting its mechanical properties. This method requires specification of the geometry of the modeled structure, the loads and pressure applied to that structure, and the elastic properties of the components. The geometry is subdivided into small regions (elements) and the differential equations governing the deformation of solids are numerically solved. Computed quantities include local deformations in response to the applied loads, as well as the corresponding stresses.
Three steps are needed to accomplish these tasks:
1) Definition of the geometry of the column constitutive parts;
2) Establishement of the laws which govern the behavior of every part of the lumbar spinal unit;
3) Evaluation of the model by performing a series of numerical computations.
2.2. Geometric Model
Concerning the first step, the total obtained meshed model representing the lumbar is given in Figure 1. The geometry of the vertebrae have been obtained by segmentation of series of images using SliceOMatic 4.3. In the Finite Element software Ansys, several procedures are developed to make the change of reference and to bring the vertebrae in the same reference. The volumes of each vertebra are assembled.
The discs have been added to fit between two vertebrae as well as the totality of the physiological ligaments. Moreover, as it is suggested in most described physiological cases, the spinal cord is added to fit into the medullar canal which is between the first and the second
Figure 1. The obtained meshed lumbar unity.
lumbar vertebra. The annulus fibers and the cartilaginous plates are respectively shown in the same figure. Details of the meshing of vertebrae are shown in Figure 2(a). As we can see, the vertebrae is reproduced by using the scanned images.
The appropriate meshing for the disc is shown in Figure 2(b) showing the geometry of the nucleus and the annulus. Finally, the fitted spinal cord between the first and the second lumbar vertebra is shown in Figure 2(c).
2.3. Material Properties
The second step consists in modeling the vertebrae as an elastic orthotropic structure with Young’s Modulii and Poisson’s ratios obtained through the mineral bone density. The maximun and the minimum values of the the young modulii (Ex, Ey, Ez) taken in the three directions (x, y, z) are summeraized in the Table 1. The same table gives the values of the shear modulus (Gxy, Gyz, Gxz) and
(a)
(b)
(c)
Figure 2. (a) The meshed geometry of L1; (b) Geometry of the nucleus and annulus; (c) The spinal cord fitted in the medullar canal.
the three Poisson ratios (vxy, vxz, vzx). For the intervertebral disc which is divided into the annulus and pulposus part, Young’s Modulii and Poisson’s ratios were taken in the litterature [11]. Cartilages fibers are also taken into account. The material properties are specified in the same table. All these parts are simplified as a structure with an isotropic and linear elastic behavior. According to the work presented in [12], the spinal cord tissue has roughly analogous rheological properties as ligamentous tissue. An average elastic modulus of the spinal cord described in [13] is found to be equal to 1.40 MPa. An isotropic elastic behavior is also established for this physiological part.
3. RESULTS AND DISCUSSION
3.1. Boundary Conditions
The third step is performed through different loadings using the finite element method, FEM. The boundary
Table 1. Material properties for the lumbar unit L1-L5.
conditions are specified in Figure 3.
The boundary conditions for the different cases use a fixed bottom of the lumbar human unit (in the three directions) and the appropriate load for each presented case.
For the case of compression, a pressure is assigned to the top surface of the geometry to ensure the homogeinity of the applied load. The pressure is computed by dividing each load by the area located on the first lumbar L1.
Equivalent Stress
The equivalent von Mises stress given in the following equation is used to evaluate the maximum value. This value is taken as a limit to the ultimate elastic stress for most structural computation, Equation (1).
(1)
where
, 
and 
are the principal stresses.
The results will be given for the whole structure. Results obtained for the bony structure (vertebrea) and the intervetbral disc are presented here after followed by the maximun (critical values) of stresses obtained for the spinal cord.
3.2. First Load Case: Compression
Figure 4 shows the details of a subdivision of the vertebra into six different points. The three first ones are related to the posterior region (POST), respectively external (out), medium (mid) and internal (in). In the same way, three points are taken in the anterior region of the vertebra (ANT), respectively, external (out), medium (mid) and internal (in).
Figure 5(a) gives the average of the von Mises equivalent stresses 
obtained for the five lumbar vertebrae for these specific six points in the case of a load of 600 N. It is clear that the stresses are higher for more important external loads and that the most loaded vertebra is the bottom one (L5). The peaks of stresses are located at the outer posterior points.
The disc is also split into six different parts as done for each vertebra. As known, the disc is composed of nucleus pulposus and annulus fibrosis. In the normal healthy disc, the pulposus hydrated core exerts a hydrostatic pressure (pressure intradiscal (IDP) on fibers of the annulus fibrosis [12].
First Case Data Analysis
Figure 5(b) gives the average values of 
for the various points of the disc located between the 4th and the 5th vertebra (L45). The values of 
indicate that the nucleus pulpous can absorb some compressive loads and that the maximum values are obtained for the outer posterior and anterior points, which are located in the annulus region as pointed out in [14]. The same authors showed