Nel seguente report vengono introdotte diverse campagne di simulazione volte all’omogeneizzazione del tessuto osseo trabecolare.
Background, analizzando la struttura dell’osso trabecolare è evidente come la forma dei pori e la struttura trabecolare influenza la rigidezza e la risposta strutturale. Viene analizzato come la variazione della geometria dei pori e la loro distribuzione influenza i parametri di omogenizzazione. Inoltre vengono fatte delle analisi caso specifiche.
Risultati, si ottiene una matrice di omogenizzazione e diverse analisi sulle variazione dei suoi parametri legate alle variazioni dei pori. Il parametro più influente è la frazione di pori ma può esserci anche un grande contributo dato dalla variazione del rapporto assiale;
Conclusioni, i risultati ottenutimostrano come il rimodellamento osseo contribuisce fortemente ad aumentare la rigidezza nella direzione di carico. Una matrice di omogenizzazione per considerare insieme entrambi i costituenti può essere usata ma con l’accortezza di rispettare la frazione volumetrica e verificare che il rapporto assiale sia vicino ad 1.


Il seguente report è estratto da un progetto del corso di Meccanica Computazionale dei Tessuti e Biomateriali, tenuto per Ingegneria Medica presso l’Università degli Studi di Roma Tor Vergata.


Indice

Introduzione

Lo scopo della seguente analisi è quello di indagare l’influenza delle proprietà microstrutturali sul comportamento macroscopico dell’osso trabecolare. Il flusso logico delle campagne di simulazioni è riportato in fig. 2.

Si parte dall’analizzare la struttura generica dell’osso trabecolare arrivando a scegliere un volume rappresentativo (RVE) idealizzato con una geometria tale da rispettare il rapporto tra i pori e la massa ossea trabecolare. Viene quindi discretizzato tale RVE e dunque il passaggio successivo prevede di ottimizzare la scelta della mesh in modo tale da avere un risultato sufficientemente accurato riducendo le richieste di sforzo computazionale. Si ottiene quindi una prima matrice di omogenizzazione che rispecchia le proprietà strutturali del RVE e
dell’osso trabecolare così idealizzato.

Vengono poi analizzate diverse proprietà geometriche e materiali seguendo le risposte dell’osso al carico,
fenomeno noto come adattamento e rimodellamento osseo. In questa campagna di simulazione viene analizzato come la variazione dei pori, della loro volume fraction e della loro distribuzione influenza il comportamento meccanico macroscopico.

Inoltre, vengono fatte ulteriori analisi andando a vedere come il rapporto assiale, la variazione dei
parametri materiali dei singoli elementi, quali matrice ossea e contenuto dei pori, influenzano la struttura. Vengono poi fatte ulteriori analisi andando a variare la finestra di interesse e considerando un caso paziente specifico.

Tessuto osseo

Il tessuto osseo è il tessuto biologico che forma le ossa. È caratterizzato da una notevole durezza e resistenza, dovute alla mineralizzazione della matrice extracellulare. È un tessuto in continua evoluzione, continuamente rinnovato e rimodellato per tutta la durata della vita. Costituisce l’impalcatura dello scheletro, da inserzione a muscoli e tendini e al suo interno contiene midollo osseo.

Il tessuto osseo è costituito da cellule e da una matrice extracellulare, organica ed inorganica. La parte organica contiene fibre collagene e sostanza amorfa formata da glicoproteine e proteoglicani.

Figure 1. Idealizzazione del RVE. (a) Ricostruzione tridimensionale; (b) Slice di riferimento no. 650; (c) ingrandimento del volume di interesse; (d) idealizzazione geometrica

Si possono distinguere due tipologie di osso: osso compatto e osso spugnoso [1]. L’osso compatto appare, all’esame macroscopico, come una singola massa solida. L’osso spugnoso invece ha un aspetto alveolare ed è costituito da sottili trabecole, un’insieme di piccole piastre ed aste che si anastomizzano in una rete tridimensionale. All’interno delle maglie che si formano con le spine trabecolari vi è accolto il midollo osseo. L’osso trabecolare costituisce il 20% della massa ossea complessiva ed è situato all’interno delle strutture ossee. La sua distribuzione varia fortemente da osso ad osso.

