COMPUTER SIMULATION OF THE ELASTIC MODULUS OF HIP HEAD BONE TISSUE DURING THE HIP JOINT LOADING VARIATIONS
Yu.V.Akulich *, A.S.Denisov**, R.M.Podgayets*, A.Yu.Akulich**
*Department of Theoretical Mechanics, Perm State Technical University, 29a, Komsomolsky Prospect, 614600, Perm, Russia
**Department of orthopaedy and traumatology, Perm State Medical Academy, 39, Kuibyshev Street, 614000, Perm, Russia
Abstract: On the basis of deformation model of bone tissue adaptation the computer experiment was performed for simulation the changes in elastic modulus of human hip head bone tissue under time variations of the hip joint loading. The rate of load was reduced to 5% of the physiological one, that may be considered as a model of a bed (regimen) confinement, and then it was increased step by step until the normal level. The computing showed the bone was weakened under the load reducing, i.e. its elastic modulus was decreasing. During the slow step rise of load the elastic modulus was fully recovered, and the stresses at each stage did not exceed the limiting values. It was shown that application of the total physiological load to the weakened bone results in its failure.
Key words: internal remodeling, elastic modulus, bone tissue, remodeling stimulus, hip head
Introduction
The adequate mathematical description of internal remodeling of bone tissue considered as a result of its adaptation to load changes is one of fundamental problems of the bone tissue biomechanics. Solution of this problem has important practical applications, and makes it possible to predict the evolution of different pathological processes in bone structures after surgical correction.
In modern studies of bone tissue remodeling kinetic equations based on deformation and energy stimuli are in use [1]. Cowin and Hegedus [2] assumed that the rate of density changes in bone tissue is in a direct proportion with the difference between the actual strain and the strain at homeostatic equilibrium, this assumption was verified by tests on a functionally isolated turkey ulna preparation [3]. The first description of bone tissue adaptation around femoral implant has been performed with similar equation by Van Rietbergen, Huiskes, Weinans, et al [4]. The main obstacle for using this methodology to simulation of the morphogenesis of femur consists in the difficulty to determine the spatial distribution of natural bone density. Fyhrie and Carter [5] assumed that the rate of density changes in bone tissue is in direct proportion to the difference between actual strain energy density and its mean value over the whole area, that allowed them to construct a model of morphogenesis in the proximal part of human femur that agrees with the Wolf's law. However, the authors did not make any further progress after this conclusion.
In all the mentioned kinetic equations the density of the bone tissue was the subject of internal adaptive remodeling, meanwhile if we want to determine the stress strain conditions in the bone tissue, we ought to know the spatial distribution of its elastic modulus E. Usually the elastic modulus is considered to be directly proportional to the cubed density with a coefficient A, which magnitude depends on the parameters of bone tissue structure. The factor
A variation within the boundaries of a hip head, according to Fyhrie and Carter [6], runs into 20%. Hence, the hypothesis about the factor A constancy may cause the stress and strain fields to be in error during the simulation of bone remodeling. Therefore it is reasonable to exclude the bone tissue density from the adaptation model, and include the elastic modulus in it. Akulich, Nyashin and Podgayets [7] compared the equations of adaptive remodeling of bone tissue with elastic modulus as a subject of remodeling, and suggested the following kinetic equation with a local strain stimulus
dE(x,t)/dt = C [sf(x,0- s/^Cx)] , (1)
which was used in this paper. In the equation (1) E(x,t), Sz(x,t) are the elastic modulus and the strain intensity at a point x at an instant of time t; C is a factor of the rate of remodeling;
Si (x) is the strain intensity at the point x under conditions of homeostasis, that take place in the bone tissue under physiological load (we shall call it the homeostatic strain intensity). If the bone is subjected to several alternative loads with different probability, the actual strain intensity Si(x,t) is a result of its averaging over the whole ensemble of loads.
The adaptive remodeling of elastic modulus, regarding the equation (1), runs in the following way. Under the change of external loads, strains at every point are changing too, and the remodeling stimulus is introduced as a difference of the local values of actual and homeostatic strain. This local strain stimulus results in such variation of elastic modulus at the given point that the actual strain intensity at this point tends to the homeostatic one. For example, the dropped load causes the strain to be immediately decreased, and its subsequent rise to the homeostatic strain is performed due to decrease of the elastic modulus, that is the bone is weakened. On the contrary, the load rise results in the increased strain, and the necessity of its diminishing to the homeostatic one causes the elastic modulus to be increased, that is the bone is strengthened.
Methods
The computer simulation of the elastic modulus changes in the bone tissue of human proximal femur under time variations of hip joint loading was performed by the above-mentioned local strain remodeling stimulus, with a special attention to remodeling the hip head and neck bone tissue.
The investigation of stress-strain conditions in the proximal femur was carried out by the quasi-two-dimensional finite element model with the side cortical plates [8]. According to Weinans, Huiskes and Grootenboer [9], our model took into account two forces: the force Fi applied to the hip head at angle a with the vertical and the force F2 applied to the greater trochanter at angle P with the vertical, with three load cases (Fig. 1).
As a normal physiological load per day we assumed 6000 cycles of load applications for the first load case and 2000 loading cycles each for the second and third cases. In our model the force Fi was distributed on a arc about a quarter of the hip head circle by the cosine law, and the force F2 was uniformly distributed on the greater trochanter surface.
In our computations we used the following scheme of the load time variations. As the initial state we had stress, strain and elastic modulus patterns under the prescribed above normal physiological load (about the way of their determinations see below). In the first stage the load abruptly fell off to 5% of the normal one (e.g. bed regimen after the trauma or surgical intervention), duration of the stage was 50 tentative time units (days, for example). In the next stages, from the second to the sixth, the intensity of load stepwise increased up to
10%, 20%, 40%, 80% and 100% of the physiological load, respectively, and the length of each stage was 50 "days". The assumption was made that the time variations of loading consisted only in proportional change of force magnitudes, but their inclinations to the vertical and number of each type loading cycles remained unchanged.
Results
In order to study the process of the bone tissue adaptation under loading variations, first of all one have to know the initial pattern of bone tissue elastic modulus in the proximal femur under the physiological load. For this purpose the problem was solved regarding the bone adaptation during the morphogenesis, with the energy stimulus [5]. The obtained pattern of elastic modulus (Fig. 2) corresponded well enough with the modulus field in the hip head determined experimentally by Brown, Way and Ferguson [10]. The pattern of homeostatic strain in the hip head is shown in Fig. 3.
Our calculations showed that the remodeling process with the energy stimulus under proportionally changed forces acting on the joint resulted in the same distribution of elastic modulus, the duration of remodeling was the only thing changed. Therefore in order to study the influence of loading variations on the hip head mechanical characteristics we used then only the kinetic equation (1) with the local strain stimulus of bone tissue remodeling.
The calculated elastic modulus distributions in the hip head at different instants of time under decrease and stepwise rise of the joint loading are shown in Figs 4 - 9. During the first stage after the load abruptly fell off to 5% of the physiological one, the elastic modulus was decreasing (Fig. 4). Comparing the distributions of bone tissue elastic modulus in the weakened bone and in the normal bone (Fig. 2), and bearing in mind the pattern of homeostatic strain intensity (Fig. 3), we can see the elastic modulus was most decreased in the very loaded and initially the strongest parts of the hip head. On the contrary, the drop of the elastic modulus is comparatively small in the slightly loaded areas of subhondral zone around
a
b
c
Fig. 1. The central layer of the proximal femur finite element model under three
load cases.
a) the first load case: F = 2317 N; a = 24°; F2 = 702 N; p = 280;b) the second load case: F1 = 1158 N; a = -15°; F2 = 351 N; p = -8°; c) the third load case: F1 = 1548 N; a = 56°; F2 = 459 N; p = 35°.
■ <
miE
E3 " □ 1
Fig. 2. The elastic modulus pattern in the hip head under physiological load (in 0.1 MPa).
Fig. 3. The homeostatic strain intensity pattern in the hip head.
the head, and at the end of the first stage these areas became perceptibly stronger then the others (Fig. 4). The elastic modulus continued to be slightly decreasing at the second stage under the load increased from 5% to 10% of physiological one, and only in a few areas it began to rise (Fig. 5). Such reduction of elastic modulus continued in some areas even at the third stage, when the load was 20% of the physiological (Fig. 6), and only at the fourth stage under the 40% of physiological load the elastic modulus began to grow in all areas (Fig. 7). The elastic modulus distribution at the and of the sixth stage of the computer simulation, i.e. in 50 "days" with the full load (Fig. 9), was practically identical to the initial pattern (Fig. 2).
n m
E3 a.
E3
□
- 50.000 7 .422E+0001 'e.700E+0001 5.978E+0001 5 .256E+0D01 4.533E+0001 3 . BllE+OQOl 3.089E+0001 2 .367E+OOQ1 1.644E+0001 9 .222E+ODOO
Fig. 4. The elastic modulus distribution in the hip head at the end of the first stage after 50 "days" under 5 % of the physiological load (in 0.1 MPa).
Fig. 5. The elastic modulus distribution in the hip head at the end of the second stage after 50 "days" under 10 % of the physiological load (in 0.1 MPa).
Fig. 6. The elastic modulus distribution in the hip head at the end of the third stage after 50 "days" under 20 % of the physiological load (in 0.1 MPa).
Fig. 7. The elastic modulus distribution in the hip head at the end of the fourth stage after 50 "days" under 40 % of the physiological load (in 0.1 MPa).
Fig. 8. The elastic modulus distribution in the hip head at the end of the fifth stage after 50 "days" under 80 % of the physiological load (in 0.1 MPa).
j~~1 6,
11 5®
□
= 30D.□□□ 1.5DOE+OOD2 '1.350E+00D2 1.2DOE+OOD2 1.050E+00D2 9.□□□E+OOD1 7.5DOE+OOD1 6.ODOE+OOD1 4.5DOE+OOD1 3.□□OE+OODl
1.5DOE+OOD1
□
Fig. 9. The elastic modulus distribution in the hip head at the end of the sixth stage after 50 "days" under the normal physiological load (in 0.1 MPa).
Time variations of elastic modulus and strain intensity in some finite element in the hip head are shown in Fig. 10 and 11. We can see that elastic modulus at each stage tends to some magnitude, equilibrated for the given level of load, and the strain intensity tends to the homeostatic one.
As each step of load rise is accompanied with abrupt change of strains (Fig. 11), the stresses are increasing too. It may result in excess of strength limit, and in destruction or necrosis of bone tissue. The computations showed that the stress intensity in the hip head in the beginning of each stage of the stepwise load increase did not exceed the sponge bone strength limit of 5.5 MPa. Such stresses existed only in a stronger bone in the neck of femur. During each stage the elastic modulus was increasing under the constant load, but the strains were decreasing, so the stresses in different elements altered slightly in the upwards or downwards directions. We also considered the case, when instead of the stepwise increase the total physiological load was applied to the weakened bone immediately after the first stage. Here the stress intensity in the hip head and neck exceeded the strength limit (Fig. 12), and the deformation of the hip head became critically large (Fig. 13), that practically meant its destruction.
Conclusion
The analysis carried out in the present paper showed that the proposed phenomenological model of bone tissue remodeling describes the clinical facts qualitatively correct. In order to achieve more exact quantitative correspondence, additional investigations of the bone remodeling are required. Firstly, we assumed that the bone tissue remodeling began immediately after the loading of the joint was changed and the current strain intensity at every point deviated from the homeostatic one at this point. However, the clinical practice validates the existence of some "dead" zone, within which the bone tissue is irresponsive to variations of load and strain at least for a while. Secondly, the computations of the remodeling were performed with some tentative time units. For their reference to real time scale, the justified specification of the remodeling rate factors is required, and these factors undoubtedly must be different for the cases of bone tissue being hardened or weakened. Thirdly, we considered some arbitrary scheme of the load rise after that the bone was weakened. However, to solve the real problem of rehabilitation after the surgical intervention it is necessary find the optimal safe scheme of the load rise, during which the complete recovery of the bone tissue properties would be achieved within the shortest possible time. At last, we must remember that all the considered above characteristics are undoubtedly very individual.
2.0 n------
1.6
1.2
0.8
0.4
0.0 ------
0 50 100 150 200 250 300
Fig. 10. The time variations of elastic modulus in one of the finite elements in the hip head (the ratio of the current elastic modulus with the nominal for the sponge bone tissue magnitude of 500 MPa is represented).
0.005
0.004
0.003
0.002
0.001
0.000 ------
0 50 100 150 200 250 300
I
L
Fig. 11. The time variations of strain intensity in one of the finite elements in the hip head (the horizontal dashed line is the homeostatic strain at the given finite element).
■5-m 4.
E33-
E3 □
= 51.000
1.234E+0000 '1.lllE+OOOO 9.874E-0001 8 .642E-0001 7.411E-0001 6.179E-D001 4 .948E-0001 3.716E-0001 2.484E-0001 1.253E-0001 2.128E-0003
Fig. 12. The stress intensity in the hip head (0.1 MPa) under the abrupt increase of the load from 5% of the
physiological to the complete load.
Fig. 13. The deformation of the hip head under the abrupt increase of the load from 5% of the physiological to
the complete load.
References
1. HUISKES R., WEINANS H., GROOTENBOER H.J., DALSTRA M., FUDALA B., SLOOF T.J. Adaptive bone remodeling theory applied to prosthetic-design analysis. J Biomech, 20: 1135 -1150, 1987.
2. COWIN S.C., HEGEDUS D.H. Bone remodeling I: theory of adaptive elasticity. J Elasticity, 6: 313-326, 1976.
3. RUBIN C.T., LANYON L.E. Osteoregulatortory nature of mechanical stimuli: function as a determinant for adaptive remodeling in bone. J Orthop Res, 5: 300-310, 1987.
4. Van RIETBERGEN B., HUISKES R., WEINANS H., SUMMER D.R., TURNER T.M., GALANTE J.O. The mechanism of bone remodeling and resorption around press-fitted THA stems. J Biomech, 26: 369-382, 1993.
5. FYHRIE D.P., CARTER D.R. A Unifying Principle Relating Stress to Trabecula: Bone Morphology, J Orthop Res, 4: 304, 1986
6. FYHRIE D.P., CARTER D.R. Femoral head apparent density distribution predicted from bone stresses. J Biomech, 23: 1-10, 1990.
7. AKULICH Yu.V., NYASHIN Yu.I., PODGAYETS R.M. The phenomenological model of adaptable adult spongy bone tissue. Russian Journal of Biomechanics, 1-2: 53-57, 1998.
8. SVESNSSON N.L., VALLIAPAN S., WOOD R.D. Stress analysis of human femur with implanted charnley prosthesis. J Biomech, 10: 581-588, 1977.
9. WEINANS H., HUISKES R., GROOTENBOER H.J. Effects of fit and bonding characteristics of femoral steams on adaptive bone remodeling. J Biomech Eng, 116: 393-400, 1994.
10. BROWN T.D., WAY M.E., FERGUSON A.B., Jr. Stress transmission anomalies in femoral heads altered by aseptic necrosis. J Biomech, 13: 687-699, 1980.
Компьютерное моделирование модуля упругости костной ткани головки бедра при изменении нагруженности тазобедренного сустава
Ю.В.Акулич, А.С.Денисов, Р.М.Подгаец, А.Ю.Акулич (Пермь, Россия)
На основе деформационной модели адаптации костной ткани проведен численный эксперимент по моделированию изменения модуля упругости губчатой костной ткани головки бедренной кости человека при изменении уровня нагруженности
тазобедренного сустава. Нагрузка снижалась до 5% от физиологической, что можно считать моделью постельного режима, а затем постепенно повышалась до нормального уровня. Численный эксперимент показал, что при снижении нагрузки кость ослабляется, что выражается в уменьшении ее модуля упругости. При медленном ступенчатом повышении нагрузки модуль упругости полностью восстанавливается, причем на каждом этапе напряжения не превышают допустимых значений. Показано, что приложение полной физиологической нагрузки к ослабленной кости приводит к ее разрушению. Библ. 10
Ключевые слова: внутренняя перестройка, упругий модуль, костная ткань, стимул перестройки, головка бедра
Received 07 October 1998