Symulacja numeryczna procesów adaptacyjnych występujących w kości piszczelowej pod wpływem zmian obciążenia cd.

2. MATERIAŁ I METODA

Symulację procesów adaptacyjnych przeprowadzono przy użyciu modeli MES (Metoda Elementów Skończonych). Niezbędne było odtworzenie doświadczenia Cartera. W tym celu stworzono dwuwymiarowy model kości udowej o geometrii zgodnej z Cartera. Ponieważ był to model płaski w celu zasymulowania trójwymiarowej struktury kości dodano side plate. Side plate jest równoległą do modelu strukturą zaczynającą się od podstawy modelu kości udowej, biegnącą przez cała długość trzonu kości i kończącą się pod krętarzem większym kości udowej. Side plate ma geometrię identyczną z modelem właściwym, również podział na siatkę elementów skończonych jest taki sam. Obie struktury (side plate i model właściwy) są ze sobą połączone trzema najbardziej skrajnymi węzłami zarówno po stronie bocznej jak i po stronie przyśrodkowej. Side plate ma różną grubość na trzech poziomach, są to następujące grubości: 2 mm, 1 mm oraz 0.5 mm. Na model właściwy kości została zadana grubość 4 mm. Przy dzieleniu side plate na siatkę elementów skończonych użyto trójwymiarowych elementów typu belkowego z dwoma stopniami swobody w każdym węźle. Natomiast podczas dzielenia na siatkę elementów skończonych modelu kości udowej użyto trójwymiarowych elementów typu powłokowego z dwoma stopniami swobody w każdym węźle. Liczba elementów całej struktury (model właściwy kości udowej +side plate) wynosi 2170, a samego modelu kości 1717. Bazując na zdjęciach rentgenowskich kości piszczelowej, stworzono model kości. Następnie analogicznie do kości udowej utworzono side plate symulujące trójwymiarową strukturę kości. W przypadku kości piszczelowej side plate rozpoczyna się u podstawy modelu, biegnie wzdłuż długości trzonu kości i kończy się pod kłykciami: bocznym i przyśrodkowym. Również w tym przypadku geometria oraz podział na siatkę elementów skończonych struktury side plate, są identyczne jak w modelu kości piszczelowej. Obie struktury (side plate i model właściwy) tak jak w kości udowej są połączone ze sobą trzema najbardziej skrajnymi węzłami zarówno po stronie bocznej jak i po stronie przyśrodkowej. Również w tym przypadku side plate ma różną grubość na trzech poziomach, są to następujące grubości: 2 mm, 1 mm oraz 0.5 mm, a model właściwy kości ma grubość 4 mm. Przy dzieleniu side plate na siatkę elementów skończonych użyto trójwymiarowych elementów typu belkowego z dwoma stopniami swobody w każdym węźle. Natomiast, podczas dzielenia na siatkę elementów skończonych modelu kości piszczelowej, użyto elementów typu powłokowego z dwoma stopniami swobody w każdym węźle. Liczba elementów całej struktury (model właściwy kości piszczelowej + side plate) wynosi 2370, a samego modelu kości 1278.

Początkowe właściwości były jednakowe dla modelu kości udowej i kości piszczelowej. Przyjęte wartości: gęstość tkanki kostnej 0,5 g/cm2, moduł Younga 356 MPa oraz współczynnik Poissona 0,2. Analizując zmiany własności mechanicznych w tkance kostnej pod wpływem obciążeń o charakterze zmiennym (chód) obciążono model siłami odpowiadającymi poszczególnym fazom jednego cyklu chodu – od momentu styku pięty z podłożem do momentu oderwania od podłoża palców stopy.







Na tak przygotowanych modelach zasymulowano przebudowę kości wg następującego algorytmu: 1.Obliczenie energii odkształcenia oraz gęstości energii odkształcenia.


gdzie: Uo – wartość energii odkształcenia
ρ - wartość gęstości tkanki kostnej

2.Obliczenia wartości naprężeń „efektywnych” dla każdego elementu skończonego w modelu wg zależności podanej przez Cartera [4]:


gdzie : E- wartość modułu Younga dla danego elementu
U – wartość gęstości energii odkształcenia.

3.Obliczenia rozkładu bodźca mechanicznego przebudowy wg zależności [4]:


gdzie: i – kolejny przypadek obciążenia
n – ilość cykli dla i-tego przypadku obciążenia
m – współczynnik wyznaczany doświadczalnie.

Wg Cartera istnieje tzw. „strefa martwa”, w której nieznaczne zmiany obciążenia w stosunku do obciążenia fizjologicznego nie powodują przebudowy struktury kości. Tutaj zmodyfikowany algorytm Cartera poprzez uwzględnienie dodatkowego parametru – gęstości kości.

 4. Obliczenie współczynnika absorpcji/resorpcji kości.

gdzie: - stała wartość bodźca
- wartość bodźca mechanicznego w danej iteracji
c – stała wartość (współczynnik nachylenia krzywej)
w – szerokość przedziału określającego „strefę martwą”

5. Obliczenie gęstości kości w zależności od jej powierzchni. 

gdzie: BS/TV - powierzchnia kości
r' - wartość średnia liniowego współczynnika absorpcji/resorpcji kości
6.Obliczenie wartości zmiany gęstości tkanki kostnej dla każdego elementu skończonego.
7.Obliczenie wartości zmian modułu sprężystości dla każdego elementu skończonego
w modelu wg zależności:


gdzie: ρ - wartość gęstości tkanki kostnej
A, B – współczynniki wyznaczane doświadczalnie.
8.Zmiana wartości modułu sprężystości w każdym elemencie skończonym modelu na wartość uzyskaną w wyniku obliczeń.
9.Zakończenie obróbki modelu i przejście do kolejnego kroku iteracji.
Modele wykonano w programie Ansys, a procedury symulacyjne stworzono programując w języku APDL.