Osso corticale e trabecolare differiscono, oltre che per la loro architettura, anche per la velocità di sviluppo, funzione, vicinanza al midollo osseo, irrorazione sanguigna, rapidità nel turnover, fratture e risposta agli stimoli. Tutti questi fattori forniscono comportamenti materiali diversi. È bene sottolineare che l’osso, nel suo complesso, è costituito sia da una parte di osso compatto che una di osso spugnoso, quindi la risposta dell’intera struttura avrà entrambi i contributi.

Tipicamente il tessuto osseo spugnoso è meno resistente di quello compatto. In letteratura sono stati
condotti molti test per valutare le proprietà materiali dell’osso spugnoso ed esistonomolti risultati, alcuni anche controversi. Tipicamente vengo eseguite prove di: compressione, flessione, trazione uniassiale, misure ultrasoniche, calcolo a ritroso tramite simulazioni FEM 2D/3D e microindentazione [2]. Il range per ilmodulo elastico va da 076 GPa a 20 GPa. È ormai accettato che ilmodulo dell’osso trabecolare è il 20  30% inferiore a quello dell’osso compatto.

Inoltre, è noto in letteratura che il parametro fondamentale e più influente nel legare la microstruttura alla
risposta di rigidezza è la frazione volumetrica di pori.[3]

Ci sono però anche altri aspetti. Bisogna tenere conto che l’osso è un tessuto biologico e quindi soggetto al rimodellamento. In particolare è noto che l’osso risponde molto bene a stimoli di compressione andando ad aumentare la sua sezione trasversa ed inspessendo lo strato di osso compatto. Anche nell’osso spugnoso si vedono alcuni cambiamenti tra cui l’allineamento delle trabecole lungo le linee di forza. [4]

Struttura multiscala

In questa analisi viene presa come riferimento una vertebra lombare, L3. Le considerazioni sono comunque estensibili a qualunque altro osso trabecolare ma i parametri materiali fanno riferimento a dati selezionati per vertebre lombari. Le vertebre lombari sono solo alcune delle 33-34 vertebre che formano la colonna vertebrale umana. Il corpo vertebrale è formato principalmente da osso trabecolare con le singole trabecole allineate lungo le linee di forza. Questa disposizione fornisce un contributo importante per la distribuzione dei carichi e il sostegno del busto. Circa tre quarti dell’intero carico strutturale ricadono sulla colonna anteriore, ovvero sul corpo vertebrale, dischi intervertebrali e sulle placche terminali delle vertebre. La struttura ossea trabecolare può subire diversi cambiamenti con l’invecchiamento o a causa di alcune patologie [4] come la perdita di densità minerale, cambiamenti morfologici come un aumento della frazione volumetrica di pori, ispessimento o assottigliamento delle trabecole e perdita di connettività. Questi cambiamenti influenzano notevolmente la risposta meccanica del tessuto.

Osso trabecolare
Figure 2. Flusso logico della campagna di simulazione


L’osso è composto da diverse fasi. Contiene una fase organica, una inorganica e acqua. La fase organica è
composta da collagene e proteine mentre la fase inorganicada fosfato di calcio, simile all’idrossiapatite [5]. Alla nanoscala si trovano le molecole di collagene e i cristalli di idrossiapatite. Si passa poi ad una scala leggermente maggiore dove delle fibrille mineralizzate si uniscono a formare una singola lamella di spessore alcuni μm. Le singole lamelle si impacchettano a formare delle trabecole con uno spessore medio di 50μm e una lunghezza di circa 1mm. Siamo alla microscala. Si passa poi ad una scala intermedia dove possiamo considerare il volume di interesse di questa analisi dove troviamo sia trabecole che i pori riempiti di midollo fino ad arrivare alla macroscala dove si guarda all’osso nel suo insieme, le singole proprietà materiali e la struttura forniscono una risposta struttura ben definita.
L’obiettivo della seguente analisi è quello di portare in conto le differenze alla microscala, considerando un volume rappresentativo ben definito, e omogenizzare le proprietà meccaniche ottenendo un’unica matrice di rigidezza omogenizzata. Tale matrice omogenizzata descrive il comportamento materiale complessivo, tenendo conto sia delle singole proprietà materiali che della loro distribuzione spaziale e geometrica.

