**Estimating Lower Limb Skeletal Loading**

### Timo Rantalainen and Adam Kãodowski

*Department of Mechanical Engineering, Lappeenranta University of Technology, Department of Health Sciences, University of Jyväskylä Finland* 

### **1. Introduction**

Osteoporosis, accidents and subsequent bone fractures cause suffering on an individual level, as well as an economical burden to the society (Ortiz-Luna et al., 2009; Stevens & Olson, 2000). It has been estimated that, in Finland alone, between 30,000 to 40,000 osteoporosis-related fractures occur annually and that 400,000 Finnish people have osteoporosis (Duodecim, 2008). There are a few potential ways of preventing bone fracture, i.e. strengthening bones and/or preventing falls (Ortiz-Luna et al., 2009; Stevens & Olson, 2000). In order to withstand prevalent loading without breaking; while remaining relatively light in weight to allow for locomotion, bones have the ability to adapt their structure to functional loading (Frost, 2000; 2003; Sievänen, 2005). It has been demonstrated that physical activity affects the weight bearing skeleton more than the non-weight bearing one (Mikkola et al., 2008), and it may therefore be argued that, the skeleton is loaded mainly by locomotory actions that impart strains on bones. Bones are loaded in daily activities by muscles accelerating and decelerating body segments and resisting the pull of gravity (Burr et al., 1996). Since falling is the single most significant bone fracture risk factor (Järvinen et al., 2008) and up to 90% of fractures are caused by falls (Cummings & Melton, 2002; Stevens & Olson, 2000; Wagner et al., 2009), exercise can be viewed as a potential intervention for fracture prevention. Exercise seemingly has a potential of both reducing the fall rate and also increasing bone strength. In agreement, exercise interventions have been shown to successfully decrease the fall rate (Kemmler et al., 2010; Korpelainen et al., 2006), to strengthen the bones and to decrease the fracture rate (Korpelainen et al., 2006; Sinaki et al., 2002).

It is quite obvious, that in some cases even the strongest bone could not withstand the loading associated with falling. Therefore, it is reasonable to question whether strengthening bones makes a significant contribution in preventing fractures. A prospective study has shown that people with higher calcaneal bone mineral density had lower fracture rate even while the fall rates between the groups were relatively similar (Cheng et al., 1997). The differences between the groups in bone mineral density were relatively large (~10%). In agreement, prospective studies have shown that increasing the amount of bone mineral with drugs by about 10% is effective in decreasing fracture rates (Cummings et al., 2009). On the other hand, the increments in bone mineral amount associated with year long exercise interventions are relatively modest, i.e. in the order of 1 - 2% (Nikander et al., 2010). There are examples of fractures for which this kind of relatively modest bone strength increase does play an important role. One of the most convincing examples is vertebral fractures. It has been shown in a prospective study that vertebral bone mineral density (represents

estimated by calculating an osteogenic index, OI, (Fig. 2) (adapted from (Turner, 1998; Turner

Estimating Lower Limb Skeletal Loading 245

*i* is the number of loading cycles, and *t* is the interval between loading expressed in hours.

Fig. 2. Osteogenicity of loading is linearly related to strain magnitude and rate, while loading

However, for human applications, bone strains are generally not known. Furthermore, measuring bone strains in vivo is invasive and limited to superficial bone sites (Hoshaw et al., 1997; Lanyon et al., 1975; Milgrom et al., 2000), thus not feasible in practice. Consequently, indirect estimates have to be used in practical applications. Since, intuitively there is an association between ground reaction force and bone strains, it has been suggested that osteogenic index may be calculated from reaction forces (Turner & Robling, 2003; Whalen et al., 1988). Subsequently, reaction force based osteogenic index has been shown to be able to successfully differentiate two exercise regimes producing differing bone responses from each

Although, using reaction force as a substitute for bone strains is tempting there are some shortcomings of this approach which will be demonstrated later under heading 3. There are at least a few examples of osteogenic exercise interventions reported, where bone response has been observed on one bone of the lower limb (and spine), while other lower limb bones have not responded (e.g. Vainionpää et al. showed increase in femoral neck bone strength following a high impact exercise intervention in post menopausal women, while tibial mid-shaft bone strength was unaffected by the same intervention (Vainionpää, 2007; Vainionpää et al., 2005). While at first glance it may appear that the aforementioned is unexpected, the phenomenon may be explained by acknowledging that bone response is site specific. High impact exercises are characterized as high impact by the nature of causing high reaction forces. High reaction forces occur in exercises in which the lower limb is mainly compressed (e.g. above 10 times body weight has been measured under one foot in triple jump 2nd jump eccentric phase (Perttunen et al., 2000)), while the leg is straight. This kind of posture enables high reaction forces, since the inertia of the body is efficiently transmitted through hip, knee and ankle joint without the joints giving. If force is not transmitted through the joint centres, the joints give and reaction forces above 6 times body weight are rarely recorded (Ishikawa et al., 2005) even under two feet. In such situation the kinetic energy of the body is damped by muscles and tendons. Bone loading condition (compression or bending) has probably a major effect on bone strains in long bone mid-shaft (Heinonen et al., 1996). Consider a 41 cm long hollow

cycle dependency and interval between loadings dependency are non-linear.

other in a prospective human study (von Stengel et al., 2007; 2005).

*ln*(*i* + 1)(1 − *e*

interval.

<sup>−</sup>*t*/6) (1)

*<sup>δ</sup><sup>t</sup>* denotes the mechanical strain rate in [Hz],

(b) Sensitivity as a function of loading

*OI* = *�*

Where *�* is the magnitude of mechanical strain, *δ�*

(a) Osteogenicity as a function of

*δ� δt*

& Robling, 2003)) Eq. 1.

repetitions.

bone strength indirectly) may be increased with back strengthening exercise program and that the higher vertebral bone mineral density is also associated with lower incidence of vertebral fractures (Sinaki et al., 2002). Furthermore, exercise induced bone gains may be more appropriately situated in order to maintain bone strength compared to drug induced bone gains. In case of drug therapy, the bone mineral amount increment is distributed rather evenly through the bone1, while physical activity associated bone increments are situated mechanically advantageously (Leppänen et al., 2010; Ma et al., 2009; Rubin & Lanyon, 1984). In pure compressive or pure tensile loading, the geometry of the bone makes no difference to bone strength (excluding slender structures such as the ones found in trabeculae, which are susceptible to buckling). The only important cross-section parameter is its area. In bending and torsion, however, the geometry of the bone makes a big difference. As a rule of thumb, the further away from the bending/rotation axis the material is situated, the more significant contribution it makes to the strength of the bone (further details given under section 2.1) (Fig. 1). Coincidentally, the toughest part of bone is the external shield, which is also subjected to the highest loads.

Fig. 1. A) Stress (*σ*) distribution in a bone subjected to bending moment (M), B) Hollow oval model of tibia mid-shaft cross-section, C) Possible change of tibia mid-shaft cross-section due to exercise.

Since exercise interventions are effective and feasible, how does one optimize the osteogenic loading of the intervention? Strains (Fig. 4) are one of the most important stimulants to bone adaptation Lanyon (1987); Turner (1998), and therefore designing osteogenic interventions could benefit from the knowledge of bone strains at different cross-sections in a wide range of exercises. General rules for bone adaptation are relatively well established, i.e.; bones respond increasingly with increasing 1) strain rate, 2) strain magnitude and 3) number of loading cycles (Lanyon, 1987; Turner, 1998; Whalen et al., 1988). In addition, loading direction affects the response, i.e.; unusual loading direction causes a larger response (Rubin & Lanyon, 1984). In bone adaptation, also the interval between loadings plays an important role. Bone responsiveness to loading decreases relatively rapidly, and dividing the loading bout into smaller bouts separated by several hours is more osteogenic than doing the exercise all in one bout (Umemura et al., 2002). If bone strains are known, the osteogenicity of loading may be

<sup>1</sup> At the moment of writing this chapter, we were unaware of a published article we could cite to confirm this statement. However, the FREEDOM study(Cummings et al., 2009) did include computed tomography measurements, the results of which were presented at least at the ECCEO-IOF11 congress held March 23-26 2011 in Valencia, Spain. Although not mentioned in the presentation, one of the authors of the present chapter was in the audience and queried the presenter (prof. Genant) afterwards whether the bone geometry had changed as a result of the intervention. As one might have expected for an anti-resorptive drug, the geometry had not changed according to the presenter.

2 Will-be-set-by-IN-TECH

bone strength indirectly) may be increased with back strengthening exercise program and that the higher vertebral bone mineral density is also associated with lower incidence of vertebral fractures (Sinaki et al., 2002). Furthermore, exercise induced bone gains may be more appropriately situated in order to maintain bone strength compared to drug induced bone gains. In case of drug therapy, the bone mineral amount increment is distributed rather evenly through the bone1, while physical activity associated bone increments are situated mechanically advantageously (Leppänen et al., 2010; Ma et al., 2009; Rubin & Lanyon, 1984). In pure compressive or pure tensile loading, the geometry of the bone makes no difference to bone strength (excluding slender structures such as the ones found in trabeculae, which are susceptible to buckling). The only important cross-section parameter is its area. In bending and torsion, however, the geometry of the bone makes a big difference. As a rule of thumb, the further away from the bending/rotation axis the material is situated, the more significant contribution it makes to the strength of the bone (further details given under section 2.1) (Fig. 1). Coincidentally, the toughest part of bone is the external shield, which is also subjected to

Fig. 1. A) Stress (*σ*) distribution in a bone subjected to bending moment (M), B) Hollow oval model of tibia mid-shaft cross-section, C) Possible change of tibia mid-shaft cross-section due

Since exercise interventions are effective and feasible, how does one optimize the osteogenic loading of the intervention? Strains (Fig. 4) are one of the most important stimulants to bone adaptation Lanyon (1987); Turner (1998), and therefore designing osteogenic interventions could benefit from the knowledge of bone strains at different cross-sections in a wide range of exercises. General rules for bone adaptation are relatively well established, i.e.; bones respond increasingly with increasing 1) strain rate, 2) strain magnitude and 3) number of loading cycles (Lanyon, 1987; Turner, 1998; Whalen et al., 1988). In addition, loading direction affects the response, i.e.; unusual loading direction causes a larger response (Rubin & Lanyon, 1984). In bone adaptation, also the interval between loadings plays an important role. Bone responsiveness to loading decreases relatively rapidly, and dividing the loading bout into smaller bouts separated by several hours is more osteogenic than doing the exercise all in one bout (Umemura et al., 2002). If bone strains are known, the osteogenicity of loading may be

<sup>1</sup> At the moment of writing this chapter, we were unaware of a published article we could cite to confirm this statement. However, the FREEDOM study(Cummings et al., 2009) did include computed tomography measurements, the results of which were presented at least at the ECCEO-IOF11 congress held March 23-26 2011 in Valencia, Spain. Although not mentioned in the presentation, one of the authors of the present chapter was in the audience and queried the presenter (prof. Genant) afterwards whether the bone geometry had changed as a result of the intervention. As one might have expected

for an anti-resorptive drug, the geometry had not changed according to the presenter.

the highest loads.

to exercise.

estimated by calculating an osteogenic index, OI, (Fig. 2) (adapted from (Turner, 1998; Turner & Robling, 2003)) Eq. 1.

$$OII = \epsilon \frac{\delta \varepsilon}{\delta t} \ln(i+1)(1 - e^{-t/6}) \tag{1}$$

Where *�* is the magnitude of mechanical strain, *δ� <sup>δ</sup><sup>t</sup>* denotes the mechanical strain rate in [Hz], *i* is the number of loading cycles, and *t* is the interval between loading expressed in hours.

Fig. 2. Osteogenicity of loading is linearly related to strain magnitude and rate, while loading cycle dependency and interval between loadings dependency are non-linear.

However, for human applications, bone strains are generally not known. Furthermore, measuring bone strains in vivo is invasive and limited to superficial bone sites (Hoshaw et al., 1997; Lanyon et al., 1975; Milgrom et al., 2000), thus not feasible in practice. Consequently, indirect estimates have to be used in practical applications. Since, intuitively there is an association between ground reaction force and bone strains, it has been suggested that osteogenic index may be calculated from reaction forces (Turner & Robling, 2003; Whalen et al., 1988). Subsequently, reaction force based osteogenic index has been shown to be able to successfully differentiate two exercise regimes producing differing bone responses from each other in a prospective human study (von Stengel et al., 2007; 2005).

