Bone as a microcontinuum 1 2 Josef Rosenberg & Robert Cimrman & Ludˇ ek Hynˇ cík University of West Bohemia in Plzeˇ n Department of Mechanics & New Technology Research Centre Univerzitní 22, 301 14 Plzeˇ n How to treat the microstructure? • homogenization • theory of mixtures, of composites • microcontinuum theories 1 Presentation for the conference Výpoˇ ctová Mechanika 2001 , Neˇ ctiny, 29.-31. October 2001. 2 Typeset by ConT EXt ( http://www.pragma-ade.nl ). 1
Basic kinematics • Continuum "points" can translate, but also rotate and deform → micromorphic continuum. Position within a particle given by x ′ = x + ξ , y ′ = y + η . • • Special types: − microstretch continuum: rotation + volume change, − micropolar continuum: rotation only. Figure 1 Coordinates within particles. 2
General balance equations The balance of forces and balance of stress moments equations: ,k + ρf l = 0 , ,k + t ml − s ml + ρl lm = 0 . t kl m klm (1) stress tensor in a particle, t ′ kl = t ′ lk , t ′ kl . . . s lm . . . micro-stress average — stress tensor of the macrovolume averaged across the volume (symmetric), t kl . . . stress tensor of the macrovolume averaged across the surface (non-symmetric), m klm . . . the first stress moment — moment of the forces acting on the surface of the macrovolume with respect to its centre of gravity, l lm . . . the first body moment of the volume forces with respect to the centre of gravity of the macrovolume, f l . . . averaged volume force. � � k d s ′ = t kl n k d S , k d s ′ = m klm n k d S . t ′ kl n ′ ξ ′ m t ′ kl n ′ Some defining relations: d S d S 3
Special types • microstretch continuum — 7 degrees of freedom m klm = 1 3 m k δ lm − 1 2 e lmr m k r , (2) l kl = 1 3 lδ lm − 1 2 e klr l r . (3) • micropolar continuum — 6 degrees of freedom m k = 0 , l = 0 . (4) 4
Micropolar continuum - the boundary value problem • Basic equations: ,k + ρf l = 0 , l,k + e lmn t mn + ρl l = 0 , t kl m k (5) t kl = ρ ∂ Ψ ∂y k ∂ Ψ ∂y k ∂x K χ lL , m kl = ρ 0 ∂x K χ lL , ∂ Γ LK ∂ Ψ KL k = ∂η l k = ∂ξ l Γ KL = 1 Ψ KL = y k MN χ kM χ kN , where χ l ∂ξ k , χ l ,K χ kL , 2 e ∂η k . K For the isotropic continuum holds (denoting γ ij = φ i,j , ε kl = ∂u l ∂x k + e lkm φ m ): • t kl = λε m m kl = αγ m m δ kl + ( µ + κ ) ε kl + µε lk , m δ kl + βγ kl + γγ lk . (6) � t kl n k = ˆ � The boundary conditions: u k = ˆ u k t l • on ∂ Ω 1 , on ∂ Ω 2 . m kl n k = ˆ φ k = ˆ φ k m l 5
Variational formulation The solution is the stationary point of the potential (see [8] ) � � m + ( µ + κ ) ε kl + µε kl � Π( u, φ ) = 1 λδ kl ε m ε kl d x 2 Ω � � m + βγ kl + γγ lk � � + 1 αδ kl γ m u i n j + g ij ) τ ij d x γ lk d x + (ˆ 2 Ω ∂ Ω 2 � � � � (ˆ ρ ˆ ρ ˆ φ k n l + γ kl ) m kl d x − l l φ l d x + f i u i d x − τ i u i d x − ˆ ∂ Ω 2 Ω ∂ Ω 1 Ω with the constraints ε kl = ∂u l γ kl = ∂φ k − φ i u j = γ ij on ∂ Ω 2 . ∂x k + e lkm φ m , − u i n j = g ij on ∂ Ω 2 , ∂x l , The weak solution of the problem at page 5 satisfies (we omit loading terms here) � Π( u, φ ; δu ) = 0 → τ kl δ u ε kl d Ω = 0 , (7) Ω � Π( u, φ ; δφ ) = 0 → ( τ kl δ φ ε kl + m kl δ φ φ l,k ) d Ω = 0 . (8) Ω 6
FE discretization Denote: 1 ≡ [1 , 1 , 1 | 0 , 0 , 0 | 0 , 0 , 0] T , J . . . a permutation matrix, G , ν strain operators. t e = ( λ 11 T + ( µ + κ ) I + µ J ) [ G + | ν ] d e = D 1 Bd e , m e = ( α 11 T + β J + γ I ) G + φ e = D 2 G + φ e . � �� � � �� � D 1 D 2 Discrete balance equations for one element: � � � � � | ξ q = � | ξ q · d e = [ A e , B e ] d e = 0 , G + T t e J 0 W G + T D 1 B J 0 W U e ≡ ( ← Eq. 7 ) q q � � � � ( ν T t e + G + T m e ) J 0 W � | ξ q = � ν T D 1 B J 0 W | ξ q · d e φ e ≡ q q � � | ξ q · φ e = [ C e , D e ] d e + E e φ e = 0 . + � G + T D 2 G + J 0 W ( ← Eq. 8 ) q ⇒ Linear system with indefinite matrix: � A e � f e � � u e � � B e = . (9) φ e g e C e D e + E e 7
Analytical verification I The analytical solution is known in some cases (cf. [3] , results taken from [8] ) , e.g.: • a plane with a hole loaded in tension, • compute the stress concentration factor on the boundary of the hole. mesh microrotations Figure 2 Plane with a hole. 8
Analytical verification II R = radius of the hole (macroscopic characteristic length) [ m ] c [ m ] = characteristic length of the microstructure K = stress concentration factor Theory: • linear elasticity: red curve ( K = 3 ) • micropolar elasticity: green curve Numerical values: • linear elasticity: magenta curve • micropolar elasticity: blue curve • adjusted (shifted by LE numeric − LE theory): cyan curve Stress concentration( R/c ). Figure 3 9
Femur bone with nail — motivation Figure 4 Example of a fixation of a bone. 10
Femur bone with nail — material data λ [Pa] µ [Pa] κ [Pa] α [N] β [N] γ [N] set 1 . 8 · 10 10 − 1 . 468 · 10 10 3 . 837 · 10 10 − 120 120 240 MP1 1 . 8 · 10 10 − 1 . 468 · 10 10 3 . 837 · 10 10 − 12000 12000 24000 MP2 1 . 8 · 10 10 4 . 5 · 10 9 LE — — — — Table 1 Material data. • Equivalent LE set was obtained using λ E = λ M , µ E = µ M + κ/ 2 ( → E = 1 . 26 · 10 10 [Pa], ν = 0 . 4 ). Material data of the steel nail: E = 2 . 1 · 10 11 [Pa], ν = 0 . 3 . • • Characteristic lengths of the microstructure: − MP1: c = 0 . 1283 [mm] − MP2: c = 1 . 283 [mm] • Characteristic length of the macrostructure = radius of the hole. • LE set was used in PAM-Crash code for verification of our solver — the results are denoted as "PC". 11
Femur bone with nail — loads • Two kinds of loading: bending and torsion. • Observed micropolar effect: decrease of stress on the femur–nail interface bending torsion Figure 5 Original (white) + deformed femur mesh (magnified displace- ments), LE set used for the bone. 12
Femur bone with nail — evaluation lines The nail was considered to be fixed to the bone—nomovementbetweenthetwoma- terials was allowed. The stress was evalu- ated along these lines on the surface of the hole drilled into the bone: t 22 [kPa], torsion case. Figure 6 13
Femur bone with nail — stress along the lines • Bending load: different behaviour (tension-compression) of middle and "non- middle" rows of elements ⇒ separate plots. • Torsion load: no such phenomenon. t 33 [kPa] (bending) t 22 [kPa] (torsion) Figure 7 Stress along the lines, MP2 set used for the bone. 14
Femur bone with nail — example Ia • We plot "averaged" stress along the front and back lines of Figure at page 13 . • The "averaging" = the least squares fitting of stress in the elements of Figure 7 ). middle element row upper element row t 33 along the lines, bending. Figure 8 15
Femur bone with nail — example Ib • The bending case — fitting with the second order polynomial. • The torsion case — fitting with the third order polynomial. t 22 along the lines, torsion. Bending, middle element row, MP1 set. Figure 9 Averaging example + torsion case results. 16
Femur bone with nail — example IIa • Dependence of stress on c : l t varied in range � 0 . 2 , 2 � [mm] while keeping the other parameters constant. This resulted in c variation in range � 0 . 1283 , 1 . 283 � [mm]. • Stress was evaluated in 6 selected elements (“left” end of the hole (the lowest x coordinate), see Figure 7 , Table 2 . • Note the difference between middle and non-middle elements in the bending case. element 5786 4236 4351 6103 6050 6123 line front front front back back back row upper middle lower upper middle lower Table 2 Selected elements. 17
Femur bone with nail — example IIb t 33 ( c ) , bending t 22 ( c ) , torsion Dependence on c in the selected elements. Figure 10 18
Femur bone with nail — example IIc t 33 ( c ) , bending t 22 ( c ) , torsion Dependence on c in element 4236. Figure 11 19
Conclusion • Linear micropolar elasticity was introduced. • Presented examples showed a strong influence of the microstructural parameters on the stress. • Further work: − micropolar anisotropic continuum − micromorphic continuum − material parameter identification 20
Recommend
More recommend