Analisi computazionale dello stress meccanico nella diverticolosi del colon

Modello del colon

Il software ad elementi finiti Abaqus (versione 6.13, SIMULIA, Providence, RI) è stato utilizzato per sviluppare i modelli e condurre le simulazioni associate. Tutte le simulazioni sono state eseguite su una Dell Workstation T7810 con due processori Intel Xeon e 32 GByte di memoria. Nel nostro studio precedente del comportamento del tessuto suino, abbiamo stabilito un modello costitutivo materiale anisotropo per le regioni discendente e spirale del colon suino basato su prove simultanee di inflazione ed estensione15. Questo modello ha dimostrato in diversi studi precedenti di catturare adeguatamente il comportamento anisotropo del colon19,34,35,36. I parametri del materiale stabiliti nel nostro studio precedente sono stati mantenuti in quanto, a nostra conoscenza, è l’unico studio che modella un colon di grandi animali sotto inflazione ed estensione, i due principali modi di deformazione del colon. Inoltre, è stato stabilito che il colon suino è rilevante per lo studio della diverticolosi33,37. Il modello di materiale selezionato e i parametri sono quindi scelte adeguate per questo studio. Poiché l’anatomia del colon a spirale non è rilevante per gli esseri umani, i risultati dai campioni del colon discendente sono stati utilizzati qui per caratterizzare il materiale del tessuto del colon. Abbiamo riportato nel nostro lavoro precedente che c’è una variazione trascurabile nel diametro localmente lungo la lunghezza di un segmento del colon discendente suino, quindi il colon normale è stato modellato come un segmento cilindrico (Fig. 1a). Il raggio interno ed esterno sono stati impostati a \({R}_{i}=11,5\) mm e \({R}_{o}=12,7\) mm, rispettivamente, in base ai valori medi misurati nel nostro lavoro precedente. La lunghezza assiale è stata impostata a \(L=10\) cm in modo che sia diverse volte maggiore del raggio. In particolare, il tessuto del colon è stato assunto come un materiale iperelastico omogeneo, anisotropo e incomprimibile il cui comportamento meccanico è descritto dalla seguente funzione di energia di deformazione anisotropa \(\bar{W}:\)15,35

$$\bar{W}={C}_{10}({I}_{1}-3)+\frac{k}_{1}^{l}}{k}_{2}^{l}}+\frac{k}_{1}^{s}{k}_{2}^{s}}$$
(1)

Il primo termine rappresenta una risposta Neo-Hookean che caratterizza il comportamento dei costituenti non collageni del tessuto del colon. La quantità \({I}_{1}={{lambda}^{2}+{lambda}^{2}+frac{1}{{lambda}^{2}{lambda}^{2}{lambda}^{2}}) rappresenta la prima invariante del tensore di deformazione con \({\lambda}_{z}) e \({\lambda}_{\theta}) che si riferiscono agli stiramenti assiali e circonferenziali, rispettivamente. I termini successivi caratterizzano il contributo delle fibre di collagene, con \({k}_{1}} e \({k}_{2}}) che sono una misura della rigidità delle fibre. In particolare, il secondo termine (con apice \(l\)) caratterizza il comportamento delle fibre di collagene allineate lungo la direzione longitudinale (assiale), mentre il terzo termine (con apice \(s\)) caratterizza la risposta delle fibre di collagene disperse in due direzioni simmetriche preferenziali rispetto alla direzione circonferenziale. La quantità \({I}_{4}\ è indicata come pseudo-invariante e caratterizza la risposta meccanica delle fibre lungo le direzioni preferenziali. Per il secondo termine, \({I}_{4}^{l}={\lambda }_{z}^{2}}}) poiché tiene conto delle fibre lungo la direzione longitudinale. Per il terzo termine, \({I}_{4}^{s}}) è espresso come:

$${I}_{4}^{s}={lambda }_{theta }^{2}{cos }^{2}{gamma }^{s}+{lambda }_{z}^{2}{sin }^{2}{gamma }^{s}$$
(2)