Although, using reaction force as a substitute for bone strains is tempting there are some shortcomings of this approach which will be demonstrated later under heading 3. There are at least a few examples of osteogenic exercise interventions reported, where bone response has been observed on one bone of the lower limb (and spine), while other lower limb bones have not responded (e.g. Vainionpää et al. showed increase in femoral neck bone strength following a high impact exercise intervention in post menopausal women, while tibial mid-shaft bone strength was unaffected by the same intervention (Vainionpää, 2007; Vainionpää et al., 2005). While at first glance it may appear that the aforementioned is unexpected, the phenomenon may be explained by acknowledging that bone response is site specific. High impact exercises are characterized as high impact by the nature of causing high reaction forces. High reaction forces occur in exercises in which the lower limb is mainly compressed (e.g. above 10 times body weight has been measured under one foot in triple jump 2nd jump eccentric phase (Perttunen et al., 2000)), while the leg is straight. This kind of posture enables high reaction forces, since the inertia of the body is efficiently transmitted through hip, knee and ankle joint without the joints giving. If force is not transmitted through the joint centres, the joints give and reaction forces above 6 times body weight are rarely recorded (Ishikawa et al., 2005) even under two feet. In such situation the kinetic energy of the body is damped by muscles and tendons. Bone loading condition (compression or bending) has probably a major effect on bone strains in long bone mid-shaft (Heinonen et al., 1996). Consider a 41 cm long hollow

**2. Stress and strain**

(Currey, 2001; Turner & Burr, 1993).

Fig. 3. Schematic illustration of a load-deformation curve.

Bending stress in a bone cross-section in a given loading configuration is determined by the second moment of area (also known as cross-sectional moment of inertia) of the cross-section as shown in (Fig. 1 A). Stress (indicated by the shaded area) is distributed unevenly throughout the cross-section, being by definition zero at the neutral axis and increasing its magnitude linearly while moving away from the neutral axis. The stress (*σ*) at a given point

**2.1 Example of how geometry affects bone strength**

It needs to be mentioned that loading, forces, stresses and strains are used somewhat interchangeably in the present text, while generally, forces, stresses and strains should be differentiated. Therefore, the remainder of the paragraph tries to explain what are the differences between forces, stresses and strains and why they should not be mixed. The relationship between the applied load (=force) and relative deformation (Fig. 4) caused by the load can be plotted as a load-deformation curve. The curve may be divided into two parts, the linear elastic region and the non-linear plastic region. The stiffness of the structure is the slope of the elastic region, yield is the load at transition from the linear to elastic region and toughness is the area under curve up to yield or failure (Fig. 3) (Currey, 2001; Martin & Burr, 1989; Turner & Burr, 1993). The load-deformation curve may be converted to stress-strain curve by accounting for the geometry of the structure and the loading situation. The strain by definition is the deformation divided by the original dimension. If the loading force is compressive (or tensile), then the stress is simply force divided by the cross-sectional area of the structure. In bending and torsion, the stress varies throughout a cross-section. The slope of the stress-strain curve in turn is the Young's modulus, i.e. stiffness, of the material. The strength of the structure is defined as the load at which the structure either yields (yield strength) or breaks (breaking strength, which equals ultimate strength in bone). As was the case with stiffness, the strength may be reported as a material property (yield or breaking stress) or structural property (yield or breaking load). Because of the anisotropic structure of bone material, Young's modulus and breaking stress are different in longitudinal and transverse directions. Furthermore, bone ultimate stress is different in tension, compression and shear. While looking at different types of loading: bending, compression and shear, it becomes clear that the geometry of bone is of utmost importance in bending and shear loading

Estimating Lower Limb Skeletal Loading 247

cylinder with an inner radius of 6 mm and an outer radius of 12 mm (corresponds to the dimensions of a rather typical tibia). Should the other end of this cylinder be rigidly fixed and either an axial force (compression) or a perpendicular force (bending moment) with the same magnitude be applied on the other end of the cylinder, the latter would result in 50x greater maximal stress (or strain) at the outer surface of the mid-shaft (furhter details on how stress may be estimated in section 2) (Rantalainen et al., 2010). A force during a high impact exercise, which is compressive for tibia will most likely cause bending in femoral neck and therefore it is plausible that a difference may be observed between femoral neck and tibial mid-shaft bone adaptation (or lack of) following the same high impact exercise intervention.

### **1.1 Estimating lower limb loading with reaction forces**

Several indirect methods have been used in estimating skeletal loading ranging from qualitative analyses to energy expenditure measurements. However, as per the title of the book, only mechanical loading estimates will be considered in the present text2. One relatively popular method to estimate loading, and thus strains, caused by locomotory actions on the bones is to examine the ground reaction forces or accelerations of the center of mass registered during these actions Jämsä et al. (2006); Turner & Robling (2003); Vainionpää (2007); Vainionpää et al. (2006). The advantage of measuring reaction forces or accelerations of the center of mass is the relative simplicity of the measurement. However, it is rarely brought into question if estimating skeletal loading from ground reaction forces is reasonable. Joint angles as well as muscle activity have a great influence on the loading of different bones and should be considered. In addition, diverse bone geometry and mineral content also have a significant influence on bone strains. An equal amount of force applied to different bodies will lead to different bone strains if the mineral amount and/or geometry of the loaded bones differ. The next evident step is thus to include the moment arms and joint angles into the calculations in inverse dynamics. This step, however, comes with an increase in the complexity of measurements as kinematics (i.e. movements of body parts in space) need to be recorded in addition to the aforementioned reaction forces. Furthermore the point of force application and direction of force needs to also be determined (detailed description of inverse dynamics is outside of the scope of the present chapter and detailed explanations can be found from other biomechanical textbooks, e.g. (Robertson, 2004)). However, even with inverse dynamics, the whole loading environment is not determined, since only external forces have been considered thus far. Muscles cause the movement of body parts in relation to each other and have relatively short moment arms compared to external forces. Therefore, the muscle (i.e. internal) forces are much higher than the measured reaction forces and neglecting these internal forces affects the loading estimate. This last step of estimating the loading environment is by far the most complex and includes arguably the biggest sources of error. Many muscles cross a given joint causing the same movement of the joint (e.g. soleus and gastrocnemii muscles plantar flex the ankle) making calculating the internal forces an overdetermined (redundant) problem, which means that there are infinite number of possible solutions of muscle forces3.

<sup>2</sup> And are thought of being the most meaningful estimates by the authors in any case.

<sup>3</sup> We have ignored at least one resourceful way of estimating skeletal loading. A walking simulator using cadaveric specimens has been built (Peterman et al., 2001), which produces reasonable strain results. However, the problem of redundant muscle forces applies also to that application.

### **2. Stress and strain**

4 Will-be-set-by-IN-TECH

cylinder with an inner radius of 6 mm and an outer radius of 12 mm (corresponds to the dimensions of a rather typical tibia). Should the other end of this cylinder be rigidly fixed and either an axial force (compression) or a perpendicular force (bending moment) with the same magnitude be applied on the other end of the cylinder, the latter would result in 50x greater maximal stress (or strain) at the outer surface of the mid-shaft (furhter details on how stress may be estimated in section 2) (Rantalainen et al., 2010). A force during a high impact exercise, which is compressive for tibia will most likely cause bending in femoral neck and therefore it is plausible that a difference may be observed between femoral neck and tibial mid-shaft bone

Several indirect methods have been used in estimating skeletal loading ranging from qualitative analyses to energy expenditure measurements. However, as per the title of the book, only mechanical loading estimates will be considered in the present text2. One relatively popular method to estimate loading, and thus strains, caused by locomotory actions on the bones is to examine the ground reaction forces or accelerations of the center of mass registered during these actions Jämsä et al. (2006); Turner & Robling (2003); Vainionpää (2007); Vainionpää et al. (2006). The advantage of measuring reaction forces or accelerations of the center of mass is the relative simplicity of the measurement. However, it is rarely brought into question if estimating skeletal loading from ground reaction forces is reasonable. Joint angles as well as muscle activity have a great influence on the loading of different bones and should be considered. In addition, diverse bone geometry and mineral content also have a significant influence on bone strains. An equal amount of force applied to different bodies will lead to different bone strains if the mineral amount and/or geometry of the loaded bones differ. The next evident step is thus to include the moment arms and joint angles into the calculations in inverse dynamics. This step, however, comes with an increase in the complexity of measurements as kinematics (i.e. movements of body parts in space) need to be recorded in addition to the aforementioned reaction forces. Furthermore the point of force application and direction of force needs to also be determined (detailed description of inverse dynamics is outside of the scope of the present chapter and detailed explanations can be found from other biomechanical textbooks, e.g. (Robertson, 2004)). However, even with inverse dynamics, the whole loading environment is not determined, since only external forces have been considered thus far. Muscles cause the movement of body parts in relation to each other and have relatively short moment arms compared to external forces. Therefore, the muscle (i.e. internal) forces are much higher than the measured reaction forces and neglecting these internal forces affects the loading estimate. This last step of estimating the loading environment is by far the most complex and includes arguably the biggest sources of error. Many muscles cross a given joint causing the same movement of the joint (e.g. soleus and gastrocnemii muscles plantar flex the ankle) making calculating the internal forces an overdetermined (redundant) problem, which means that there are infinite number of possible

adaptation (or lack of) following the same high impact exercise intervention.

<sup>2</sup> And are thought of being the most meaningful estimates by the authors in any case.

However, the problem of redundant muscle forces applies also to that application.

<sup>3</sup> We have ignored at least one resourceful way of estimating skeletal loading. A walking simulator using cadaveric specimens has been built (Peterman et al., 2001), which produces reasonable strain results.

**1.1 Estimating lower limb loading with reaction forces**

solutions of muscle forces3.

It needs to be mentioned that loading, forces, stresses and strains are used somewhat interchangeably in the present text, while generally, forces, stresses and strains should be differentiated. Therefore, the remainder of the paragraph tries to explain what are the differences between forces, stresses and strains and why they should not be mixed. The relationship between the applied load (=force) and relative deformation (Fig. 4) caused by the load can be plotted as a load-deformation curve. The curve may be divided into two parts, the linear elastic region and the non-linear plastic region. The stiffness of the structure is the slope of the elastic region, yield is the load at transition from the linear to elastic region and toughness is the area under curve up to yield or failure (Fig. 3) (Currey, 2001; Martin & Burr, 1989; Turner & Burr, 1993). The load-deformation curve may be converted to stress-strain curve by accounting for the geometry of the structure and the loading situation. The strain by definition is the deformation divided by the original dimension. If the loading force is compressive (or tensile), then the stress is simply force divided by the cross-sectional area of the structure. In bending and torsion, the stress varies throughout a cross-section. The slope of the stress-strain curve in turn is the Young's modulus, i.e. stiffness, of the material. The strength of the structure is defined as the load at which the structure either yields (yield strength) or breaks (breaking strength, which equals ultimate strength in bone). As was the case with stiffness, the strength may be reported as a material property (yield or breaking stress) or structural property (yield or breaking load). Because of the anisotropic structure of bone material, Young's modulus and breaking stress are different in longitudinal and transverse directions. Furthermore, bone ultimate stress is different in tension, compression and shear. While looking at different types of loading: bending, compression and shear, it becomes clear that the geometry of bone is of utmost importance in bending and shear loading (Currey, 2001; Turner & Burr, 1993).

Fig. 3. Schematic illustration of a load-deformation curve.

### **2.1 Example of how geometry affects bone strength**

Bending stress in a bone cross-section in a given loading configuration is determined by the second moment of area (also known as cross-sectional moment of inertia) of the cross-section as shown in (Fig. 1 A). Stress (indicated by the shaded area) is distributed unevenly throughout the cross-section, being by definition zero at the neutral axis and increasing its magnitude linearly while moving away from the neutral axis. The stress (*σ*) at a given point

Fig. 4. Strain describes the relative deformation of a structure. If the rectangle (representing a solid structure) on the left is subjected to a tension force as indicated by the arrows on the narrower rectangle (representing the same structure under tension), the rectangle will deform. If the structure is incompressible, the width will change roughly half of the relative change in height but to opposite direction (remember that any given real life structure is three dimensional and volume is maintained if the structure is incompressible, i.e. Poisson's ratio = 0.5), i.e. the rectangle becomes longer but narrower. If the force was compressive, the rectangle would become shorter but thicker. In addition to being subjected to both tensile and compressive strain at the same point, due to being incompressible, the opposite forms of

Estimating Lower Limb Skeletal Loading 249

strain also cause shear strain at the same point. Unlike the example given, bone is

strains.

compressible to some extent (cortical bone Poisson's ratio is somewhere between 0.2 to 0.4). Nevertheless, any given loaded point in bone will nevertheless always develop all forms of

Furthermore, if one were to estimate the magnitude of bone stresses and strains, just calculating based on external forces, (i.e. joint moments), and neglect the internal forces (i.e. muscle forces), the stresses and strains could potentially be greatly miss-characterized (an example is given under section 4). Again, using tibial mid-shaft as an example consider standing on the balls of the feet. If the fact that the moments are caused by muscles were neglected, one would include only ankle joint moment and the weight of the body above tibial mid-shaft in calculating the tibial mid-shaft stresses and strains. However, since the ankle joint moment is in fact produced by the ankle extensors, the compressive and bending forces cause by the muscle activity should be taken into consideration instead of joint moments. Due to the ratio of contact point to ankle joint moment arm and ankle extensor muscles' insertion to ankle joint moment arm, the force required from the ankle extensor musculature is in the order of two to three times body weight (Finni et al., 1998; Ishikawa et al., 2005; Komi et al., 1992). Furthermore, bending moments and compressive forces load the bone in a very different manner and neglecting the muscles as the cause of the joint moment has the potential of greatly influencing the result (Fig. 6). Knowing bone stresses and strains in a single point is useful only to a limited extent and for most practical applications, the stresses and strains need to be estimated across the bone cross-section and at several points through the length of the

within the cross-section is determined by multiplying the bending moment (M) with the distance coordinate (y) from the neutral axis (dashed line) and dividing the product with the second moment of area (I) of the cross-section in relation to the neutral axis (i.e. *σ* = *<sup>M</sup> <sup>I</sup> y*). This leads to the conclusion that the bending stress will be positive (tensile) on one side of the neutral axis and negative (compressive) on the other. In case of symmetrical cross-sections the stress magnitudes on both extremes will be equal. In case of a real bone it is rarely the case. Let's simplify the bone mid-shaft cross-section to a hollow oval with external radii *ae*, *be* and internal *ai*, *bi* as depicted in (Fig. 1 B). The area of an oval is equal to *πab*, where a and b are the radii of the oval and thus the area of an oval ring is given by *π*(*aebe* − *aibi*). If the area is increased by 2% in a way observed in association to exercise interventions, only *be* will be increased as shown in figure (Fig. 1 C). Therefore the following relation holds:

$$(1.02 \cdot \pi (a\_{\varepsilon} b\_{\varepsilon} - a\_{i} b\_{i}) = \pi (a\_{\varepsilon} b\_{\varepsilon - \text{new}} - a\_{i} b\_{i}) \tag{2}$$

and can be transformed to compute the corresponding new external radius, *be*−*new*, with the means of following equation:

$$b\_{\varepsilon-n\varepsilon w} = \frac{1.02 \cdot (a\_{\varepsilon}b\_{\varepsilon} - a\_{l}b\_{l}) + a\_{l}b\_{l}}{a\_{\varepsilon}}.\tag{3}$$

The cross-sectional second moment of area of a solid ellipse,*I*, is equal to:

$$I = \frac{\pi}{4}ab^3\tag{4}$$

and thus for an elliptic ring, the second moment of area will be:

$$I = \frac{\pi}{4} (a\_{\varepsilon} b\_{\varepsilon}^3 - a\_i b\_i^3). \tag{5}$$

Lets assume *ai* = 0.5, *bi* = 0.5 and *ae* = 1.0, *be* = 1.5, which if the units were in cm, corresponds to young adult men tibial mid-shaft. In this particular case the cross-sectional inertia prior to the intervention would have been 2.60*cm*<sup>4</sup> and after 2.74*cm*4, which is a 5% increase. Accordingly the change would have lead to a 3% decrease in maximal stress. For further details on bone mechanics, please refer to pertinent literature such as prof. Currey's book (Currey, 2002).

### **3. Reaction forces and inverse dynamics in estimating skeletal loading**

Lower limb skeletal loading (bone stresses) may be estimated with a combination of kinetic and kinematic measurements (Anderson et al., 1996). Kinetic measurements are typically the measurement of vertical and horizontal ground reaction forces e.g. during locomotion, whereas the kinematic measurements (displacements, velocities and accelerations of segments and respective angular equivalents of joints) are derived from the motion captured from the subject. The measured kinematics is combined with the measured kinetics in inverse dynamics to calculate joint moments. Judging by the apparent similarity of the patterns of ankle moment and tibial mid-shaft maximal principal strain in walking (Fig. 5), simply looking at joint moments appears to give a reasonable estimate of tibial mid-shaft bone strains. However, it needs to be remembered that tibia is also affected by the knee joint moment and the patterns of the knee joint moments differ quite markedly from the patterns observed at the ankle joint.

6 Will-be-set-by-IN-TECH

within the cross-section is determined by multiplying the bending moment (M) with the distance coordinate (y) from the neutral axis (dashed line) and dividing the product with the second moment of area (I) of the cross-section in relation to the neutral axis (i.e. *σ* = *<sup>M</sup>*

This leads to the conclusion that the bending stress will be positive (tensile) on one side of the neutral axis and negative (compressive) on the other. In case of symmetrical cross-sections the stress magnitudes on both extremes will be equal. In case of a real bone it is rarely the case. Let's simplify the bone mid-shaft cross-section to a hollow oval with external radii *ae*, *be* and internal *ai*, *bi* as depicted in (Fig. 1 B). The area of an oval is equal to *πab*, where a and b are the radii of the oval and thus the area of an oval ring is given by *π*(*aebe* − *aibi*). If the area is increased by 2% in a way observed in association to exercise interventions, only *be* will be

and can be transformed to compute the corresponding new external radius, *be*−*new*, with the

*be*−*new* <sup>=</sup> 1.02 · (*aebe* <sup>−</sup> *aibi*) + *aibi*

*<sup>I</sup>* <sup>=</sup> *<sup>π</sup>*

<sup>4</sup> (*aeb*<sup>3</sup>

Lets assume *ai* = 0.5, *bi* = 0.5 and *ae* = 1.0, *be* = 1.5, which if the units were in cm, corresponds to young adult men tibial mid-shaft. In this particular case the cross-sectional inertia prior to the intervention would have been 2.60*cm*<sup>4</sup> and after 2.74*cm*4, which is a 5% increase. Accordingly the change would have lead to a 3% decrease in maximal stress. For further details on bone mechanics, please refer to pertinent literature such as prof. Currey's

Lower limb skeletal loading (bone stresses) may be estimated with a combination of kinetic and kinematic measurements (Anderson et al., 1996). Kinetic measurements are typically the measurement of vertical and horizontal ground reaction forces e.g. during locomotion, whereas the kinematic measurements (displacements, velocities and accelerations of segments and respective angular equivalents of joints) are derived from the motion captured from the subject. The measured kinematics is combined with the measured kinetics in inverse dynamics to calculate joint moments. Judging by the apparent similarity of the patterns of ankle moment and tibial mid-shaft maximal principal strain in walking (Fig. 5), simply looking at joint moments appears to give a reasonable estimate of tibial mid-shaft bone strains. However, it needs to be remembered that tibia is also affected by the knee joint moment and the patterns of the knee joint moments differ quite markedly from the patterns observed at

*ae*

*<sup>e</sup>* <sup>−</sup> *aib*<sup>3</sup>

1.02 · *π*(*aebe* − *aibi*) = *π*(*aebe*−*new* − *aibi*) (2)

increased as shown in figure (Fig. 1 C). Therefore the following relation holds:

The cross-sectional second moment of area of a solid ellipse,*I*, is equal to:

*<sup>I</sup>* <sup>=</sup> *<sup>π</sup>*

**3. Reaction forces and inverse dynamics in estimating skeletal loading**

and thus for an elliptic ring, the second moment of area will be:

means of following equation:

book (Currey, 2002).

the ankle joint.

*<sup>I</sup> y*).

. (3)

<sup>4</sup> *ab*<sup>3</sup> (4)

*<sup>i</sup>* ). (5)

Fig. 4. Strain describes the relative deformation of a structure. If the rectangle (representing a solid structure) on the left is subjected to a tension force as indicated by the arrows on the narrower rectangle (representing the same structure under tension), the rectangle will deform. If the structure is incompressible, the width will change roughly half of the relative change in height but to opposite direction (remember that any given real life structure is three dimensional and volume is maintained if the structure is incompressible, i.e. Poisson's ratio = 0.5), i.e. the rectangle becomes longer but narrower. If the force was compressive, the rectangle would become shorter but thicker. In addition to being subjected to both tensile and compressive strain at the same point, due to being incompressible, the opposite forms of strain also cause shear strain at the same point. Unlike the example given, bone is compressible to some extent (cortical bone Poisson's ratio is somewhere between 0.2 to 0.4). Nevertheless, any given loaded point in bone will nevertheless always develop all forms of strains.

Furthermore, if one were to estimate the magnitude of bone stresses and strains, just calculating based on external forces, (i.e. joint moments), and neglect the internal forces (i.e. muscle forces), the stresses and strains could potentially be greatly miss-characterized (an example is given under section 4). Again, using tibial mid-shaft as an example consider standing on the balls of the feet. If the fact that the moments are caused by muscles were neglected, one would include only ankle joint moment and the weight of the body above tibial mid-shaft in calculating the tibial mid-shaft stresses and strains. However, since the ankle joint moment is in fact produced by the ankle extensors, the compressive and bending forces cause by the muscle activity should be taken into consideration instead of joint moments. Due to the ratio of contact point to ankle joint moment arm and ankle extensor muscles' insertion to ankle joint moment arm, the force required from the ankle extensor musculature is in the order of two to three times body weight (Finni et al., 1998; Ishikawa et al., 2005; Komi et al., 1992). Furthermore, bending moments and compressive forces load the bone in a very different manner and neglecting the muscles as the cause of the joint moment has the potential of greatly influencing the result (Fig. 6). Knowing bone stresses and strains in a single point is useful only to a limited extent and for most practical applications, the stresses and strains need to be estimated across the bone cross-section and at several points through the length of the

Fig. 6. Freebody diagrams for estimating tibial mid-shaft bone loading while standing on one leg on the ball of the foot. A) Forces accounted for when estimating tibial loading with just the ground reaction force, i.e. ground reaction force (GRF) and its counter force (weight of the subject). B) Estimating tibial loading after inverse dynamics. Compressive loading is accounted for by the GRF and bending is accounted for by the ankle joint moment. It is quite easy to imagine a situation, in which tibial loading is not proportional to ground reaction forces. Consider the situation in pane A and then imagine one standing with the ankle fully plantar flexed. In both cases the GRFs are the same, while ankle joint moment is changed (if one were a ballerina, the ankle joint moment could actually be made to drop to zero). C) Tibial loading is estimated from joint contact forces comprising reaction forces and muscular

Estimating Lower Limb Skeletal Loading 251

Lumped mass formulation is another method that can be used to describe mechanical flexibility (Hughes, 1979; Hughes & Winget, 1980; Huston, 1981; 1991; Kamman & Huston, 2001; Raman-Nair & Baddour, 2003). In this formulation, a flexible body is replaced by a set of point masses connected via springs. Using a sufficient amount of springs and masses allows for reasonably accurate mass distribution, within inhomogeneous bodies as well. Similar to the linear theory of elastodynamics, the lumped mass approach does not require that any changes be made in the standard rigid multibody solver. In contrast to the rigid body representation, the performance decays with the rise of discretization precision. Therefore,

The third major method for introducing flexibility into the multibody formulation is the floating frame of reference (Shabana, 2005). Formulation relies on coordinate partitioning so that one set of coordinates is used to describe the flexible body's reference frame in the global coordinate system and another set of coordinates is used to describe the deformation of the body in the local frame of reference. Originally, the deformation of the body was described in a similar fashion to the finite element method, resulting in a remarkable increase in regards to computational effort. However, Shabana and Wehage (Shabana, 1983) have developed a solution to this problem by replacing the full finite element models of flexible structures with deformation modes description. This allowed for a reduction in the amount of deformation

this method can only be practically applied to beam structures.

forces.

Fig. 5. Schematic illustration of tibial mid-shaft anteromedial aspect maximal principal strain (Lanyon et al., 1975) (upper pane) and ankle sagittal plane joint moment (Silder et al., 2008) during a walking cycle.

bone. Therefore, it quickly becomes evident that building models for individual applications, solving the muscle forces (which are redundant and therefore need to be estimated with some kind of optimization criteria) and calculating the stresses and strains is not feasible and a more general modelling approach needs to be chosen. Flexible multibody dynamics methods can be used as such a generalized modelling method.

### **3.1 Flexible multibody dynamics**

Flexibility can be accounted for in a number of ways in multibody applications. The linear theory of elastodynamics can be considered a traditional approach to account for flexibility. This approach relies on the assumption of small deformations in the flexible bodies. Thus, rigid body simulation is decoupled with the deformation computation. The rigid body simulation is performed to obtain the external, as well as internal, forces acting upon each of the bodies. These forces are later imposed on the finite element model of the body for which deformations, stresses and strains can be obtained (Erdman & Sandor, 1972; Lowen & Chassapis, 1986; Lowen & Jandrasits, 1972; Turcic & Midha, 1984; Turcic et al., 1984; Winfrey, 1971). Standard multibody and finite element solvers can be used for this method, which is a great benefit. Additionally, a considerable increase in computational speed can be achieved if the force application points are known a priori and do not vary over the duration of the simulation. In such cases, the use of linear finite element analysis is reduced to a single computation of the full model. Strains, stresses and deformations can then be computed for each time step as post-processing. Conversely, the flexibility of the bodies does not affect the multibody simulation behaviour, which is the main disadvantage, especially in case of considerable deformations.

8 Will-be-set-by-IN-TECH

Fig. 5. Schematic illustration of tibial mid-shaft anteromedial aspect maximal principal strain (Lanyon et al., 1975) (upper pane) and ankle sagittal plane joint moment (Silder et al., 2008)

bone. Therefore, it quickly becomes evident that building models for individual applications, solving the muscle forces (which are redundant and therefore need to be estimated with some kind of optimization criteria) and calculating the stresses and strains is not feasible and a more general modelling approach needs to be chosen. Flexible multibody dynamics methods can

Flexibility can be accounted for in a number of ways in multibody applications. The linear theory of elastodynamics can be considered a traditional approach to account for flexibility. This approach relies on the assumption of small deformations in the flexible bodies. Thus, rigid body simulation is decoupled with the deformation computation. The rigid body simulation is performed to obtain the external, as well as internal, forces acting upon each of the bodies. These forces are later imposed on the finite element model of the body for which deformations, stresses and strains can be obtained (Erdman & Sandor, 1972; Lowen & Chassapis, 1986; Lowen & Jandrasits, 1972; Turcic & Midha, 1984; Turcic et al., 1984; Winfrey, 1971). Standard multibody and finite element solvers can be used for this method, which is a great benefit. Additionally, a considerable increase in computational speed can be achieved if the force application points are known a priori and do not vary over the duration of the simulation. In such cases, the use of linear finite element analysis is reduced to a single computation of the full model. Strains, stresses and deformations can then be computed for each time step as post-processing. Conversely, the flexibility of the bodies does not affect the multibody simulation behaviour, which is the main disadvantage, especially in case of

during a walking cycle.

**3.1 Flexible multibody dynamics**

considerable deformations.

be used as such a generalized modelling method.

Fig. 6. Freebody diagrams for estimating tibial mid-shaft bone loading while standing on one leg on the ball of the foot. A) Forces accounted for when estimating tibial loading with just the ground reaction force, i.e. ground reaction force (GRF) and its counter force (weight of the subject). B) Estimating tibial loading after inverse dynamics. Compressive loading is accounted for by the GRF and bending is accounted for by the ankle joint moment. It is quite easy to imagine a situation, in which tibial loading is not proportional to ground reaction forces. Consider the situation in pane A and then imagine one standing with the ankle fully plantar flexed. In both cases the GRFs are the same, while ankle joint moment is changed (if one were a ballerina, the ankle joint moment could actually be made to drop to zero). C) Tibial loading is estimated from joint contact forces comprising reaction forces and muscular forces.