Omogenizzazione

L’obiettivo è quindi quello di ottenere delle proprietà materiali alla macroscala a partire dalle proprietà alla microscala come le rigidezza dei duemateriali principali, matrice ossea e contenuto dei pori.

Ad una scala sufficientemente grande rispetto le singole trabecole è possibile considerare i pori, ripieni di midollo osseo, e le trabecole ossee come due materiali diversi tra loro ma entrambi eterogenei e dal comportamento isotropo. L’isootropia è una prima approssimazione che però ci permette di procedere inmodo più efficiente nell’omogenizzazione andando a vedere come la geometria, la disposizione e la frazione volumetrica dei pori influenzi il comportamento macroscopico.

Con la procedura di omogenizzazione andremo a considerare una media delle proprietà meccaniche di osso e matrice ma questa media viene intesa in senso più ampio e il significato può variare a seconda dell’approccio. I contributi delle singole proprietà materiali vengono influenzati dalla geometria della struttura e dalle loro distribuzioni.

Volume di interesse

Il primo passo, per poter procedere ad un’omogenizzazione, è quello di definire un volume rappresentativo (RVE) tale da essere rappresentativo dal punto di vista statistico della distribuzione geometrica dei costituenti, quindi della microscala. Per farlo il codice calcola automaticamente il lato nel range di qualche mm a partire da alcune considerazioni. Viene rispettato il diametro medio dei pori [6] e la frazione volumetrica [2].

Sulla base di queste considerazioni vengono inseriti 25 pori equidistanti e tali che la loro distanza dal bordo sia lametà della distanza tra ogni poro. Questo consente di avere un RVE simmetrico e rappresentativo di una selezione delle trabecole. Maggiori informazioni sulla procedura di automazione del calcolo del RVE sono nella sezione 8.2.

Inoltre, per fare le opportune considerazioni nell’analizzare le tensioni medie, i pori vengono considerati pieni di un materiale a bassa rigidezza rappresentativo del midollo osseo. Questo serve sia a considerare un fenomeno più rappresentativo della realtà sia a permettere una corretta applicazione del teorema dellamedia delle tensioni.

Stime analitiche

Per confrontare i dati dell’omogenizzazione vengono considerate le stime analitiche di Voigt e Reuss. L’approssimazione di Voigt considera la matrice di rigidezza come unamedia pesata per le volume fraction delle rigidezze dei costituenti.

\begin{equation}
\mathbb{C}_{\text {Voigt }}=\left(1-v_{f}\right) \mathbb{C}_{b}+v_{f} \mathbb{C}_{m}
\end{equation}

Le stime di Reuss sono analoghe ma in cedevolezza:

\begin{equation}
\mathrm{S}_{\text {Voigt }}=\left(1-v_{f}\right) \mathbb{S}_{b}+v_{f} \mathbb{S}_{m}
\end{equation}

Vengono considerate queste due stime per via del teorema di Hill che garantisce, sotto le opportune ipotesi, che la matrice di rigidezza sia compresa tra esse [7]. Reuss costituisce il limite inferiore e Voigt il limite superiore dell’omogenizzazione.

I due materiali vengono considerati singolarmente isotropi, quindi tali da avere un comportamento materiale descritto dalla matrice di cedevolezza nell’equazione 3.

\begin{equation}
[\mathbb{S}]_{\text {iso }}=\left[\begin{array}{cccccc}
\frac{1}{E} & -\frac{v}{E} & -\frac{v}{E} & 0 & 0 & 0 \\
-\frac{v}{E} & \frac{1}{E} & -\frac{v}{E} & 0 & 0 & 0 \\
-\frac{v}{E} & -\frac{v}{E} & \frac{1}{E} & 0 & 0 & 0 \\
0 & 0 & 0 & \frac{2(1+v)}{E} & 0 & 0 \\
0 & 0 & 0 & 0 & \frac{2(1+v)}{E} & 0 \\
0 & 0 & 0 & 0 & 0 & \frac{2(1+v)}{E}
\end{array}\right]
\end{equation}

