The model shown in Fig. (1) was developed using MSC MARC^® 2003 software.

Fig. (1) Set-up for FE pull out test (a) and alternative loading condition (b).

The bone consists of cortical bone (E =11 GPa, ν=0.33, 5 mm in thickness) and trabecular bone. Two different trabecular bone densities were simulated (Table 1), and their respective mechanical properties were obtained from the relations [14Doblarè M, Garcìa JM. Anisotropic bone remodelling model based on a continuum damage-repair theory J Biomech 2002; 35: 1-17., 16Fyhrie DP, Vashishth D. Bone stiffness predicts strength similarly for human vertebral cancellous bone in compression and for cortical bone in tension Bone 2000; 26: 169-73.]:

Table 1

Mechanical Properties of Trabecular Bone as a Function of Bone Density

where:

ρ is bone density (g/cm³);

E is Young’s modulus (MPa);

σ_U is the ultimate tensile stress (MPa).

The bone density values were chosen to simulate recurrent clinical situations (adult patient with mild bone resorbing) and significantly different Young’s modulus (1:2 ratio).

With regard to the screws, four different constitutive materials were simulated in order to study both currently used materials such as stainless steel (screw models n. 9-10) and titanium Ti6Al4V (screw models n. 3-8), and more innovative solutions such as low stiffness titanium (Ti12Mo5Ta, screw model n.2) and PMMA reinforced by an inner Ti12Mo5T cylinder (screw model n.1); the last two materials have been chosen in order to better approach trabecular bone’s Young modulus.

The general geometry of the simulated screws and the mesh details are shown in Fig. (2) (symbols as in Table 2); they were generated from an effectively produced geometry (screw model n. 9), varying all parameters, one by one.

Table 2

Geometric Parameters of Screws

Fig. (2) Screw parameters and detail of mesh.

All the screw threads were simulated with symmetrical flanks, fillet radii were all equal, and the pitch was constant along the screw axis. The axisymmetric model entails shorter calculation times.

The constraints and boundary conditions are shown in Fig. (1): they correspond to the pull out test (Fig. 1a) or to an alternative condition, with a different position of the constraint (Fig. 1b); this alternative condition should represent a more realistic working condition. According to Saint Venant’s, the second condition is equivalent to the first one, however the small displacement hypothesis cannot be assumed a priori in this case, due to the low stiffness of trabecular bone. For both configurations, the immediate post-operative condition and that of complete osseointegration were simulated: the first condition was created through a ‘touch’ contact between the screw and trabecular bone (whatever displacement is allowed with the exception of penetration), while the second was obtained through a ‘glue’ contact (the two components are completely bounded to each other).

The mean number of elements per model was 7500.

The maximum allowed force was calculated by iteration (non-linear case), specifying that the maximum stress on the bone must not exceed its ultimate tensile stress. This hypothesis is quite restrictive: actually bone is not a fragile material, tension peaks produce localized plastic deformations and stresses are redistributed over a wider area. However bone is likely to be stressed by dynamic loads and fatigue damage may occur when stress peaks produce a crack which progressively propagates during loading, until structural failure occurs.

The finite element model with ‘touch’ boundary condition was validated by means of experimental tests. All tests were performed on an hydraulic INSTRON 8872 test machine, equipped with a ±5 kN load cell.

Pull out tests [19Cheung LK, Zhang Q. Stability consideration for internal maxillary distractors J Craniomaxillofac Surg 2003; 31: 142-8.] were conducted replacing bone by polyurethane foam, as reported in the literature [20Palissery V, Taylor M, Browne M. Fatigue characterization of a polymer foam to use as a cancellous bone analog material in the assessment of orthopaedic devices J Mater Sci Mater Med 2004; 15: 61-7.].

The mechanical properties of polyurethane foam were determined by ten compression and ten tensile static tests, where a linearly increasing displacement was applied. The nominal stress was calculated from the force, divided by the initial specimen section; the nominal strain was calculated from the current displacement, divided by the initial specimen height.