Lumped mass formulation is another method that can be used to describe mechanical flexibility (Hughes, 1979; Hughes & Winget, 1980; Huston, 1981; 1991; Kamman & Huston, 2001; Raman-Nair & Baddour, 2003). In this formulation, a flexible body is replaced by a set of point masses connected via springs. Using a sufficient amount of springs and masses allows for reasonably accurate mass distribution, within inhomogeneous bodies as well. Similar to the linear theory of elastodynamics, the lumped mass approach does not require that any changes be made in the standard rigid multibody solver. In contrast to the rigid body representation, the performance decays with the rise of discretization precision. Therefore, this method can only be practically applied to beam structures.

The third major method for introducing flexibility into the multibody formulation is the floating frame of reference (Shabana, 2005). Formulation relies on coordinate partitioning so that one set of coordinates is used to describe the flexible body's reference frame in the global coordinate system and another set of coordinates is used to describe the deformation of the body in the local frame of reference. Originally, the deformation of the body was described in a similar fashion to the finite element method, resulting in a remarkable increase in regards to computational effort. However, Shabana and Wehage (Shabana, 1983) have developed a solution to this problem by replacing the full finite element models of flexible structures with deformation modes description. This allowed for a reduction in the amount of deformation

Fig. 7. Model diagram with tibia cross-section at which the stresses and strains will be