Questematrici vengono calcolate e per ogni iterazione del processo di omogenizzazione viene verificato che il risultato sia compreso in tali stime tramite la funzione VoigtReussTest[]. In particolare, passando la variabile verbOutput==True è possibile estrarne i valori:

\begin{equation}
\mathbb{C}_{\text {Voigt }}=\left(\begin{array}{ccc}
2019.340 & 865.404 & 0 \\
865.404 & 2019.340 & 0 \\
0 & 0 & 576.969
\end{array}\right)
\end{equation}
\begin{equation}
\mathbb{C}_{\text {Reuss }}=\mathbb{S}_{\text {Reuss }}^{-1}=\left(\begin{array}{ccc}
0.2262 & 0.0399 & 0 \\
0.0399 & 0.2262 & 0 \\
0 & 0 & 0.09316
\end{array}\right)
\end{equation}

Omogenizzazione numerica

Viene quindi strutturata una campagna di simulazioni tale da analizzare il RVE ed estrarre una matrice di rigidezza omogenizzata. Ulteriori considerazioni vengono poi fatte sul parametro C(om) 11 . L’approccio all’omogenizzazione segue condizioni al bordo omogenee di tipo spostamento. In particolare, viene applicato un campo di spostamento alla macroscala supponendo che sia rappresentativo di quello alla microscala. Quindi si ottiene un certo campo di tensione alla microscala, funzione del punto. Dunque è possibile considerare il campo di tensione alla macroscala come lamedia di quello alla microscala. In particolare, queste condizioni sono garantite dall’Average Strain Theorem. Quindi le tensioni medie misurate nel RVE sono definite come:

\begin{equation}
\langle\sigma\rangle=\frac{1}{\|\tilde{\Omega}\|} \int_{\Omega} \varepsilon d \tilde{\Omega}
\end{equation}

dove ˜Ω è l’area del RVE. Applicando un campo di deformazione alla macroscala si ottengono le deformazioni alla microscala da cui è possibile estrarre le tensioni semplicemente risolvendo il problema dell’equilibrio elastico:

\begin{equation}
\sigma=\langle\tilde{\sigma}\rangle=\mathbb{C}<\tilde{\varepsilon}\rangle=\mathbb{C} \varepsilon
\end{equation}

In particolare, applicando condizioni di deformazione unitaria su “1,”2 e 12 si ottengono rispettivamente prima, seconda e terza colonna di C. Quindi è possibile farlo in tre processi carico separati e considerare il contributo complessivo in quanto siamo in elasticità lineare sotto l’ipotesi di piccole deformazioni per cui vale la sovrapposizione degli effetti. La campagna di simulazione numerica viene risolta tramite il metodo agli elementi finiti utilizzando elementi triangolari a tre nodi. Questo viene implementato tramite il software Mathematica sfruttando il tool agli elementi finiti AceFEM. Viene sfruttata la procedura TEST[], scritta ad hoc, che permette di applicare in sequenza le tre condizioni di spostamento e ottenere così la matrice di rigidezza:

\begin{equation}
\mathbb{C}^{(o m)}=\left(\begin{array}{ccc}
756.78 & 80.90 & 0.00 \\
80.90 & 756.78 & 0.00 \\
0.00 & 0.00 & 40.87
\end{array}\right)
\end{equation}

Da questa matrice è possibile innanzitutto verificare che sono soddisfatti i limiti inferiore e superiore (Reuss e Voigt). Inoltre, si osserva un comportamento isotropo e come l’elevata pore fraction e quindi la grande presenza di materiale a bassa rigidezza va ad abbassare la rigidezza della matrice ossea. Ulteriori informazioni sull’automazione e sulla convergenza di questo risultato sono presenti nella sezione dell’analisi a convergenza §8.4.