Young’s modulus was determined through five dynamic tests, applying a pulsating sinusoidal load (Fig. 3, left).

Fig. (3) Curves of compression (upper) and tensile (lower) tests on polyurethane foam: nominal stress (σnom) versus nominal strain (εnom).

Young’s modulus was calculated as 6 MPa (s.d. = ± 1.21 MPa), while ultimate tensile stress was 0.42 MPa (s.d. = ± 0.084 MPa); these data were input into the numeric model to be validated. Compression tests produced a progressive packing of the foam (reaching up to 60% height reduction) which does not break but shows buckling.

The pull out tests required specific equipment: two holes drilled in the opposite sides of a rectangular metal profile (obtained from an extruded bar) and a polyurethane block was placed within the box profile (Fig. 4).

Fig. (4) Set-up of pull out test and load vs displacement curve.

A bolt screwed into the lower hole was then held by the lower jaw of the machine (the bolt was left free to move in the cross-sectional plane, in order to avoid bending moments acting on the screw). A screw was inserted through the upper hole and inserted into the polyurethane foam. The metal profile simulated the cortical layer of the FE model, and prevented the polyurethane foam deforming when the screw was pulled. A preparatory hole (6.5 mm in diameter) was drilled in the foam before inserting the screw, which was always set at the same height. The tests were conducted on a standard metrical screw (M10, UNI 4536) whose geometry was known in detail, allowing an accurate FE model to be constructed. The lower end of the shaft was threaded with 17 threads.

The screw was pulled at a speed of 2 mm/min. The pull out test was repeated 23 times and the mean force determined was 120 N (s.d. = ±14 N). A typical force versus displacement curve is shown in Fig. (4).

Ten different screw models were realized (Table 2). Each model was implemented with four different pitches (1.4 mm, 2.4 mm, 3.4 mm, 4.4 mm) and two different values of bone density (0.35 g/cm³; 0.47 g/cm³). Overall, each screw model was analyzed in eight FE models. The numerical results were summarized through a multifactor model [22Milton JS, Arnold JC. Introduction to Probability and Statistics. New York: McGraw-Hill 1995; pp. 604-54.], examining both the effects of each single factor and the interaction between the pitch and bone density.

where:

Y = pull-out limit load; d = bone density; p = pitch μ = mean

β = linear effect of bone density

γ_L = linear effect of the pitch; γ_Q = square effect of the pitch

γ_C = cube effect of the pitch

βγ_L = effect of interaction between the pitch and bone density

βγ_Q = effect of interaction between the bone density and the square effect of the pitch

βγ_C =effect of interaction between the bone density and the cube effect of the pitch

Up to the third power of p could be considered because four different pitches were simulated; the only linear effect of bone density could be considered because only two density values have been simulated.

Eq. (3) can be reformulated in order to isolate the contribution of the screw diameter:

or in order to isolate the contribution of the screw pitch:

The coefficients of this numerical model were calculated for every screw, solving a mathematical system with eight parameters (the effects) and eight equations (pull out forces, calculated from the FE models).

In order to validate FE models, it is necessary to select a failure criterion. The experimental tests showed that the polyurethane foam exhibited fragile behavior in the tensile tests, while compressive tests were much less critical. The maximum (positive) main stress failure criterion was chosen with an ultimate tensile stress of 0.42 MPa.

As shown in Fig. (5), the experimentally evaluated pull out force (120 N) generates ultimate tensile stresses at the last thread, close to the fillet. The FE model is therefore validated because it correctly predicts a localized failure when the pull out force is applied.

Fig. (5) Maximum principal stresses (left and center) and equivalent Von Mises stresses (right).

Actually this validation is quite specific: it has been demostrated that FE model correctly predicts the failure load which is the object of this paper; a more detailed validation would be required if other information were sought from the finite element model (stiffness behavior, stress distribution, etc.).