*Fg* <sup>=</sup> *<sup>m</sup>*

Estimating Lower Limb Skeletal Loading 253

forces, *GRFf* , as well as vertical component of the ground reaction force, *GRFn* can be solved

Σ*MiD* = 0 ⇔ *Fg* · (*Ax* − *Dx*) − *Fs* · (*Ay* − *Dy*) = 0

where, *Ax*, *Ay*, *Dx*, *Dy* are x and y coordinates of hip joint and foot-ground contact point

Knowing the value of *Fg* from equation 10, value of *GRFn* can be obtained from equation 8, yielding *GRFn* = 392.43*N*. Substituting numerical values to equation 9, *Fs* can be solved as:

From equation 7 value of *GRFf* can be expressed as: *GRFf* = 19.28*N*. Knowing all the external forces applied to the model in the equilibrium, internal forces and moments can be computed.

*Ay* − *Dy*

<sup>2</sup> · 9.81066 *<sup>N</sup>*

<sup>=</sup> 392.43*<sup>N</sup>* · (−35.33*mm* <sup>−</sup> (−81.47*mm*))

<sup>⇔</sup> *Fs* <sup>=</sup> *Fg* · (*Ax* <sup>−</sup> *Dx*)

correspondingly. Substituting to equation 6 numerical values *Fg* can be obtained:

<sup>2</sup> · *<sup>g</sup>* <sup>=</sup> <sup>80</sup>*kg*

*Fg* <sup>=</sup> *<sup>m</sup>*

*Fs* <sup>=</sup> *Fg* · (*Ax* <sup>−</sup> *Dx*) *Ay* − *Dy*

<sup>2</sup> · *<sup>g</sup>* (6)

*kg* . The unknown stabilization, *Fs*, and friction

*kg* <sup>=</sup> 392.43*<sup>N</sup>* (10)

<sup>−</sup>24.84*mm* <sup>−</sup> (−964.11*mm*) <sup>=</sup> 19.28*<sup>N</sup>* (11)

(9)

Σ*Fix* = 0 ⇔ *GRFf* − *Fs* = 0 ⇔ *GRFf* = *Fs* (7)

Σ*Fiy* = 0 ⇔ *GRFn* − *Fg* = 0 ⇔ *GRFn* = *Fg* (8)

computed with the use of ground reaction forces only

where *g* is the gravity constant equal to 9.81066 *<sup>N</sup>*

from the equilibrium equations:

coordinates (in the range of thousands) to a reasonable amount of modal coordinates, making this method an effective compromise between accuracy and computational effort.

### **3.2 Modelling results obtained with flexible multibody dynamics in human locomotion**

There is a relative paucity of dynamic bone loading modelling results in human locomotion in the present literature. One notable exception is the results from Anderson et al. (1996), who modelled femoral neck maximal stresses in several gym exercises. However, Anderson et al. (1996) reported bone stresses from only one location (i.e. the femoral neck). The floating frame of reference method with modal coordinates has been successfully applied in modelling tibial and femoral bone strains during walking (Klodowski et al., 2011; Nazer et al., 2011; Nazer, Klodowski, Rantalainen, Heinonen, Sievänen & Mikkola, 2008; Nazer, Rantalainen, Heinonen, Sievänen & Mikkola, 2008). When looking at femoral neck strain estimates and tibial mid-shaft strain estimates it may be seen that the timing of peak strains differ between these two bone sites. Remembering from the introduction, that literature reports femoral neck adaptation while no adaptation is seen at tibial mid-shaft following an intervention, it certainly appears that modelling based approaches of estimating skeletal loading may prove useful in the future.

### **4. Numerical example**

To illustrate the differences between presented modelling techniques a simple numerical example will be presented. In the example axial stresses and strains at the mid-shaft cross-section of the tibia of a subject standing on the balls of his feet will be estimated. To make the example more realistic, the dimensions of the subject were taken from the LifeMod4 database and correspond to a male whose height is 180cm, weight 80kg and age 24 years. In the example the subject extends his ankle joints by 20◦, and tilts his body forward by 10◦. In this condition the center of mass will be slightly ahead of the foot contact, which will add extra torque at joints to make the example more interesting. Only a right half of a lower body will be considered, meaning: pelvis, right upper leg, right lower leg and a right foot. To make simulation models comparable with hand calculations, which are based solely on ground reaction forces, location of the cross-section at which stresses and strains are computed is preserved together with it's geometry described by a polygon. All models share also the material model for the bone which for simplicity is isotropic and homogeneous, with Young's modulus of 17GPa (Dong & Guo, 2004) and Poisson's ratio equal to 0.3 (Reilly & Burstein, 1975).

### **4.1 Approach based on ground reaction force (model 1)**

Let's start the analysis with the simplest 2D (saggital plane) model, which would be the most intuitive first choice (Fig. 7).

In the model the mass is concentrated in a point at the hip joint, represented by body gravity force, *Fg*. This location is close to the real center of mass for the human. At the same location stabilization force *Fs* is applied to constrain the model from falling. At the foot-ground contact point, ground reaction force is applied, for simplicity divided into normal, *GRFn*, and tangential components, *GRFf* which correspond to static friction. Knowing the mass of the subject, *m* = 80*kg*, we can assume for symmetric case, that the load carried by each leg is equal, thus:

<sup>4</sup> LifeMod plug-in for MSC Adams, version 2008.0.0, LifeModdeler, Inc., San Clemente, California, USA

10 Will-be-set-by-IN-TECH

coordinates (in the range of thousands) to a reasonable amount of modal coordinates, making

To illustrate the differences between presented modelling techniques a simple numerical example will be presented. In the example axial stresses and strains at the mid-shaft cross-section of the tibia of a subject standing on the balls of his feet will be estimated. To make the example more realistic, the dimensions of the subject were taken from the LifeMod4 database and correspond to a male whose height is 180cm, weight 80kg and age 24 years. In the example the subject extends his ankle joints by 20◦, and tilts his body forward by 10◦. In this condition the center of mass will be slightly ahead of the foot contact, which will add extra torque at joints to make the example more interesting. Only a right half of a lower body will be considered, meaning: pelvis, right upper leg, right lower leg and a right foot. To make simulation models comparable with hand calculations, which are based solely on ground reaction forces, location of the cross-section at which stresses and strains are computed is preserved together with it's geometry described by a polygon. All models share also the material model for the bone which for simplicity is isotropic and homogeneous, with Young's modulus of 17GPa (Dong & Guo, 2004) and Poisson's ratio equal to 0.3 (Reilly & Burstein,

Let's start the analysis with the simplest 2D (saggital plane) model, which would be the most

In the model the mass is concentrated in a point at the hip joint, represented by body gravity force, *Fg*. This location is close to the real center of mass for the human. At the same location stabilization force *Fs* is applied to constrain the model from falling. At the foot-ground contact point, ground reaction force is applied, for simplicity divided into normal, *GRFn*, and tangential components, *GRFf* which correspond to static friction. Knowing the mass of the subject, *m* = 80*kg*, we can assume for symmetric case, that the load carried by each leg is

<sup>4</sup> LifeMod plug-in for MSC Adams, version 2008.0.0, LifeModdeler, Inc., San Clemente, California, USA

this method an effective compromise between accuracy and computational effort.

**4. Numerical example**

intuitive first choice (Fig. 7).

**4.1 Approach based on ground reaction force (model 1)**

1975).

equal, thus:

**3.2 Modelling results obtained with flexible multibody dynamics in human locomotion** There is a relative paucity of dynamic bone loading modelling results in human locomotion in the present literature. One notable exception is the results from Anderson et al. (1996), who modelled femoral neck maximal stresses in several gym exercises. However, Anderson et al. (1996) reported bone stresses from only one location (i.e. the femoral neck). The floating frame of reference method with modal coordinates has been successfully applied in modelling tibial and femoral bone strains during walking (Klodowski et al., 2011; Nazer et al., 2011; Nazer, Klodowski, Rantalainen, Heinonen, Sievänen & Mikkola, 2008; Nazer, Rantalainen, Heinonen, Sievänen & Mikkola, 2008). When looking at femoral neck strain estimates and tibial mid-shaft strain estimates it may be seen that the timing of peak strains differ between these two bone sites. Remembering from the introduction, that literature reports femoral neck adaptation while no adaptation is seen at tibial mid-shaft following an intervention, it certainly appears that modelling based approaches of estimating skeletal loading may prove useful in the future.

Fig. 7. Model diagram with tibia cross-section at which the stresses and strains will be computed with the use of ground reaction forces only

$$F\_{\mathcal{S}} = \frac{m}{2} \cdot g \tag{6}$$

where *g* is the gravity constant equal to 9.81066 *<sup>N</sup> kg* . The unknown stabilization, *Fs*, and friction forces, *GRFf* , as well as vertical component of the ground reaction force, *GRFn* can be solved from the equilibrium equations:

$$
\Sigma F\_{\rm lx} = 0 \Leftrightarrow GRF\_f - F\_{\rm s} = 0 \Leftrightarrow GRF\_f = F\_{\rm s} \tag{7}
$$

$$
\Sigma F\_{\rm liy} = 0 \Leftrightarrow GRF\_{\rm li} - F\_{\mathcal{S}} = 0 \Leftrightarrow GRF\_{\rm li} = F\_{\mathcal{S}} \tag{8}
$$

$$\begin{aligned} \Sigma M\_{\text{iD}} = 0 &\Leftrightarrow F\_{\text{\(\)}} \cdot (A\_{\text{\(\)}} - D\_{\text{\(\)}}) - F\_{\text{\(\)}} \cdot (A\_{\text{\(\)}} - D\_{\text{\(\)}}) = 0 \\ &\Leftrightarrow F\_{\text{\(\)}} = \frac{F\_{\text{\(\)}} \cdot (A\_{\text{\(\)}} - D\_{\text{\(\)}})}{A\_{\text{\(\)}} - D\_{\text{\(\)}}} \end{aligned} \tag{9}$$

$$\begin{array}{c} \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \quad \cdot \\ \text{where,} \quad A\_{\mathcal{X}}, A\_{\mathcal{Y}}, D\_{\mathcal{X}}, D\_{\mathcal{Y}} \text{ are x and y coordinates of hip joint and foot-ground contact points,} \\ \text{corresponding. Substituting to equation 6 numerical values } F\_{\mathcal{S}} \text{ can be obtained:} \end{array}$$

$$F\_{\mathcal{S}} = \frac{m}{2} \cdot \text{g} = \frac{80 \text{kg}}{2} \cdot 9.81066 \frac{\text{N}}{\text{kg}} = 392.43 \text{N} \tag{10}$$

Knowing the value of *Fg* from equation 10, value of *GRFn* can be obtained from equation 8, yielding *GRFn* = 392.43*N*. Substituting numerical values to equation 9, *Fs* can be solved as:

$$F\_{\rm S} = \frac{F\_{\rm g} \cdot (A\_{\rm x} - D\_{\rm x})}{A\_{\rm y} - D\_{\rm y}} = \frac{392.43N \cdot (-35.33mm - (-81.47mm))}{-24.84mm - (-964.11mm)} = 19.28N\tag{11}$$

From equation 7 value of *GRFf* can be expressed as: *GRFf* = 19.28*N*. Knowing all the external forces applied to the model in the equilibrium, internal forces and moments can be computed.

**Point x [mm] ¯ z [mm] ¯** 13.11 -4.1 10.29 -5.59 7.12 -9.23 0.67 -10.06 -4.84 -10.35 -12.45 -3.75 -13.61 3.68 -8.18 13.6 5.76 5.26 9.31 1.04

Estimating Lower Limb Skeletal Loading 255

Table 1. Coordinates of the vertices of the polygon defining the cross-section of the bone

where *n* is the number of points, point *n* + 1 refers to point 1 as the polygon must be closed. Note that the equation 15 holds when the points are ordered in clockwise direction, for the anticlockwise order of points the minus sign has to be added. After computations the area of

point P can be obtained from the geometrical dependencies depicted in figure 8B. Note that

*<sup>x</sup>* is a coordinate of the point P with respect to point P' not a distance (absolute value of a



*Px* = −159.55*mm* + 225.76*mm* · sin(10◦) − 1.84*mm* · cos(10◦) = −122.16*mm* (19)

*Py* = −840.22*mm* + 225.76*mm* · cos(10◦) − (−1.84*mm*) · sin(10◦) = −617.57*mm* (20)

2

*<sup>x</sup>*˜ <sup>=</sup> *<sup>x</sup>*¯ <sup>−</sup> *<sup>P</sup>*¯

Substituting numerical values to equations 21 and 22 yields: *Iz*¯ = 15008*mm*4. Finally, having all cross-sectional parameters, and most importantly coordinates of point P in global