NomeValoreReference
Eb3.8 GPaCowin [2]
νb0.3Dalstra et al. [8],Wirtz et al. [9]
Em15 KPaJansen et al. [10]
νm0.1
C(b)116730.77MPa§4.3 *
C(m)110.15Mpa§4.3 *
C(om)11756.78Mpa§4.3 *
L *1.324mm§4.1
Ta. 1. Parametri materiali (*Parametri calcolati)

Sebbene questo primo risultato di omogenizzazione possa sembrare sufficiente a considerare matrice ossea e pori come un’unica entità con tali proprietà materiali, è bene verificare cosa succede se le proprietà del singolo componente o della sua disposizione spaziale variano. I passaggi seguenti cercano di studiare se questo accade e come influenza il comportamento complessivo.

Variazione della microstruttura

Il carico della colonna vertebrale ricade per la maggior parte sul corpo vertebrale e sui dischi intervertebrali. Il corpo vertebrale è formato da osso trabecolare ad alta porosità con un guscio esterno altamente denso. Tuttavia, il guscio esterno è molto sottile ed è formato da un addensamento di osso trabecolare molto compatto, istologicamente differente dall’osso corticale. Diverse analisi agli elementi finito hanno stimato che il guscio esterno contribuisce per meno del 15% alla capacità di carico globale della vertebra [4, 11].

Inoltre, diversi campioni di osso trabecolare hanno mostrato proprietà meccaniche differenti a seconda della regione dove vengono prelevati [12, 13]. Questo è segno della risposta dell’osso che, essendo un tessuto vivo e in continuo rimodellamento, è variabile e si adatta all’ambiente. L’aumento delle proprietà di resistenza a compressione nella parte centrale dell’osso non è nient’altro che la risposta al maggior carico verticale trasmesso dalla parte polposa del disco vertebrale rispetto all’adicente regione periferica di anelli fibrosi.

Sulla base di queste considerazioni sono state condotte ulteriori analisi variando la geometria del RVE e osservando come varia la risposta materiale risultante. Le variazioni considerate seguono dall’idealizzazione di alcuni fenomeni fisiologici che portano al rimodellamento osseo indotto dall’ambientemeccanico in cui si trova immerso.

Variazione della simmetria

Le trabecole quindi rispondono al carico verticale posizionandosi prevalentemente lungo le linee di forza del carico e aumentano la loro densità e il loro spessore. Per indagare la variazione delle proprietà materiali è stato considerato un RVE con pori variabili. Un modo semplice, ma accurato, di analizzare questo comportamento è quello di considerare pori ellittici. Introducendo pori ellittici la risposta dell’osso al carico verticale può essere interpretata come una riduzione del rapporto assiale e quindi un restringimento dell’ellisse lungo la direzione x.

Fig. 3. Variazione della matrice di rigidezza al variare del rapporto assiale dei pori

Le campagne di simulazione condotte mostrano diversi cambiamenti nella matrice di rigidezza e quindi nelle proprietà meccaniche. Partendo da un rapporto assiale molto basso, quindi un’ellissemolto schiacciata, si osserva un grande aumento nel parametro C22 e un lieve aumento anche del parametro C11 rispetto al caso di riferimento (§4.3). Viene completamente meno la simmetria materiale e l’isotropia della matrice di rigidezza, segno che la risposta di rigidezza lungo la direzione del semidiametro maggiore è più grande di quella nella direzione del semidiametro minore. Questo rispecchia quanto produce l’osso con la sua risposta di rimodellamento.

Fig. 4. Relazione tra C11, C22 e AxR

Con un rapporto assiale di AxR = 0.1 il coefficiente C22 arriva a 5,376 GPa contro un C11 = 1.462 GPa. Questa differenza si riduce diminuendo il rapporto assiale. I due coefficienti convergono nel caso in cui il rapporto assiale tende a 1 ovvero nel caso in cui C11 = C22 = C(om) 11 . I limiti di Reuss e Voigt sono verificati in tutti i casi. I dati sono presenti in fig. 3 e fig. 4. Nel diminuire il rapporto assiale la riduzione dei coefficienti di rigidezza segue un andamento lineare. Il coefficiente C22 risente notevolmente di questo effetto mostrando una tasso di diminuzione molto più alto del coefficiente C11 che invece segue un andamento quasi costante.