The location of ultimate tensile stress at the last thread may be unexpected, given that in the literature, the first thread is considered the most critical [21Shigley JE, Mischke CR, Budynas RG. Mechanical Engineering Design. Seventh ed.. New York: McGraw-Hill 2004; pp. 345-92.]. However, the pull out test is a static test, while the above considerations refer to fatigue life. Further, the mechanical behavior of threaded joints using polyurethane foam (or bone) differs from the behavior of metal joints: both foam and bone are porous, leading to completely different modalities of failure for tension and compression. The first “female” thread is subjected to compression stresses (contact between the screw and the bone), while the last is subjected to tensile stresses (as evidenced in Fig. 5).

The results of the multifactor analysis are illustrated in Fig. (6), where the trends of the ‘f(p)’ function (eq. 5) vs pitch are plotted for different screw models, bone densities and boundary conditions (post-operative or complete osseointegration). Table 3 shows the values of the D coefficient (eq. 4) reported for different screw models.

Table 3

D[10^-3Nm³/kg] from Eq. (5)

Fig. (6) fp=p⋅γL+βγL⋅d+p2⋅γQ+βγQ⋅d+p3⋅γC+βγC⋅dN

Fig. (7) Load distribution between screw and “female” thread for different materials.

Some considerations can be made:

• In general, the differences in the stress distribution between pull out and the alternative constraint condition (with all other conditions constant) are minimal: the position of the full constraint does not influence the stress pattern on the “female” thread. In this sense, the pull out test is representative of a wide range of operating conditions.
• The post-operative and full osseointegration conditions are very different: the latter is more sensitive to bone density (in Table 3 the D coefficient values calculated for the post-operative condition are lower than that calculated for complete osseointegration) and to thread pitch (in Fig. (6), curves B, D, F are steeper than curves A, C, E). This can be explained considering that, in full osseointegration, the screw and the bone are forced to undergo the same deformation. The post-operative condition is more critical because a lower pull out force is required. It is therefore advisable to construct non-linear numerical models where this situation is analyzed, implementing a ‘touch’ contact between the screw and bone.
• The coefficients of eq. 5 vary from model to model.
• The influence of pitch varies more for high values of bone density (Fig. 6). Effectively, greater bone stiffness implies a greater sensitivity to thread geometry because stresses are more localized.
• For a given screw, the optimal pitch (that providing maximum pull out force) is not always the maximum simulated pitch. For example, in the pull out test regarding the post-operative condition, given the screw model N. 8, the optimal pitch is between 3.4 and 4.4 mm for a bone density of 0.35 g/cm³, while the optimal pitch is between 2.4 mm and 3.4 mm for a bone density of 0.47 g/cm³.
• With regard to the interaction between the pitch and bone density, coefficients βγ_L, βγ_Q and βγ_C can assume either positive and negative values (for example, βγ_L>0 in the screw model 1, while βγ_L has a negative value in model 2). This means that in some models, the effect of bone density on pull out force is smaller for greater pitches. This result further confirms the conclusion that the optimal screw geometry will be different for different bone densities.
• Many screw models show a similar pattern of function p vs pitch (Fig. 6).
• There is a difference between load distribution in mechanical threaded joint and biomechanical threaded joint, as demonstrated by further numerical models, which implemented different values of “r”, the parameter which defines the ratio between the Young’s modulus of bone and that of the screw (Fig. 7). In the biomechanical case, r is very low because the Young’s modulus of bone is approx. 300 MPa, while the Young’s modulus of the screw is about one thousand times greater. Conversely, they have similar values in the mechanical and consequently a much more homogeneous load distribution is obtained on the second, third and fourth threads the biomechanical case.
• According to our results, screw model n. 10 is generally to be preferred; the performance of this model is moderately influenced by the thread pitch; considering the post-operative condition which is the most critical, a pitch larger than 1.4 mm should be recommended for the higher bone density.