(*x*¯*i*+<sup>1</sup> · *z*¯*<sup>i</sup>* − *x*¯*<sup>i</sup>* · *z*¯*i*+1) (15)

(*x*¯*<sup>i</sup>* + *x*¯*i*+1) · (*x*¯*i*+<sup>1</sup> · *z*¯*<sup>i</sup>* − *x*¯*<sup>i</sup>* · *z*¯*i*+1) (16)

*<sup>x</sup>*, can be obtained from:

*<sup>x</sup>* = −1.84*mm*. Global coordinates of

*<sup>x</sup>* · cos(*α*) (17)

*<sup>x</sup>* · sin(*α*) (18)

*<sup>i</sup>*+1) · (*x*˜*i*+<sup>1</sup> · *z*¯*<sup>i</sup>* − *x*˜*<sup>i</sup>* · *z*¯*i*+1); (21)

*<sup>x</sup>* (22)

Area of a polygon can be obtained using the following equation:

*<sup>A</sup>* <sup>=</sup> <sup>1</sup> 2 *n* ∑ *i*=1

the bone cross-section is *A* = 393.07*mm*2. Centroid y coordinate, *P*¯

*n* ∑ *i*=1

*Px* = *Cx* + |*CP*�

*Py* = *Cy* + |*CP*�

Second moment of area of the bone cross-section about *z*¯ axis can be defined as:

*<sup>i</sup>* + *x*˜*<sup>i</sup>* · *x*˜*i*+<sup>1</sup> + *x*˜

*P*¯ *<sup>x</sup>* <sup>=</sup> <sup>1</sup> 6*A*

*P*¯

where:

Substituting numerical values to equation 16 yields *P*¯

coordinate), thus sign convention in the equations 17 and 18.

Substituting numerical values to equations 17, 18 yields:

*n* ∑ *i*=1 (*x*˜ 2

*Iz*¯ <sup>=</sup> <sup>1</sup> 12

Fig. 8. A) Freebody diagram of the foot and lower part of tibia. B) Lower part of tibia diagram used for calculations of global coordinates of point P

The cross-section at which the stresses and strains will be computed is located at point P, as marked in figure 8. Solving of internal normal force, *Fn*, transverse force, *Ft*, and bending moment, *M* will be based on free body diagram depicted in figure 8.

Equilibrium equations for the bottom part of tibia and a foot can be expressed in the cross-section coordinate system *x*¯ − *y*¯ for simplicity:

$$\begin{aligned} \Sigma F\_{\text{if}} &= 0 \Leftrightarrow F\_{\text{f}} + \mathcal{G}RF\_{\text{f}} \cdot \cos(\alpha) - \mathcal{G}RF\_{\text{h}} \cdot \sin(\alpha) = 0 \\ &\Leftrightarrow F\_{\text{f}} = \mathcal{G}RF\_{\text{h}} \cdot \sin(\alpha) - \mathcal{G}RF\_{\text{f}} \cdot \cos(\alpha) \end{aligned} \tag{12}$$

$$\begin{aligned} \Sigma F\_{i\mathcal{G}} &= 0 \Leftrightarrow -F\_{\mathcal{U}} + GRF\_f \cdot \sin(a) + GRF\_{\mathcal{U}} \cdot \cos(a) = 0 \\ &\Leftrightarrow F\_{\mathcal{U}} = GRF\_f \cdot \sin(a) + GRF\_{\mathcal{U}} \cdot \cos(a) \end{aligned} \tag{13}$$

$$\begin{aligned} \Sigma M\_{lP} = 0 &\Leftrightarrow M - GRF\_f \cdot (P\_y - D\_y) - GRF\_{ll} \cdot (P\_x - D\_x) = 0\\ &\Leftrightarrow M = GRF\_f \cdot (P\_y - D\_y) + GRF\_{ll} \cdot (P\_x - D\_x) \end{aligned} \tag{14}$$

Note that coordinates of points P and D are expressed in the global coordinate system due to the fact, that forces *GRFf* and *GRFn* act along global coordinate system axes. The resulting moments from those forces around point P are then computed as the force multiplied the moment arm, which is the distance between contact point and point P and is measured along axis perpendicular to the force direction. In this example we assume that moment has a positive sign if it results in clockwise rotation. In order to solve the equilibrium equations location of point P is required. Point P lies at the neutral axis of the bone cross-section. Thus at this point we will compute the cross-sectional parameters: second moment of area along bending axis, centroid location, and area. Bone cross-section is taken from the generic bone model generated by LifeMod. It is described as a polygon (Fig. 7) with vertices coordinates given in table 1.

12 Will-be-set-by-IN-TECH

Fig. 8. A) Freebody diagram of the foot and lower part of tibia. B) Lower part of tibia

The cross-section at which the stresses and strains will be computed is located at point P, as marked in figure 8. Solving of internal normal force, *Fn*, transverse force, *Ft*, and bending

Equilibrium equations for the bottom part of tibia and a foot can be expressed in the

Σ*Fix*¯ = 0 ⇔ *Ft* + *GRFf* · cos(*α*) − *GRFn* · sin(*α*) = 0

Σ*Fiy*¯ = 0 ⇔ −*Fn* + *GRFf* · sin(*α*) + *GRFn* · cos(*α*) = 0

Σ*MiP* = 0 ⇔ *M* − *GRFf* · (*Py* − *Dy*) − *GRFn* · (*Px* − *Dx*) = 0

Note that coordinates of points P and D are expressed in the global coordinate system due to the fact, that forces *GRFf* and *GRFn* act along global coordinate system axes. The resulting moments from those forces around point P are then computed as the force multiplied the moment arm, which is the distance between contact point and point P and is measured along axis perpendicular to the force direction. In this example we assume that moment has a positive sign if it results in clockwise rotation. In order to solve the equilibrium equations location of point P is required. Point P lies at the neutral axis of the bone cross-section. Thus at this point we will compute the cross-sectional parameters: second moment of area along bending axis, centroid location, and area. Bone cross-section is taken from the generic bone model generated by LifeMod. It is described as a polygon (Fig. 7) with vertices coordinates

<sup>⇔</sup> *Ft* <sup>=</sup> *GRFn* · sin(*α*) <sup>−</sup> *GRFf* · cos(*α*) (12)

<sup>⇔</sup> *Fn* <sup>=</sup> *GRFf* · sin(*α*) + *GRFn* · cos(*α*) (13)

<sup>⇔</sup> *<sup>M</sup>* <sup>=</sup> *GRFf* · (*Py* <sup>−</sup> *Dy*) + *GRFn* · (*Px* <sup>−</sup> *Dx*) (14)

diagram used for calculations of global coordinates of point P

cross-section coordinate system *x*¯ − *y*¯ for simplicity:

given in table 1.

moment, *M* will be based on free body diagram depicted in figure 8.


Table 1. Coordinates of the vertices of the polygon defining the cross-section of the bone

Area of a polygon can be obtained using the following equation:

$$A = \frac{1}{2} \sum\_{i=1}^{n} (\bar{\mathbf{x}}\_{i+1} \cdot \bar{z}\_i - \bar{x}\_i \cdot \bar{z}\_{i+1}) \tag{15}$$

where *n* is the number of points, point *n* + 1 refers to point 1 as the polygon must be closed. Note that the equation 15 holds when the points are ordered in clockwise direction, for the anticlockwise order of points the minus sign has to be added. After computations the area of the bone cross-section is *A* = 393.07*mm*2. Centroid y coordinate, *P*¯ *<sup>x</sup>*, can be obtained from:

$$\bar{P}\_{\mathbf{x}} = \frac{1}{6A} \sum\_{i=1}^{n} (\bar{\mathbf{x}}\_{i} + \bar{\mathbf{x}}\_{i+1}) \cdot (\bar{\mathbf{x}}\_{i+1} \cdot \bar{z}\_{i} - \bar{\mathbf{x}}\_{i} \cdot \bar{z}\_{i+1}) \tag{16}$$

Substituting numerical values to equation 16 yields *P*¯ *<sup>x</sup>* = −1.84*mm*. Global coordinates of point P can be obtained from the geometrical dependencies depicted in figure 8B. Note that *P*¯ *<sup>x</sup>* is a coordinate of the point P with respect to point P' not a distance (absolute value of a coordinate), thus sign convention in the equations 17 and 18.

$$P\_{\mathcal{X}} = \mathbb{C}\_{\mathcal{X}} + |\mathbb{C}P^{\prime}| \cdot \sin(\mathfrak{a}) + \bar{P}\_{\mathcal{X}} \cdot \cos(\mathfrak{a}) \tag{17}$$

$$P\_y = \mathbb{C}\_y + |\mathbb{C}P'| \cdot \cos(\mathfrak{a}) - \bar{P}\_{\mathfrak{X}} \cdot \sin(\mathfrak{a}) \tag{18}$$

Substituting numerical values to equations 17, 18 yields:

$$P\_{\text{\tiny X}} = -159.55 mm + 225.76 mm \cdot \sin(10^{\circ}) - 1.84 mm \cdot \cos(10^{\circ}) = -122.16 mm \tag{19}$$

$$P\_y = -840.22 mm + 225.76 mm \cdot \cos(10^\circ) - (-1.84 mm) \cdot \sin(10^\circ) = -617.57 mm \tag{20}$$

Second moment of area of the bone cross-section about *z*¯ axis can be defined as:

$$I\_{\tilde{z}} = \frac{1}{12} \sum\_{i=1}^{n} (\tilde{\mathbf{x}}\_{i}^{2} + \tilde{\mathbf{x}}\_{i} \cdot \tilde{\mathbf{x}}\_{i+1} + \tilde{\mathbf{x}}\_{i+1}^{2}) \cdot (\tilde{\mathbf{x}}\_{i+1} \cdot \tilde{z}\_{i} - \tilde{\mathbf{x}}\_{i} \cdot \tilde{z}\_{i+1});\tag{21}$$

where:

$$
\mathfrak{X} = \mathfrak{X} - \bar{P}\_{\mathfrak{X}} \tag{22}
$$

Substituting numerical values to equations 21 and 22 yields: *Iz*¯ = 15008*mm*4. Finally, having all cross-sectional parameters, and most importantly coordinates of point P in global

**Point** σ*Fn* **[MPa]** σ*<sup>M</sup>* **[MPa]** σ **[MPa]** � [μ] -0.99 -18.89 -19.88 -1169 -0.99 -15.33 -16.32 -960 -0.99 -11.33 -12.32 -725 -0.99 -3.17 -4.16 -245 -0.99 3.79 2.80 165 -0.99 13.41 12.42 731 -0.99 14.87 13.88 817 -0.99 8.01 7.02 413 -0.99 -9.61 -10.6 -623 -0.99 -14.09 -15.08 -887

Estimating Lower Limb Skeletal Loading 257

Table 2. Axial stresses resulting from normal force, *σFn*, bending moment, *σM*, total stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the cross-section of the bone

global x axis. In order to maintain the desired body position, ankle, knee and hip joints cannot remain unconstrained, thus at each of this joints torques are applied. In the presented case this corresponds to fixing the body segments together rigidly as only static case is discussed. However in the general case, where motion would be studied, those torques would vary over time in order to reproduce subject's motion. The model used in inverse dynamics simulation

Models used in simulation can be more complex comparing to the simple mechanical calculations. For this reason natural mass distribution along with inertia properties for each body segment can be achieved. Reasonable guess of geometry and mass distribution can be obtained from LifeMod package. According to the subject specification the foot weight 2.14kg, lower leg mass is 4.19kg, upper leg's mass is 7.47kg and the upper body plus pelvis weight 52.4kg. As only half of the model is considered, the upper body and pelvis mass have to be reduced to 26.2kg to maintain symmetry conditions. Stresses and strains calculations could be performed in analogical way like in the model 1 using the multibody model just to obtain the cross-sectional bending moments and normal force. However use of flexible multibody dynamics allows to compute the stresses and strains directly at any nodal location within the flexible bone. To benefit from the flexible multibody approach, the rigid bone at which stresses and strains are to be computed needs to be replaced with it's flexible version. To obtain the most accurate results bone geometry should be reconstructed from medical images to correspond to the subject, however for this example just the simple shell geometry created by LifeMod will be used. After exporting the geometry from Adams5, it can be imported to any finite element software which is capable of generating modal neutral file. For this example

In the finite element software the initial shell geometry is converted to solid body which can be discretized with solid elements. As in the model 1, it is assumed that bone is homogeneous, thus the whole volume of the bone is meshed with 5mm in size four-node tetrahedral linear solid elements. Each of the elements is assigned Young modulus of 17GPa and Poisson's ratio of 0.3. In order to perform modal analysis rigid beams are added at the distal ends of tibia connecting surface nodes with knee and ankle joint centres respectively. Craig-Bampton (Craig & Bampton, 1968) modal reduction allows to represent large finite element models with the use of several deformation modes, thus remarkable reducing the computational effort in

<sup>5</sup> MSC Adams, version R3, MSC.Software Corporation, Santa Ana, California, USA

is presented in figure 9.

ANSYS 12 will be used.

coordinate system, equilibrium equations 12, 13 and 14, can be solved for: bending moment, *M*, normal force, *Fn* and tangential force, *Ft* as follows:

$$F\_l = 392.43N \cdot \sin(10^\circ) - 19.28N \cdot \cos(10^\circ) = 49.16N \tag{23}$$

$$F\_{\rm H} = 19.28 \text{N} \cdot \sin(10^{\circ}) + 392.43 \text{N} \cdot \cos(10^{\circ}) = 389.82 \text{N} \tag{24}$$

$$M = 19.28 \text{N} \cdot (-617.57 mm - (-840.22 mm))\tag{25}$$