Rimane evidente come al diminuire del rapporto assiale, e quindi anche della volume fraction, sia i limiti di Reuss e Voigt sia l’omogenizzazione numerica mostrino il coefficiente maggiore tendere al coefficiente C(b) 22 della sola matrice ossea pur rimanendo inferiore al suo 80%.

Lo stesso ragionamento può essere fatto nel caso contrario in cui malattie, come l’osteoporosi, portano ad una riduzione della massa ossea. Questo, rivisitato in ottica di un aumento del rapporto assiale, mostra come l’osso trabecolare tende a ridurre le sue proprietà strutturali anche in zone dove precedentemente il rimodellamento osseo potrebbe aver apportato un precedentemente miglioramento strutturale (meno pori e maggior resistenza). Cambiamenti morfologici tra cui l’assottigliamento delle trabecole, aumento dello spazio inter-trabecolare e la perdita di connettività posso portare ad un calo drastico della resistenza fino al rischio di frattura vertebrale [4].

Estensione dell’area di indagine

Un’ulteriore indagine è stata condotta andando ad analizzare se e come varia la risposta meccanica variando la finestra di interesse. Avendo idealizzato un RVE simmetrico con pori simmetrici disposti su una griglia regolare è chiaro che, aumentando semplicemente il lato, si andrebbero ad includere soltanto più pori. La volume fraction rimarrà la stessa così come la rispostameccanica. Anche se aumenta la somma delle tensioni queste vengonomediate su un area più grande e il loro rapporto, quindi la matrice di omogenizzazione, rimane costante. Lo studio si sposta allora su un RVE idealizzato in modo differente. I casi di interesse sonomostrati in fig. 5. L’analisi è divisa secondo due casistiche:

  • Analisi della rispostameccanica con una disposizione causale dei pori secondo una distribuzione precisa identificata da Doktor et al. [6].
  • Analisi della risposta meccanica con una distribuzione di pori ipotizzata a partire da un caso reale, paziente specifica [14].

Nell’analizzare i pori disposti secondo una distribuzione nota vengono considerate diverse inclusioni di raggio fisso r = {0.1, 0.15, 0.2, 0.03, 0.35}. Questi pori vengono disposti casualmente considerandone prima 100, fig. 5 (a), e poi 300, fig. 5 (b). Per l’implementazione computazionale maggiori informazioni in §8.3.4.

Si ottengono dei coefficienti della matrice di rigidezza C11, C22 superiori al caso §4.3. Questo in parte è dovuto ad una pore fraction minore, ovvero ad una presenza maggiore di matrice ossea che contribuisce ad aumentare la rigidezza complessiva. È interessante notare che il contributo della frazione volumetrica di pori si riflette in una maggior rigidezza lungo le direzioni principali anche se la struttura non sembra avere una direzione preferenziale per l’accrescimento osseo come in §4.3.

A rafforzare queste conclusioni sono anche i risultati della precedente analisi (§5.1) . Confrontandoli con questi dati si vede come il fattore dominante rimane la pore fraction ma anche come se aumenta la quantità di matrice ossea disposta lungo la direzione di carico aumenta la rigidezza in tale direzione.

Fig. 5. Variazione della regione di interesse. (a) introduzione di 100 pori con disposizione casuale; (b) introduzione di 300 pori con disposizione casuale; (c) introduzione di
28 pori su un RVE di lato 3 mm; (d) introduzione di 52 pori su un RVE di lato 5 mm.

Nell’aumentare la regione di interesse vengono considerati 300 pori e viene aumentato il lato del RVE. Questo permette di includere un numero maggiore di pori e la loro distribuzione è tale da ridurre leggermente i coefficienti della matrice di rigidezza. Si perde la simmetria esatta della matrice di rigidezza ma rimangono dei coefficienti C11 e C22 relativamente stabili. L’effetto è meno evidente nel caso di un maggior numero di inclusioni. La perdita delle proprietà di isotropia della matrice è molto evidente nei coefficienti che ci si aspetterebbe fossero nulli ma anche se molto piccoli, sono numeri diversi da zero. Questi ultimi risultati si discostano dall’omogenizzazione nel caso ideale e mostrano complessivamente una maggior rigidezza attribuibile all’aumento della massa ossea.