qui, \(\pm {\gamma }^{s}) indica l’angolo dei due orientamenti simmetrici delle fibre rispetto alla direzione circonferenziale. Cinque set di parametri del materiale, precedentemente derivati per cinque campioni di colon discendente suino, sono stati usati in questo lavoro per i parametri del materiale \({C}_{10},\,{k}_{1}^{l},\,{k}_{2}^{l},\,{k}_{1}^{s},\,{k}_{2}^{s},\) e \({\gamma }^{s}}). Questa funzione di energia di deformazione non è direttamente disponibile in Abaqus. È stata quindi integrata per la prima volta implementando una subroutine utente basata su Fortran chiamata UANISOHYPER_INV, che permette di definire il comportamento anisotropo iperelastico del materiale usando la formulazione invariante. Questa subroutine richiede la definizione delle derivate prima e seconda della funzione deformazione-energia (Eq. 1) rispetto agli invarianti scalari \({I}_{1}}) e \({I}_{4}}). Per verificare che la subroutine è stata implementata correttamente e convalidare il modello computazionale del tessuto del colon, abbiamo confrontato i valori di sollecitazione stimati dal modello della sollecitazione circonferenziale isocorica \({bar{sigma }}_{{{m}}}) e della sollecitazione assiale isocorica \({bar{sigma}_{z}^{m}) sulla superficie esterna del colon contro i valori ottenuti analiticamente sulla base della funzione energia di deformazione e calcolati come segue:

$$$${bar{sigma }{{theta }^{a}={lambda }_{theta}=frac{parziale \bar{W}}{parziale {lambda (3)

(3)
$$$${bar{sigma }_{z}^{a}={lambda }_{z}={frac{parziale
(4)

Il coefficiente di determinazione \(R^{2}) è stato calcolato come segue per quantificare la deviazione tra i valori computazionali e analitici:

$${R}^{2}(A)=1-\frac{{\sum }_{q=1}^{n}{({A}_{q}^{m}-{A}_{q}^{a})}^{2}}{{\sum }_{q=1}^{n}{({A}_{avg}^{a}-{A}_{q}^{a})}^{2}$$
(5)

dove \(A={bar{sigma}_{z},\e il pedice “avg” indica la media dei valori analitici su tutti i punti dati \(n). Un valore di \R}^{2}} vicino a 1 indica che una buona correlazione è globalmente ottenuta.

Modello del diverticolo

Per i modelli del colon malato, una struttura simile a un diverticolo è stata modellata come l’intersezione di una sfera di diametro \D\ con la forma cilindrica del colon normale a un’altezza \H\ sopra la superficie. Solo una sacca è stata considerata in questo studio per concentrarsi sull’effetto di una singola sacca senza variabili aggiuntive che sarebbero incorporate da sacche multiple (numero di sacche, posizioni lungo il colon, posizione relativa rispetto all’altra, variabilità nelle dimensioni individuali, ecc.) A nostra conoscenza, non esiste uno studio che riporti sistematicamente la geometria dei diverticoli. Viene comunemente menzionato che essi sono tipicamente tra i 5 e i 10 mm di diametro senza un’adeguata indicazione su dove viene misurato il diametro38,39. Per coprire una varietà di dimensioni della tasca, nelle simulazioni sono stati considerati tre valori per l’altezza iniziale della tasca (H) e il diametro iniziale della tasca (D): \H=2 mm, 4 mm e 6 mm, e D = 8 mm, 10 mm e 12 mm. Le combinazioni di questi valori di D e H hanno portato a 9 modelli malati con diverse dimensioni della tasca. Un raccordo rotondo di raggio \({r}_{f}=2\) mm è stato incluso all’intersezione tra il cilindro che rappresenta il colon normale e la sfera che rappresenta la sacca per rimuovere i bordi taglienti che potrebbero non essere incontrati nel tessuto reale e potrebbero alterare i valori di stress. Lo spessore del tessuto del diverticolo è stato assunto uguale a quello del colon (cioè, \({R}_o}-{R}_i}\)). Ci sono stati tentativi di caratterizzare la meccanica del colon nei diverticoli23,30,40,41. Nessuna proprietà del materiale, tuttavia, sono stati stabiliti specificamente per il tessuto diverticolo (cioè, la caratterizzazione del sacchetto isolato diverticolo). Proprietà simili al tessuto normale sono stati quindi considerati per la sezione pouch.

Condizioni limite

Le condizioni limite sono state imposte sulla base del nostro studio precedente di colon suino15. Abbiamo osservato che la pressione luminale (pressione sulla parete interna del tessuto con pressione esterna uguale a zero) valori superiori a 1,5 kPa hanno portato a danni permanenti al tessuto in condizioni passive. Nelle simulazioni di calcolo, questo valore è stato spinto ulteriormente per visualizzare la tendenza dell’evoluzione dei valori di stress e la pressione luminale fino a 5 kPa è stata imposta sulla superficie interna del colon, compresa la tasca. La pressione della parete esterna è stata mantenuta uguale a 0. Poiché stiamo lavorando con un tessuto a parete sottile (rapporto diametro interno/spessore vicino a 20), la pressione luminale potrebbe quindi essere assimilata come il delta tra la pressione interna (ad esempio, la pressione imposta all’interno del condotto del tessuto dalla materia fecale o dal gas intestinale) e la pressione esterna (ad esempio, la pressione addominale). Tipicamente, sarebbe necessaria una piccola pressione per aprire il tessuto del colon che altrimenti collasserebbe. Tale pressione (~0,2 kPa) è, tuttavia, molto piccola rispetto alla gamma considerata qui, quindi viene ignorata. La pressione è stata applicata con incrementi quasi statici. Un allungamento assiale di circa il 10% è stato tipicamente visto in vivo nel colon durante il nostro studio precedente. Così, un tratto assiale del 10% è stato imposto lungo la direzione \(Z\) prima della pressurizzazione. Abbiamo trovato valori di angolo di apertura molto piccoli nel tessuto del colon suino, e le sollecitazioni circonferenziali residue sono state quindi trascurate nel modello. Per entrambi i modelli normali e malati, abbiamo usato la simmetria rispetto al piano \(YZ\) e al piano \(XY\) e abbiamo risolto solo un quarto del modello per ottenere una maggiore efficienza, prima di ricostruire l’intero modello a scopo di visualizzazione.

Meshing

Si è fatto in ogni caso un adeguato sezionamento della geometria per assicurare che il meshing sia possibile usando elementi esaedrici lineari. La formulazione degli elementi ibridi è stata selezionata per imporre l’incomprimibilità del materiale (cioè l’elemento C3D8H dalla libreria Abaqus). Uno studio di convergenza della mesh è stato condotto per il modello normale riducendo la dimensione della mesh di riferimento a partire dalla dimensione predefinita suggerita automaticamente da Abaqus (1,3 mm) fino a quando la variazione del valore massimo delle sollecitazioni principali massime tra due mesh successive era inferiore al 10%. La dimensione convergente della mesh è stata poi utilizzata anche nei modelli malati, con alcuni ulteriori affinamenti intorno alla tasca. Un totale di 69608 elementi e 88350 nodi sono stati utilizzati per il modello normale alla dimensione della mesh più fine mantenuta per l’analisi (un quarto della geometria). Tra 69764 elementi e 80656 elementi (88605 e 102440 nodi rispettivamente) sono stati utilizzati per i nove modelli malati (quarto della geometria). Un orientamento locale del materiale è stato specificato utilizzando l’opzione di formulazione discreta disponibile in Abaqus per definire correttamente l’orientamento delle fibre in ogni elemento come richiesto nella funzione di energia di deformazione.

Analisi dei risultati

Sono stati considerati un totale di 10 modelli (uno normale e nove malati con diverse dimensioni della tasca). Ogni modello è stato simulato con 5 diversi set di parametri del materiale per un totale di 50 casi simulati. Una misura del livello di stress è necessaria per confrontare tra i modelli. Poiché non stiamo stabilendo un criterio di fallimento, ma semplicemente guardando le variazioni in condizioni diverse, qualsiasi misura come lo stress principale massimo o lo stress di von Mises sarebbe adeguato. I risultati sono mostrati qui solo in termini di stress principale massimo per evitare informazioni ridondanti, poiché le osservazioni complessive con lo stress di von Mises sono state trovate simili. L’analisi e i risultati con la tensione di von Mises sono inclusi nei dati associati a questo studio se il lettore desidera consultarli. Per ogni caso simulato, la sollecitazione principale massima (al centroide di un elemento della mesh), indicata come \({\sigma }_{MPS}\), è stata osservata in un intervallo di pressione da 0 a 5 kPa con un incremento di 0,25 kPa. Uno script Python è stato sviluppato per cercare automaticamente il valore massimo di \({sigma }_{MPS}}, designato come \({sigma }_{MPS,max}) alle varie pressioni per ogni caso simulato. Per ogni modello, il valore medio ottenuto per ogni caso simulato (cioè, per tutti e cinque i set di parametri del materiale), è stato calcolato e conservato per la successiva analisi. Questi valori sono indicati come \({sigma }{MPS,\,max}^{avg}\) e \({sigma }_{MPS,\,max}^{avg,\,n}\ per i modelli malati e normali, rispettivamente.

I valori \({sigma }_{MPS,\ max}^{avg}\), sono stati correlati, per pressioni tra 1 e 5 kPa (con incrementi di 0,25 kPa), a diversi parametri della geometria della tasca (Fig. 1c) per identificare la relazione tra la geometria della tasca e lo stress elevato: la larghezza del collo della tasca lungo l’asse x, la larghezza del collo della tasca lungo l’asse z, l’area alla base del collo della tasca e la superficie della tasca sul lato del lume. Il coefficiente di correlazione di Pearson \({r}^{2}} è stato calcolato per ogni parametro. Il p-value a due code è stato calcolato anche per valutare la significatività della correlazione.

Inoltre, abbiamo calcolato la pressione \({P}_{d}\ nei modelli malati in modo tale che per una data pressione \({P}_{n}\ nel colon normale abbiamo:

$$max ({sigma }_{MPS,max}({P}_{d}))\le \,\max ({sigma }_{MPS,\,max}^{n}({P}_{n}))$$
(6)

Il seguente fattore di pressione fu successivamente calcolato come:

$$f==sinistra(\frac{P}_{n}-{P}_{d}}{P}_{n}}} destra)\volte 100$$
(7)

Questo fattore indica la % di cui la pressione \({P}_{n}})deve essere ridotta nel colon malato in modo che \({\sigma }_{MPS,max) sia simile a \({sigma }{MPS,\,max}^{n}), cioè il valore massimo della sollecitazione principale massima.e il valore massimo della sollecitazione principale massima nel modello malato non supera il valore corrispondente nel modello normale. La pressione ridotta è quindi designata come \({P}_{d}\). In particolare, per ogni caso simulato, \({P}_{n}}) è stato ridotto di 0,25 kPa fino a quando la relazione dell’Eq. 6 è stata soddisfatta, fornendo il valore corrispondente di \({P}_{d}\ e il fattore di pressione \(f\) è stato calcolato. Questo è stato fatto per \({P}_{n}\) nell’intervallo da 2 a 5 kPa (con incrementi di 1 kPa). Il valore medio ottenuto per i modelli su tutti i casi simulati è stato calcolato e conservato per l’analisi successiva.

Infine, abbiamo stabilito la zona di influenza di una sacca, cioè la distanza intorno a una sacca dove lo stress è elevato. Per ottenere questo, abbiamo stabilito per ogni caso malato una mappa dello stress principale massimo sulla superficie laterale del lume del colon (poiché lo stress è più alto su questa superficie) lungo i nodi di magliatura situati sui percorsi longitudinali e circonferenziali che attraversano il centro della sacca. L’evoluzione dello stress principale massimo lungo questi percorsi è stata osservata in funzione delle distanze \({d}_{l}\ e \({d}_c}\ misurate rispetto al centro della tasca. Poi, abbiamo calcolato le distanze più lontane \({d}_{l}^{inf}) lungo il percorso longitudinale e \({d}_c}^{inf}) lungo il percorso circonferenziale dopo le quali lo stress principale massimo scende entro il 10% del valore corrispondente nel modello normale lungo gli stessi percorsi. Queste distanze costituiscono la nostra misura della zona di influenza. Osservando che raggiungono un valore di plateau dopo una certa pressione, i loro valori medi su tutti i casi per ogni modello a 5 kPa sono stati correlati a diversi parametri della geometria della sacca, in modo simile a \({\sigma }_{MPS,\,max}^{avg}\).

Tutti i calcoli sono stati fatti in linguaggio di programmazione Python sviluppando i notebook Jupyter42. I dati sono stati memorizzati e manipolati con la libreria Pandas di Python43. Le operazioni matematiche di base sono state fatte con la libreria NumPy44. Il coefficiente di correlazione \({r}^{2}\) e il p-value a due code sono stati calcolati con la funzione “stats” della libreria Scipy45. Tutti i grafici sono stati realizzati con la libreria Matplotlib46.

Lascia un commento

Il tuo indirizzo email non sarà pubblicato.