$$(1 + 392.43N \cdot \stackrel{\frown}{(-122.16mm - (-159.55mm))} = 18966 \text{Nmm} \tag{25}$$

Each point within the cross-section of the bone is subjected to stresses caused by the normal force, *Fn*, and bending moment, *M*. In mechanics tensile stresses are denoted by positive sign, consequently compressive strains are marked by the negative sign. Compressive axial stresses caused by the normal force are uniform within the cross-section and can be calculated using the following relationship:

$$
\sigma\_{\mathcal{F}\_u} = -\frac{F\_n}{A} \tag{26}
$$

Stresses caused by the bending moment vary within the cross-section with respect to the coordinate *x*˜. The relationship between point location and bending stresses is:

$$
\sigma\_M = -\frac{M \cdot \mathfrak{X}}{I\_{\tilde{z}}} \tag{27}
$$

Total axial strains are sum of stresses caused by the normal force and stresses caused by bending:

$$
\sigma = \sigma\_{\text{F}\_{\text{\tiny H}}} + \sigma\_{\text{M}} = -\frac{F\_{\text{H}}}{A} - \frac{M \cdot \overline{\mathfrak{X}}}{I\_{\text{\tiny Z}}} \tag{28}
$$

Substituting numerical values to equation 26 results in:

$$
\sigma\_{\text{F}\_{\text{u}}} = -\frac{389.82N}{393.07mm^2} = -0.99 MPa\tag{29}
$$

Strains are associated with stresses with a linear relationship called the Hook's law:

$$
\epsilon = \frac{\sigma}{E} \tag{30}
$$

Table 2 contains axial stresses computed at the vertices of bone's cross-section using equations 26, 27, and 28. In the last column of the table total strains are presented computed using the Hooke's law (equation 30). It can be seen that the posterior part of the bone is subjected to tensile strains while the anterior part is in compression. Bending contribution in the loading is dominating, which can be seen from analysis of table 2.

### **4.2 Inverse dynamics simulation (model 2)**

Inverse dynamics simulation is similar to the calculations presented in the previous section. The main difference between the two approaches is the mass distribution within the model and complexity of the mechanical system. Moreover the inverse dynamics model, presented in this section, is three dimensional. This allows to account for bending also in the frontal (coronal) (and transverse) plane(s). In model 2 foot is attached to the ground by a revolute joint. Foot, lower leg and pelvis are interconnected via spherical joints, and pelvis at the center of mass is attached to the ground by planar joint in the sagittal plane. Stabilization force is applied at the center of mass of the body by a constraint restricting translation along the 14 Will-be-set-by-IN-TECH

coordinate system, equilibrium equations 12, 13 and 14, can be solved for: bending moment,

Each point within the cross-section of the bone is subjected to stresses caused by the normal force, *Fn*, and bending moment, *M*. In mechanics tensile stresses are denoted by positive sign, consequently compressive strains are marked by the negative sign. Compressive axial stresses caused by the normal force are uniform within the cross-section and can be calculated using

*<sup>σ</sup>Fn* <sup>=</sup> <sup>−</sup> *Fn*

Stresses caused by the bending moment vary within the cross-section with respect to the

*<sup>σ</sup><sup>M</sup>* <sup>=</sup> <sup>−</sup> *<sup>M</sup>* · *<sup>x</sup>*˜ *Iz*¯

Total axial strains are sum of stresses caused by the normal force and stresses caused by

*<sup>σ</sup>* <sup>=</sup> *<sup>σ</sup>Fn* <sup>+</sup> *<sup>σ</sup><sup>M</sup>* <sup>=</sup> <sup>−</sup> *Fn*

Strains are associated with stresses with a linear relationship called the Hook's law:

*�* <sup>=</sup> *<sup>σ</sup>*

Table 2 contains axial stresses computed at the vertices of bone's cross-section using equations 26, 27, and 28. In the last column of the table total strains are presented computed using the Hooke's law (equation 30). It can be seen that the posterior part of the bone is subjected to tensile strains while the anterior part is in compression. Bending contribution in the loading

Inverse dynamics simulation is similar to the calculations presented in the previous section. The main difference between the two approaches is the mass distribution within the model and complexity of the mechanical system. Moreover the inverse dynamics model, presented in this section, is three dimensional. This allows to account for bending also in the frontal (coronal) (and transverse) plane(s). In model 2 foot is attached to the ground by a revolute joint. Foot, lower leg and pelvis are interconnected via spherical joints, and pelvis at the center of mass is attached to the ground by planar joint in the sagittal plane. Stabilization force is applied at the center of mass of the body by a constraint restricting translation along the

*<sup>σ</sup>Fn* <sup>=</sup> <sup>−</sup> 389.82*<sup>N</sup>*

*M* = 19.28*N* · (−617.57*mm* − (−840.22*mm*))

coordinate *x*˜. The relationship between point location and bending stresses is:

*Ft* = 392.43*N* · sin(10◦) − 19.28*N* · cos(10◦) = 49.16*N* (23)

*Fn* = 19.28*N* · sin(10◦) + 392.43*N* · cos(10◦) = 389.82*N* (24)

<sup>+</sup> 392.43*<sup>N</sup>* · (−122.16*mm* <sup>−</sup> (−159.55*mm*)) = <sup>18966</sup>*Nmm* (25)

*<sup>A</sup>* <sup>−</sup> *<sup>M</sup>* · *<sup>x</sup>*¯ *Iz*¯

393.07*mm*<sup>2</sup> <sup>=</sup> <sup>−</sup>0.99*MPa* (29)

*<sup>E</sup>* (30)

*<sup>A</sup>* (26)

(27)

(28)

*M*, normal force, *Fn* and tangential force, *Ft* as follows:

Substituting numerical values to equation 26 results in:

is dominating, which can be seen from analysis of table 2.

**4.2 Inverse dynamics simulation (model 2)**

the following relationship:

bending:


Table 2. Axial stresses resulting from normal force, *σFn*, bending moment, *σM*, total stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the cross-section of the bone

global x axis. In order to maintain the desired body position, ankle, knee and hip joints cannot remain unconstrained, thus at each of this joints torques are applied. In the presented case this corresponds to fixing the body segments together rigidly as only static case is discussed. However in the general case, where motion would be studied, those torques would vary over time in order to reproduce subject's motion. The model used in inverse dynamics simulation is presented in figure 9.

Models used in simulation can be more complex comparing to the simple mechanical calculations. For this reason natural mass distribution along with inertia properties for each body segment can be achieved. Reasonable guess of geometry and mass distribution can be obtained from LifeMod package. According to the subject specification the foot weight 2.14kg, lower leg mass is 4.19kg, upper leg's mass is 7.47kg and the upper body plus pelvis weight 52.4kg. As only half of the model is considered, the upper body and pelvis mass have to be reduced to 26.2kg to maintain symmetry conditions. Stresses and strains calculations could be performed in analogical way like in the model 1 using the multibody model just to obtain the cross-sectional bending moments and normal force. However use of flexible multibody dynamics allows to compute the stresses and strains directly at any nodal location within the flexible bone. To benefit from the flexible multibody approach, the rigid bone at which stresses and strains are to be computed needs to be replaced with it's flexible version. To obtain the most accurate results bone geometry should be reconstructed from medical images to correspond to the subject, however for this example just the simple shell geometry created by LifeMod will be used. After exporting the geometry from Adams5, it can be imported to any finite element software which is capable of generating modal neutral file. For this example ANSYS 12 will be used.

In the finite element software the initial shell geometry is converted to solid body which can be discretized with solid elements. As in the model 1, it is assumed that bone is homogeneous, thus the whole volume of the bone is meshed with 5mm in size four-node tetrahedral linear solid elements. Each of the elements is assigned Young modulus of 17GPa and Poisson's ratio of 0.3. In order to perform modal analysis rigid beams are added at the distal ends of tibia connecting surface nodes with knee and ankle joint centres respectively. Craig-Bampton (Craig & Bampton, 1968) modal reduction allows to represent large finite element models with the use of several deformation modes, thus remarkable reducing the computational effort in

<sup>5</sup> MSC Adams, version R3, MSC.Software Corporation, Santa Ana, California, USA

**Point** σ **[MPa]** � [μ] 1 -13.04 -755 2 -12.14 -713 3 -8.38 -491 4 -4.05 -245 5 2.44 145 6 7.64 447 7 8.65 508 8 4.13 245 9 -7.83 -450 10 -10.56 -617

Estimating Lower Limb Skeletal Loading 259

Table 3. Total axial stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the

Forward dynamics simulation with the use of flexible multibody dynamics is based on the model presented in inverse dynamics example. However in the model 3, torque actuators at the joints are removed, and in turn muscle actuators are added. Muscle actuators are represented by straight line of action between muscle insertion points. The contraction of each muscle is controlled by the proportional-integral-derivative (PID) controller, which aims at replicating muscle contraction patterns obtained in the inverse dynamics. The model is actuated by 16 muscles: tibialis anterior, soleus, gastrocnemius medialis, gastrocnemius lateralis, rectus femoris, vastus lateralis, vastus medialis, biceps femoris (two actuators), illiacus, semitendinosus, gluteus medialis (two actuators), gluteus maximus (two actuators), adductor magnus. Muscle placement is depicted in figure 10. Adductor magnus is not visible

Forward dynamics simulation requires significant amount of input parameters, among them muscle contraction patterns, that are used as a reference signal for muscle controllers. The easiest way to obtain muscle contraction splines is the inverse dynamics simulation. During the inverse dynamics model is driven by the motion captured from a subject. In the example case, there is no motion, so during the inverse dynamics the model is simply fixed in desired orientation. After gathering muscle contraction splines, which in this case are just initial muscle length values, the forward dynamics simulation can be performed. External constrains are removed, meaning that the joints are free to move if no muscle force is produced. However support hinge at the foot, planar constraint at the pelvis and translation along x axis constrain

Simulated ground reaction forces are: *GRFn* = 392.22*N* and *GRFf* = 20.57*N* for the vertical and frictional components respectively. Stresses and strains at the mid-shaft cross-section of

It is worth noting that the loading pattern obtained from model 3 is completely different to the ones obtained using models 1 and 2. Model 3 predicts that the posterior part of tibia is subjected to compression while the anterior part is in tension. The magnitude of strains and stresses is around 7.62 times smaller for model 3 compared to model 2. The differences can be explained by analysing differences in mechanics of model 2 and 3. The first one behaves like a rigid structure, where sagittal ankle torque is carried as bending moment by the tibia. Model

cross-section of the bone obtained from inverse dynamics model

**4.3 Forward dynamics simulation (model 3)**

in the picture as it is on the inferior site of the femur.

at the center of mass are preserved.

the tibia are presented in table 4.

**4.4 Summary of the numerical example**

Fig. 9. Illustration of the three-dimensional model used in the inverse dynamics simulation. The model consists of four segments: foot, lower and upper leg, and a pelvis.

the simulation software. Knee and ankle joint center location nodes are used as boundary nodes in the modal analysis. For the current example 42 deformation modes will be computed. In the simulation however only 36 modes are usable, due to the fact that first six deformation modes represent rigid body motion, which is accounted for already in multibody formulation. Ground reaction forces obtained from the model are *GRFn* = 392.27*N* and *GRFf* = 20.50*N*. Stress and strain results of inverse dynamics simulation with flexible tibia are grouped in table 3. Comparing table 2 containing results from the simple model with the current results it can be seen that stresses and strains are on average 1.43 times higher in case of a simple model. At the same time, both models predict same deformation mode of the bone, which is bending in the saggital plane causing compressive strains at the anterior site of the bone and tensile strains on the posterior side.

16 Will-be-set-by-IN-TECH

Fig. 9. Illustration of the three-dimensional model used in the inverse dynamics simulation.

the simulation software. Knee and ankle joint center location nodes are used as boundary nodes in the modal analysis. For the current example 42 deformation modes will be computed. In the simulation however only 36 modes are usable, due to the fact that first six deformation modes represent rigid body motion, which is accounted for already in multibody formulation. Ground reaction forces obtained from the model are *GRFn* = 392.27*N* and *GRFf* = 20.50*N*. Stress and strain results of inverse dynamics simulation with flexible tibia are grouped in table 3. Comparing table 2 containing results from the simple model with the current results it can be seen that stresses and strains are on average 1.43 times higher in case of a simple model. At the same time, both models predict same deformation mode of the bone, which is bending in the saggital plane causing compressive strains at the anterior site of the bone and

The model consists of four segments: foot, lower and upper leg, and a pelvis.

tensile strains on the posterior side.


Table 3. Total axial stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the cross-section of the bone obtained from inverse dynamics model

### **4.3 Forward dynamics simulation (model 3)**

Forward dynamics simulation with the use of flexible multibody dynamics is based on the model presented in inverse dynamics example. However in the model 3, torque actuators at the joints are removed, and in turn muscle actuators are added. Muscle actuators are represented by straight line of action between muscle insertion points. The contraction of each muscle is controlled by the proportional-integral-derivative (PID) controller, which aims at replicating muscle contraction patterns obtained in the inverse dynamics. The model is actuated by 16 muscles: tibialis anterior, soleus, gastrocnemius medialis, gastrocnemius lateralis, rectus femoris, vastus lateralis, vastus medialis, biceps femoris (two actuators), illiacus, semitendinosus, gluteus medialis (two actuators), gluteus maximus (two actuators), adductor magnus. Muscle placement is depicted in figure 10. Adductor magnus is not visible in the picture as it is on the inferior site of the femur.

Forward dynamics simulation requires significant amount of input parameters, among them muscle contraction patterns, that are used as a reference signal for muscle controllers. The easiest way to obtain muscle contraction splines is the inverse dynamics simulation. During the inverse dynamics model is driven by the motion captured from a subject. In the example case, there is no motion, so during the inverse dynamics the model is simply fixed in desired orientation. After gathering muscle contraction splines, which in this case are just initial muscle length values, the forward dynamics simulation can be performed. External constrains are removed, meaning that the joints are free to move if no muscle force is produced. However support hinge at the foot, planar constraint at the pelvis and translation along x axis constrain at the center of mass are preserved.

Simulated ground reaction forces are: *GRFn* = 392.22*N* and *GRFf* = 20.57*N* for the vertical and frictional components respectively. Stresses and strains at the mid-shaft cross-section of the tibia are presented in table 4.

### **4.4 Summary of the numerical example**

It is worth noting that the loading pattern obtained from model 3 is completely different to the ones obtained using models 1 and 2. Model 3 predicts that the posterior part of tibia is subjected to compression while the anterior part is in tension. The magnitude of strains and stresses is around 7.62 times smaller for model 3 compared to model 2. The differences can be explained by analysing differences in mechanics of model 2 and 3. The first one behaves like a rigid structure, where sagittal ankle torque is carried as bending moment by the tibia. Model

**Point** σ **[MPa]** � [μ] 1 1.08 60 2 0.64 37 3 -0.89 -49 4 -2.42 -137 5 -4.73 -275 6 -5.69 -333 7 -5.26 -309 8 -2.75 -164 9 0.19 6 10 0.57 31

Estimating Lower Limb Skeletal Loading 261

Table 4. Total axial stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the

stress shielding for the bone Sverdlova & Witzel (2010). When it comes to the external forces

Estimating lower limb skeletal loading may be useful e.g. in association to planning osteogenic exercise interventions. Although estimating skeletal loading from ground reaction forces is tempting due to the simplicity of the measurement, more information may be gained from more elaborate approaches. As was demonstrated by the numerical example, accounting for mass distribution in the model (as shown in model 2) can lead to 143% difference in strains. Further difference of up to 19 times can be noted when internal forces are considered in addition to external forces (i.e. joint torque actuators are replaced with muscle models). Even though all of the methods predicted almost identical stabilization force and produced the same ground reaction forces, internal loading was very different. For example in implant design only relatively realistic strain estimates may be expected to be useful. Implant loosening is a well known phenomenon (van Loon et al., 1999) and designing the implant in such a way that appropriate loading levels are attained can be therefore argued to be quite relevant. Furthermore, if a specific (for example clinically relevant) bone site is to be targeted with an exercise intervention, it may be necessary to the estimate the loading in a site specific

Anderson, D. D., Hillberry, B. M., Teegarden, D., Proulx, W. R., Weaver, C. M. & Yoshikawa,

the hip joint and femoral neck, *Journal of Applied Biomechanics* 12: 292–312. Burr, D. B., Milgrom, C., Fyhrie, D., Forwood, M., Nyska, M., Finestone, A., Hoshaw, S., Saiag,

Cheng, S., Suominen, H., Sakari-Rantala, R., Laukkanen, P., Avikainen, V. & Heikkinen, E.

follow-up study in elderly people., *J Bone Miner Res* 12(7): 1075–1082.

URL: *http://dx.doi.org/10.1359/jbmr.1997.12.7.1075*

W. (1996). Biomechanical analysis of an exercise program for forces and stresses in

E. & Simkin, A. (1996). In vivo measurement of human tibial strains during vigorous

(1997). Calcaneal bone mineral density predicts fracture occurrence: a five-year

cross-section of the bone obtained from forward dynamics model

they were predicted correctly by all three models.

activity, *Bone* 18(5): 405–410.

**5. Conclusion**

manner.

**6. References**

Fig. 10. Illustration of the three-dimensional model used in the forward dynamics simulation. The model consists of four segments: foot, lower and upper leg, and a pelvis; and fifteen muscles.

1 behaves the same way. Comparing the bending stresses at the considered cross-section to normal stresses one can notice that the bending stresses are up to 19 times higher than the normal stresses. In contrast, in model 3 tibia, foot and muscles create a mechanism. Due to the fact that ankle joint is almost frictionless it carries almost no torque at all, and the contact force is transmitted through Achilles tendon, and through ankle joint just as an axial and tangential force. Tibia is then subjected to only transverse and normal forces at the knee and ankle joints, and the bending moment comes mostly from the activity of soleus (545.72N), vastus lateralis (67.6N), vastus medialis (2.39N) and rectus femoris (0.08N) muscles. All of those muscles while contracting create bending moment which creates compressive strains at the posterior site of tibia. Tensile strains from bending moment actually reduce the compressive loading from axial force at the anterior site of tibia. This example demonstrates 1) that internal forces should not be neglected in estimating skeletal loads and 2) that muscle contraction can provide 18 Will-be-set-by-IN-TECH