Questi risultati numerici sono influenzati dalla distribuzione spaziale dei pori e differenti disposizioni possono portare a matrici di rigidezza leggermente diverse pur rimanendo sotto i limiti della perdita di isotropia osservabile nel caso §5.1.

Analisi caso specifica

Un ulteriore analisi di interesse, simile al caso precedente, considera una disposizione di pori identificata a partire direttamente dal set di immagini Micro CT menzionate. Prendendo una slice di riferimento e selezionandone un determinato RVE è possibile estrarre una determinata struttura ossea fig. 5 (c) e (d). Sono stati considerati dei pori circolari disegnati manualmente tali da essere tangenti alla matrice ossea e da soddisfare la pore fraction compatibili con i dati disponibili in letteratura (circa 60%).

Il primo risultato evidente di questi test è l’aumento di rigidezza chiaramente dovuto alla minor pore fraction considerata, ovvero all’aumento di massa ossea. I coefficienti C11 e C22 pur essendo numericamente diversi rimangono molto simili ma più grandi dei casi precedenti ma inferiori ai casi limite riscontrati con la riduzione del rapporto assiale.

Fig. 6. Variazione di C11 al variare di Eb. Il coefficiente segue un andamento quasi
lineare con un incremento costante al crescere di Eb.
Fig. 7. Variazione di C11 al variare di Em, quest’ultimo rappresentato in scala
logaritmica. Il valore di C11 rimane quasi costante, si osservi che la scala è molto
grande.

Parametri materiali

Per il primo studio di omogenizzazione sono stati considerati dei parametri materiali medi tra quanto identificato in letteratura. Tuttavia è noto anche che il modulo elastico delle strutture ossee può variare in relazione all’età o a particolari patologie. Diversi risultati mostrano che persone giovani hanno moduli di rigidezza più elevati delle persone anziane [2].

È stato analizzato sia il contributo del modulo elastico della massa ossea, Eb, che dei pori, Em. In particolare, l’aumento del modulo elastico della massa ossea conferisce sempre più rigidezza e un incremento costante del parametri C11. Al variare di Eb da meno di 2 GPa a poco più di 6 GPa il coefficiente C11 passa da poco meno di 300 a 1000. Nel mezzo si trova il caso affrontato nel primo studio (4.3). Può anche essere fatto un fitting polinomiale soddisfatto dalla relazione:

\begin{equation}
\mathbb{C}_{11}=\mathbb{C}_{11}\left(E_{b}\right)=0.184+0.151 \cdot E_{b}
\end{equation}

Al contrario, variando ilmodulo dei pori, Em, non si osservano grandi cambiamenti. Variando da 0.15 MPa fino a 1000 MPa si osserva una variazione inferiore all’1 %. Questo rispecchia il fatto che il maggior contributo di rigidezza è dato dalla matrice ossea. Allo stesso tempo ci dice che, a livello numerico, è possibile utilizzare anche unmodulo leggermente più alto di quanto ipotizzato in modo da evitare numeri troppo piccoli che possono dare problemi computazionali andando oltre la precisione numerica del calcolatore.

Conclusioni

Alla luce dei risultati ottenuti è chiaro che l’elemento predominante nell’influenzare la rigidezza complessiva è la frazione di pori. Nel fare l’omogenizzazione si può tenere conto della frazione di pori e ottenere un risultato numerico attendibile.

È bene però considerare che c’è una forte influenza anche del rapporto assiale nei pori. Essendo l’osso una struttura viva e in continuo rimodellamento è bene tenere conto anche delle direzioni preferenziali di carico lungo le quali si potrebbe osservare una riduzione del rapporto assiale. Trascurare questo effetto e considerare l’osso trabecolare come isotropo potrebbe portare ad una matrice di omogenizzazione con notevoli errori numerici.