Fig. 10. Illustration of the three-dimensional model used in the forward dynamics simulation. The model consists of four segments: foot, lower and upper leg, and a pelvis;

1 behaves the same way. Comparing the bending stresses at the considered cross-section to normal stresses one can notice that the bending stresses are up to 19 times higher than the normal stresses. In contrast, in model 3 tibia, foot and muscles create a mechanism. Due to the fact that ankle joint is almost frictionless it carries almost no torque at all, and the contact force is transmitted through Achilles tendon, and through ankle joint just as an axial and tangential force. Tibia is then subjected to only transverse and normal forces at the knee and ankle joints, and the bending moment comes mostly from the activity of soleus (545.72N), vastus lateralis (67.6N), vastus medialis (2.39N) and rectus femoris (0.08N) muscles. All of those muscles while contracting create bending moment which creates compressive strains at the posterior site of tibia. Tensile strains from bending moment actually reduce the compressive loading from axial force at the anterior site of tibia. This example demonstrates 1) that internal forces should not be neglected in estimating skeletal loads and 2) that muscle contraction can provide

and fifteen muscles.


Table 4. Total axial stresses, *σ*, and total strains, *�*, at the vertices of the polygon defining the cross-section of the bone obtained from forward dynamics model

stress shielding for the bone Sverdlova & Witzel (2010). When it comes to the external forces they were predicted correctly by all three models.

### **5. Conclusion**

Estimating lower limb skeletal loading may be useful e.g. in association to planning osteogenic exercise interventions. Although estimating skeletal loading from ground reaction forces is tempting due to the simplicity of the measurement, more information may be gained from more elaborate approaches. As was demonstrated by the numerical example, accounting for mass distribution in the model (as shown in model 2) can lead to 143% difference in strains. Further difference of up to 19 times can be noted when internal forces are considered in addition to external forces (i.e. joint torque actuators are replaced with muscle models). Even though all of the methods predicted almost identical stabilization force and produced the same ground reaction forces, internal loading was very different. For example in implant design only relatively realistic strain estimates may be expected to be useful. Implant loosening is a well known phenomenon (van Loon et al., 1999) and designing the implant in such a way that appropriate loading levels are attained can be therefore argued to be quite relevant. Furthermore, if a specific (for example clinically relevant) bone site is to be targeted with an exercise intervention, it may be necessary to the estimate the loading in a site specific manner.

### **6. References**


Huston, R. L. (1991). Computer methods in flexible multibody dynamics, *Int. J. Numer. Meth.*

Estimating Lower Limb Skeletal Loading 263

Ishikawa, M., Komi, P. V., Grey, M. J., Lepola, V. & Bruggemann, G. P. (2005). Muscle-tendon

Jämsä, T., Vainionpää, A., Korpelainen, R., Vihriälä, E. & Leppäluoto, J. (2006). Effect of daily physical activity on proximal femur, *Clin.Biomech.(Bristol, Avon)* 21(1): 1–7. Järvinen, T. L., Sievänen, H., Khan, K. M., Heinonen, A. & Kannus, P. (2008). Shifting the

Kamman, J. W. & Huston, R. L. (2001). Multibody dynamics modeling of variable length cable

Kemmler, W., von Stengel, S., Engelke, K., Häberle, L. & Kalender, W. A. (2010). Exercise

Klodowski, A., Rantalainen, T., Mikkola, A., Heinonen, A. & Sievänen, H. (2011). Flexible

Komi, P. V., Fukashiro, S. & Järvinen, M. (1992). Biomechanical loading of achilles tendon

Korpelainen, R., Keinanen-Kiukaanniemi, S., Heikkinen, J., Vaananen, K. & Korpelainen, J.

Lanyon, L. E. (1987). Functional strain in bone tissue as an objective, and controlling stimulus for adaptive bone remodelling, *Journal of Biomechanics* 20(11-12): 1083–1093. Lanyon, L. E., Hampson, W. G., Goodship, A. E. & Shah, J. S. (1975). Bone deformation

Leppänen, O. V., Sievänen, H., Jokihaara, J., Pajamäki, I., Kannus, P., Cooper, D. M. & Järvinen,

Lowen, G. & Chassapis, C. (1986). The elastic behavior of linkages: An update, *Mechanism and*

Lowen, G. G. & Jandrasits, W. G. (1972). Survey of investigations into the dynamic behavior

during normal locomotion, *Clinics in sports medicine* 11(3): 521–531.

interaction and elastic energy usage in human walking, *Journal of Applied Physiology*

focus in fracture prevention from osteoporosis to falls, *BMJ (Clinical research ed.)*

effects on bone mineral density, falls, coronary risk factors, and health care costs in older women: the randomized controlled senior fitness and prevention (sefip) study.,

multibody approach in forward dynamic simulation of locomotive strains in human skeleton with flexible lower body bones, *Multibody System Dynamics* 25: 395–409.

(2006). Effect of exercise on extraskeletal risk factors for hip fractures in elderly women with low bmd: a population-based randomized controlled trial, *Journal of*

recorded in vivo from strain gauges attached to the human tibial shaft, *Acta*

T. L. N. (2010). The effects of loading and estrogen on rat bone growth., *J Appl Physiol*

URL: *http://www.sciencedirect.com/science/article/B6V46-482GN2K-NK/2/4ab4249996dbf*

of mechanisms containing links with distributed mass and elasticity, *Mechanism and*

URL: *http://www.sciencedirect.com/science/article/B6V46-482GN85-T/2/037ade3ebbe8696*

*Engng.* 32(8): 1657–1668.

99(2): 603–608.

336(7636): 124–126.

*Arch Intern Med* 170(2): 179–185.

URL: *http://dx.doi.org/10.1002/nme.1620320809*

systems, *Multibody System Dynamics* 5: 211–221. URL: *http://dx.doi.org/10.1023/A:1011489801339*

URL: *http://dx.doi.org/10.1001/archinternmed.2009.499*

URL: *http://dx.doi.org/10.1007/s11044-010-9240-9*

*Bone and Mineral Research* 21(5): 772–779.

*Orthopaedica Scandinavica* 46(2): 256–268.

URL: *http://dx.doi.org/10.1152/japplphysiol.00989.2009*

108(6): 1737–1744.

*1be022bcf551d34f9cb*

*34d8d47cb9979d1fa*

*Machine Theory* 21(1): 33 – 42.

*Machine Theory* 7(1): 3 – 17.


URL: *http://dx.doi.org/10.1056/NEJMoa0809493*


URL: *http://dx.doi.org/10.1002/jbmr.5650110410*


URL: *http://dx.doi.org/10.1002/nme.1620151210*

Huston, R. L. (1981). Multi-body dynamics including the effects of flexibility and compliance, *Computers & Structures* 14(5-6): 443 – 451. URL: *http://www.sciencedirect.com/science/article/B6V28-47X8NJK-11V/2/5a675c9ea3c09 53ea466e0caaabb4861*

20 Will-be-set-by-IN-TECH

Craig, R. R. & Bampton, M. C. C. (1968). Coupling of substructures for dynamic analysis, *American Institute of Aeronautics and Astronautics Journal* 6(7): 1313–1319. Cummings, S. R., Martin, J. S., McClung, M. R., Siris, E. S., Eastell, R., Reid, I. R., Delmas, P.,

Cummings, S. R. & Melton, L. J. (2002). Epidemiology and outcomes of osteoporotic fractures,

Currey, J. D. (2001). Bone strength: what are we trying to measure?, *Calcified tissue international*

Currey, J. D. (2002). *Bones: structure and mechanics*, Princeton University Press, Princeton,

Dong, X. N. & Guo, X. E. (2004). The dependence of transversely isotropic elasticity of human femoral cortical bone on porosity, *Journal of Biomechanics* 37(8): 1281–1287. Duodecim, S. L. (2008). Osteoporoosi, *Technical report*, Suomalainen Lääkäriseura Duodecim.

Erdman, A. G. & Sandor, G. N. (1972). Kineto-elastodynamics–areview of the state of the art

Finni, T., Komi, P. V. & Lukkariniemi, J. (1998). Achilles tendon loading during walking:

Frost, H. M. (2000). Muscle, bone, and the utah paradigm: a 1999 overview, *Medicine and*

Frost, H. M. (2003). Bone's mechanostat: a 2003 update, *The anatomical record.Part A, Discoveries*

Heinonen, A., Sievänen, H., Kannus, P., Oja, P. & Vuori, I. (1996). Effects of unilateral

Hoshaw, S. J., Fyhrie, D. P., Takano, Y., Burr, D. B. & Milgrom, C. (1997). A method suitable for

Hughes, P. C. (1979). Dynamics of a chain of flexible bodies, *Journal of the Astronautical Sciences*

Hughes, T. J. R. & Winget, J. (1980). Finite rotation effects in numerical integration of rate

Huston, R. L. (1981). Multi-body dynamics including the effects of flexibility and compliance,

*in molecular, cellular, and evolutionary biology* 275(2): 1081–1101.

URL: *http://www.sciencedirect.com/science/article/B6V46-482GN85-V/2/ad26341e4743520*

application of a novel optic fiber technique, *European journal of applied physiology and*

strength training and detraining on bone mineral mass and estimated mechanical characteristics of the upper limb bones in young women., *J Bone Miner Res*

in vivo measurement of bone strain in humans, *Journal of Biomechanics* 30(5): 521–524.

constitutive equations arising in large-deformation analysis, *Int. J. Numer. Meth.*

URL: *http://www.sciencedirect.com/science/article/B6V28-47X8NJK-11V/2/5a675c9ea3c09*

and trends, *Mechanism and Machine Theory* 7(1): 19 – 33.

361(8): 756–765.

68(4): 205–210.

United Kingdom.

*3d04178013c171e71*

11(4): 490–501.

27: 359–380.

*Engng.* 15(12): 1862–1867.

*53ea466e0caaabb4861*

*Lancet* 359(9319): 1761–1767.

URL: *http://www.kaypahoito.fi/*

*occupational physiology* 77(3): 289–291.

*science in sports and exercise* 32(5): 911–917.

URL: *http://dx.doi.org/10.1002/jbmr.5650110410*

URL: *http://dx.doi.org/10.1002/nme.1620151210*

*Computers & Structures* 14(5-6): 443 – 451.

URL: *http://dx.doi.org/10.1056/NEJMoa0809493*

Zoog, H. B., Austin, M., Wang, A., Kutilek, S., Adami, S., Zanchetta, J., Libanati, C., Siddhanti, S., Christiansen, C. & Trial, F. R. E. E. D. O. M. (2009). Denosumab for prevention of fractures in postmenopausal women with osteoporosis., *N Engl J Med* Huston, R. L. (1991). Computer methods in flexible multibody dynamics, *Int. J. Numer. Meth. Engng.* 32(8): 1657–1668.

```
URL: http://dx.doi.org/10.1002/nme.1620320809
```

URL: *http://dx.doi.org/10.1001/archinternmed.2009.499*


URL: *http://dx.doi.org/10.1152/japplphysiol.00989.2009*

Lowen, G. & Chassapis, C. (1986). The elastic behavior of linkages: An update, *Mechanism and Machine Theory* 21(1): 33 – 42.

URL: *http://www.sciencedirect.com/science/article/B6V46-482GN2K-NK/2/4ab4249996dbf 1be022bcf551d34f9cb*

Lowen, G. G. & Jandrasits, W. G. (1972). Survey of investigations into the dynamic behavior of mechanisms containing links with distributed mass and elasticity, *Mechanism and Machine Theory* 7(1): 3 – 17.

URL: *http://www.sciencedirect.com/science/article/B6V46-482GN85-T/2/037ade3ebbe8696 34d8d47cb9979d1fa*

Robertson, D. G. E. (2004). *Research methods in biomechanics*, Human Kinetics, Champaign, IL :.

Estimating Lower Limb Skeletal Loading 265

Rubin, C. T. & Lanyon, L. E. (1984). Regulation of bone formation by applied dynamic loads,

Shabana, A. A. (2005). *Dynamics of multibody systems*, Cambridge University Press, Cambridge

Shabana, A. A.; Wehage, R. A. (1983). Coordinate reduction technique for transient analysis

Sievänen, H. (2005). Hormonal influences on the muscle-bone feedback system: a perspective,

Silder, A., Heiderscheit, B. & Thelen, D. G. (2008). Active and passive contributions to joint kinetics during walking in older adults, *Journal of Biomechanics* 41(7): 1520–1527. Sinaki, M., Itoi, E., Wahner, H. W., Wollan, P., Gelzcer, R., Mullan, B. P., Collins, D. A. &

Stevens, J. A. & Olson, S. (2000). Reducing falls and resulting hip fractures among older

Sverdlova, N. S. & Witzel, U. (2010). Principles of determination and verification of muscle

Turcic, D. A. & Midha, A. (1984). Dynamic analysis of elastic mechanism systems. part i: Applications, *Journal of Dynamic Systems, Measurement, and Control* 106(4): 249–254.

Turcic, D. A., Midha, A. & Bosnik, J. R. (1984). Dynamic analysis of elastic mechanism systems.

Turner, C. H. (1998). Three rules for bone adaptation to mechanical stimuli, *Bone*

Turner, C. H. & Burr, D. B. (1993). Basic biomechanical measurements of bone: a tutorial, *Bone*

Turner, C. H. & Robling, A. G. (2003). Designing exercise regimens to increase bone strength,

Umemura, Y., Sogo, N. & Honda, A. (2002). Effects of intervals between jumps or bouts on osteogenic response to loading, *Journal of Applied Physiology* 93(4): 1345–1348. Vainionpää, A. (2007). Bone adaptation to impact loading : significance of loading intensity. Vainionpää, A., Korpelainen, R., Leppäluoto, J. & Jämsä, T. (2005). Effects of high-impact

Vainionpää, A., Korpelainen, R., Vihriälä, E., Rinta-Paavola, A., Leppäluoto, J. & Jämsä, T.

van Loon, C. J., de Waal Malefijt, M. C., Buma, P., Verdonschot, N. & Veth, R. P. (1999). Femoral

exercise on bone mineral density: a randomized controlled trial in premenopausal

(2006). Intensity of exercise is associated with bone density change in premenopausal

bone loss in total knee arthroplasty. a review, *Acta Orthopaedica Belgica* 65(2): 154–163.

URL: *http://libtest.csc.fi:8991/F/?func=direct&doc\_number=002760407&local\_base=fin01*

of spatial substructures with large angular rotations, *Journal of Structural Mechanics*

Hodgson, S. F. (2002). Stronger back muscles reduce the incidence of vertebral fractures: a prospective 10 year follow-up of postmenopausal women., *Bone*

forces in the human musculoskeletal system: Muscle forces to minimise bending

part ii: Experimental results, *Journal of Dynamic Systems, Measurement, and Control*

URL: *http://www.loc.gov/catdir/toc/ecip0414/2004002355.html*

:.

11(3): 401–431.

30(6): 836–841.

106(4): 255–260.

23(5): 399–407.

14(4): 595–608.

*Journal of Bone and Joint Surgery. American volume* 66(3): 397–402.

*Journal of Musculoskeletal & Neuronal Interactions* 5(3): 255–261.

women, *MMWR Recomm Rep.* 49(RR-2): 3–12.

stress, *Journal of Biomechanics* 43(3): 387–396.

URL: *http://link.aip.org/link/?JDS/106/249/1*

URL: *http://link.aip.org/link/?JDS/106/255/1*

*Exercise and sport sciences reviews* 31(1): 45–50.

women, *Osteoporosis International* 16(2): 191–197.

women, *Osteoporosis International* 17(3): 455–463.


URL: *http://dx.doi.org/10.1080/10255842.2010.488223*

	- URL: *http://dx.doi.org/10.1007/s11044-008-9120-8*

22 Will-be-set-by-IN-TECH

Ma, H., Leskinen, T., Alen, M., Cheng, S., Sipila, S., Heinonen, A., Kaprio, J., Suominen, H. &

Martin, R. B. & Burr, D. B. (1989). *Structure, function, and adaptation of compact bone*, Raven

Mikkola, T. M., Sipilä, S., Rantanen, T., Sievänen, H., Suominen, H., Kaprio, J., Koskenvuo,

Milgrom, C., Finestone, A., Levi, Y., Simkin, A., Ekenman, I., Mendelson, S., Millgram, M.,

tibial strains than running?, *British journal of sports medicine* 34(3): 195–199. Nazer, R. A., Klodowski, A., Rantalainen, T., Heinonen, A., Sievanen, H. & Mikkola, A. (2011).

Nazer, R. A., Klodowski, A., Rantalainen, T., Heinonen, A., Sievänen, H. & Mikkola, A.

Nazer, R. A., Rantalainen, T., Heinonen, A., Sievänen, H. & Mikkola, A. (2008). Flexible

Nikander, R., Sievänen, H., Heinonen, A., Daly, R. M., Uusi-Rasi, K. & Kannus, P. (2010).

Ortiz-Luna, G., Garcia-Hernandez, P. & Tamayo-Orozco, J. A. (2009). Treatment options for

Perttunen, J. O., Kyröläinen, H., Komi, P. V. & Heinonen, A. (2000). Biomechanical loading in

Peterman, M. M., Hamel, A. J., Cavanagh, P. R., Piazza, S. J. & Sharkey, N. A. (2001). In

Raman-Nair, W. & Baddour, R. (2003). Three-dimensional dynamics of a flexible marine riser undergoing large elastic deformations, *Multibody System Dynamics* 10: 393–423.

Rantalainen, T., Nikander, R., Heinonen, A., Suominen, H. & Sievanen, H. (2010).

Reilly, D. T. & Burstein, A. H. (1975). The elastic and ultimate properties of compact bone

A twin study, *Journal of Bone and Mineral Research* .

*Journal of Bone and Mineral Research* 23(4): 492–498.

URL: *http://dx.doi.org/10.1080/10255842.2010.488223*

*Multibody System Dynamics* 20: 287–306.

*of Biomechanics* 41(5): 1036–1043.

*Biomechanics* 34(5): 693–698.

URL: *http://dx.doi.org/10.1007/s11044-008-9120-8*

URL: *http://dx.doi.org/10.1186/1741-7015-8-47*

optimising bone strength throughout life., *BMC Med* 8: 47.

the triple jump, *Journal of sports sciences* 18(5): 363–370.

URL: *http://dx.doi.org/10.1023/A:1026213630987*

*Calcified tissue international* 86(6): 447–454.

tissue, *Journal of Biomechanics* 8(6): 393–405.

Press, New York.

*Biomed Engin* p. 1.

1: S114–25.

Kujala, U. M. (2009). Long-term leisure time physical activity and properties of bone:

M., Kauppinen, M. & Heinonen, A. (2008). Genetic and environmental influence on structural strength of weight-bearing and non-weight-bearing bone: a twin study,

Nyska, M., Benjuya, N. & Burr, D. (2000). Do high impact exercises produce higher

A full body musculoskeletal model based on flexible multibody simulation approach utilised in bone strain anaflysis during human locomotion., *Comput Methods Biomech*

(2008). Analysis of dynamic strains in tibia during human locomotion based on flexible multibody approach integrated with magnetic resonance imaging technique,

multibody simulation approach in the analysis of tibial strain during walking, *Journal*

Targeted exercise against osteoporosis: A systematic review and meta-analysis for

osteoporosis and decision making criteria: 2009, *Salud publica de Mexico* 51 Suppl

vitro modeling of human tibial strains during exercise in micro-gravity, *Journal of*

Direction-specific diaphyseal geometry and mineral mass distribution of tibia and fibula: A pqct study of female athletes representing different exercise loading types,


URL: *http://libtest.csc.fi:8991/F/?func=direct&doc\_number=002760407&local\_base=fin01*


URL: *http://link.aip.org/link/?JDS/106/255/1*


**0**

**12**

*Brazil*

**Physical Insights Into Dynamic Similarity**

Comparative zoologists, evolutionary biologists, experimental biologists, and mechanical engineers make wonderful generalizations about the movements of different-sized running bipeds and quadrupeds (Heglund et al., 1974; McMahon, 1975; Alexander, 1976; Hoyt & Taylor, 1981; Garland, 1982; McMahon et al., 1987; Rome et al., 1988; Gatesy & Biewener, 1991; Farley et al.,1993; Cynthia & Farley, 1998; Bullimore & Burn, 2006), flying birds and swimming fish (Hill, 1950; Alexander, 2003; Taylor et al., 2003). The dynamic similarity across body mass and taxa of animals maintaining a certain gait in locomotion has been thoroughly investigated (Alexander, 1976, 1985, 1989, 2005; Alexander & Bennet-Clark, 1976; Alexander & Jayes, 1983; Marden & Allen, 2002; Biewener, 2005; Bejan & Marden, 2006a, b; Bullimore & Donelan, 2008). Based on integrative approach to animal locomotion (e.g., review by Dickinson et al., 2000) and using simple physical ideas (e.g., reviews by Lin, 1982 and Alexander, 2003), many cited above researches have demonstrated the importance of scaling biomechanics via reliably established scaling relations for gait characteristics with changes of speed and body mass. Among well known empirical findings of the dynamic similarity in animals observed across body mass are scaling relations established for stride (or stroke) speed and/or frequency (Hill, 1950; Heglund et al., 1974 Greenewalt, 1975; Garland, 1982; Heglund & Taylor, 1988; Gatesy & Biewener, 1991; Farley et al., 1993; Bullen & McKenzie, 2002; Bejan & Marden, 2006a), duty factor and relative stride length (Alexander & Jayes, 1983; Gatesy & Biewener, 1991), body force output (Alexander, 1985; Marden & Allen, 2002; Bejan & Marden, 2006a), limb stress and stiffness (Rubin & Lanyon, 1984; McMahon & Cheng, 1990; Farley & Gonzalez, 1996;

Although the concept of *mechanical similarity*, well known in analytical mechanics (e.g., Duncan, 1953; Landau & Lifshitz, 1976), was accurately re-formulated in application to the *dynamic similarity* (Alexander, 1976, 1989, 2003, 2005; Bullimore & Donelan, 2008), its exploration in biomechanics is often controversial. For example, Hill's seminal observation of dynamic similarity through the optimal speeds of birds in flight gaits (Hill, 1950) was found in sharp disagreement (McMahon, 1975) with that revealed through transient (trot-to-gallop) speeds in quadrupeds (Heglund et al., 1974). Then, comparing anatomic consequences of McMahon's elastic similarity (McMahon, 1975) for the stride length, Alexander noted

**1. Introduction**

Biewener, 2005; Bullimore & Burn, 2006).

**in Animal Locomotion. I. Theoretical**

**Principles and Concepts**

*Departamento de Física, Universidade Federal*

Valery B. Kokshenev

*de Minas Gerais, Belo Horizonte*


URL: *http://link.aip.org/link/?MSE/93/268/1*