Un forte peso è dato anche daimoduli di rigidezza dei singoli costituenti ma il peso maggiore è dato dalla matrice ossea che si occupa di conferire rigidezza alla struttura globale. È possibile evitare uno studio caso specifico e considerare parametrimedi e l’errore rimane in dei limiti accettabili. In alcuni casi però questo potrebbe sottostimare la rigidezza.

Possibili sviluppi futuri possono essere quelli di analizzare l’influenza di ulteriori parametri come lo spessore trabecolare o la connettività. Inoltre, si potrebbe passare ad unmodello di osso trabecolare più avanzato considerando l’interno dei pori non più come un materiale fortemente elastico ma ripieno di fluido con possibilità di scorrere e andando a considerare unmodello poroelastico.

Reference

  1. Monesi V. Istologia. Piccin; 2012.
  2. Cowin S. Bone Mechanics Handbook. CRC Press; 2001.
  3. Nazarian vSDZDMRSBD A. Bone volume fraction explains the variation in strength and stiffness of cancellous bone affected bymetastatic cancer and osteoporosis. Calcified tissue international 007;83:368– 379.
  4. Ferguson SJ, Steffen T. Biomechanics of the aging spine. European spine journal 2003;12:S97–S103.
  5. Hamed E, Jasiuk I, Yoo A. Multi-scale modelling of elastic moduli of trabecular bone. Soc Interface
    2012;9.
  6. Doktor T, Valach J, Kytyr D, Jiroušek O. Pore Size Distribution of Human Trabecular Bone – Comparison of IntrusionMeasurements with Image Analysis p. 115–118.
  7. Aboudi J, Arnold SM, Bednarcyk BA. Chapter 3 Fundamentals of the Mechanics of Multiphase Materials. In: Aboudi J, Arnold SM, Bednarcyk BA, editors. Micromechanics of Composite Materials Oxford: Butterworth-Heinemann; 2013.p. 87–145. https://www.sciencedirect.com/science/ article/pii/B9780123970350000033.
  8. DalstraM, Huiskes R, Odgaard A, Van Erning L. Mechanical and textural properties of pelvic trabecular
    bone. Journal of biomechanics 1993;26:523–35. https://doi.org/10.1016/0021-9290(93)90014-6.
  9. Wirtz DC, Schiffers N, Pandorf T, Radermacher K, Weichert D, Forst R. Critical evaluation of known bone material properties to realize anisotropic FEsimulation of the proximal femur. Journal of
    biomechanics 2000;33:1325–30. https://doi.org/10.1016/s0021-9290(00)00069-5.
  10. Jansen LE, Birch NP, Schiffman JD, Crosby AJ, Peyton SR. Mechanics of intact bone marrow. Journal
    of the mechanical behavior of biomedical materials 2015;50:299–307. https://doi.org/10.1016/j.
    jmbbm.2015.06.023.
  11. Silva MJ, Keaveny TM, HayesWC. Load sharing between the shell and centrum in the lumbar vertebral body. Spine 1997;22,2:140–50.
  12. Keller TS, Hansson TH, Abram AC, Spengler DM,PanjabiM. Regional variations in the compressive properties of lumbar vertebral trabeculae. Effects of disc degeneration. Spine 1989;14:1012–1019.
  13. Keller TS, Ziv ME I, Spengler D. Interdependence of lumbar disc and subdiscal bone properties: a report of the normal and degenerated spine. Journal of spinal disorders 1993;6:106–113.
  14. Beller G, Burkhart M, Felsenberg D, Gowin W, Hege HC, Koller B, et al. Vertebral Body Data Set ESA29-99-L3 2005;http://bone3d.zib.de/data/2005/ESA29-99-L3/.
  15. Mastrofini A. Image processing volumetrico su osso trabecolare. website 2021 ; https://alessandromastrofini.it/image-volumetric-bone/.
  16. Cowin S. Integrated Bone Tissue Physiology. Bone Mechanics Handbook 2001;1:4–5.
  17. Mathworks. Find circles using circular Hough transform 2012;https://it.mathworks.com/help/images/
    ref/imfindcircles.html/.


Scarica il report

Omogeneizzazione numerica di materiali microstrutturati

Accedi alla repository