Model de colon
Programul de elemente finite Abaqus (versiunea 6.13, SIMULIA, Providence, RI) a fost utilizat pentru a dezvolta modelele și a efectua simulările asociate. Toate simulările au fost efectuate pe o stație de lucru Dell T7810 cu două procesoare Intel Xeon și 32 GBytes de memorie. În studiul nostru anterior privind comportamentul țesutului porcin, am stabilit un model constitutiv al materialelor anizotrope pentru regiunile descendentă și spirală ale colonului porcin pe baza unor teste simultane de umflare și extensie15. Acest model s-a dovedit în mai multe studii anterioare că surprinde în mod adecvat comportamentul anizotropic al colonului19,34,34,35,36. Parametrii materialelor stabiliți în studiul nostru anterior au fost reținuți deoarece, după cunoștințele noastre, acesta este singurul studiu care modelează un colon de animal mare în condiții de umflare și extensie, cele două moduri majore de deformare a colonului. În plus, s-a stabilit că colonul porcin este relevant pentru studiul diverticulozei33,37. Modelul de material și parametrii selectați sunt, prin urmare, alegeri adecvate pentru acest studiu. Având în vedere că anatomia colonului spiralat nu este relevantă pentru om, rezultatele obținute din probele de colon descendent au fost utilizate aici pentru caracterizarea materialului tisular al colonului. Am raportat în lucrarea noastră anterioară că există o variație neglijabilă a diametrului la nivel local pe lungimea unui segment de colon descendent de porc, astfel încât colonul normal a fost modelat ca un segment cilindric (Fig. 1a). Raza interioară și cea exterioară au fost stabilite la \({R}_{i}=11,5\) mm și, respectiv, \({R}_{o}=12,7\) mm, pe baza valorilor medii măsurate în lucrările noastre anterioare. Lungimea axială a fost stabilită la \(L=10\) cm, astfel încât să fie de câteva ori mai mare decât raza. În mod specific, s-a presupus că țesutul colonului este un material hiperelastic omogen, anizotrop și incompresibil, al cărui comportament mecanic este descris de următoarea funcție anizotropă de energie de deformare \(\bar{W}:\)15,35
Primul termen reprezintă un Neo-Hookean care caracterizează comportamentul constituenților necolagenari ai țesutului colonului. Cantitatea \({I}_{1}={\lambda }_{z}^{2}+{\lambda }_{\theta }^{2}+\frac{1}{{{\lambda }_{z}^{2}{\lambda }_{\theta }^{2}}}\) reprezintă prima invariant al tensorului de deformare cu \({\lambda }_{z}\) și \({\lambda }_{\theta }\) referindu-se la întinderile axiale și circumferențiale, respectiv. Termenii următori caracterizează contribuția fibrelor de colagen, \({k}_{1}\) și \({k}_{2}\) fiind o măsură a rigidității fibrelor. Mai exact, al doilea termen (cu superscriptul \(l\)) caracterizează comportamentul fibrelor de colagen aliniate de-a lungul direcției longitudinale (axiale), în timp ce al treilea termen (cu superscriptul \(s\)) caracterizează răspunsul fibrelor de colagen dispersate în două direcții preferențiale simetrice în raport cu direcția circumferențială. Cantitatea \({I}_{4}\) este denumită pseudo-invariantă și caracterizează răspunsul mecanic al fibrelor de-a lungul direcțiilor preferențiale. Pentru cel de-al doilea termen, \({I}_{4}^{l}={\lambda }_{z}^{2}\), deoarece reprezintă fibrele de-a lungul direcției longitudinale. Pentru cel de-al treilea termen, \({I}_{4}^{s}\) se exprimă astfel:
Aici, \(\pm {\gamma }^{s}\) indică unghiul dintre cele două unghiuri de orientare simetrică a fibrelor în raport cu direcția circumferențială. Cinci seturi de parametri de material, derivate anterior pentru cinci probe de colon descendent de porc, au fost utilizate în această lucrare pentru parametrii de material \({C}_{10},\,{k}_{1}^{l},\,{k}_{2}^{l},\,{k}_{1}^{s},\,{k}_{2}^{s},\) și \({\gamma }^{s}\). Această funcție de energie de deformare nu este disponibilă direct în Abaqus. Prin urmare, a fost integrată pentru prima dată prin implementarea unei subrutine de utilizator bazată pe Fortran, numită UANISOHYPER_INV, care permite definirea comportamentului material hiperelastic anizotropic hiperelastic folosind formularea invariantă. Această subrutină a necesitat definirea derivatelor prima și a doua ale funcției energie-deformație (Ecuația 1) în raport cu invarianții scalari \({I}_{1}\) și \({I}_{4}\). Pentru a verifica dacă subrutina a fost implementată corect și pentru a valida modelul de calcul al țesutului colonului, am comparat valorile estimate de model ale tensiunii izocorice circumferențiale \({\bar{\sigma }}}_{\theta }^{m}\) și ale tensiunii izocorice axiale \({\bar{\sigma }}_{z}^{m}\) pe suprafața exterioară a colonului cu valorile obținute analitic pe baza funcției de energie a deformației și calculate după cum urmează:
Coeficientul de determinare \({R}^{2}\) a fost calculat după cum urmează pentru a cuantifica abaterea dintre valorile de calcul și cele analitice:
unde \(A={\bar{\sigma }}_{z},\,{\bar{\sigma }}_{\theta }\), iar indicele „\(avg\)” indică media valorilor analitice pe toate punctele de date \(n\). O valoare a lui \({R}^{2}\) apropiată de 1 indică faptul că se obține o bună corelație la nivel global.
Modelul diverticulului
Pentru modelele de colon bolnav, o structură asemănătoare unui diverticul a fost modelată ca fiind intersecția unei sfere cu diametrul \(D\) cu forma cilindrică a colonului normal la o înălțime \(H\) deasupra suprafeței. În acest studiu a fost luată în considerare o singură pungă pentru a se concentra asupra efectului unei singure pungi fără variabilitățile suplimentare care ar fi încorporate de mai multe pungi (numărul de pungi, pozițiile de-a lungul colonului, poziția relativă una față de cealaltă, variabilitatea dimensiunii individuale etc.). După cunoștințele noastre, nu există niciun studiu care să raporteze sistematic geometria diverticulilor. Se menționează în mod obișnuit că acestea au de obicei un diametru cuprins între 5 și 10 mm, fără a se indica în mod corespunzător unde este măsurat diametrul38,39. Pentru a acoperi o varietate de dimensiuni ale pungilor, în cadrul simulărilor au fost luate în considerare trei valori pentru înălțimea inițială a pungii \(H\)și diametrul inițial al pungii \(D\): \(H=2\\) mm, 4 mm și 6 mm, iar \(D\) = 8 mm, 10 mm și 12 mm. Combinațiile acestor valori ale lui D și H au condus la 9 modele bolnave cu diferite dimensiuni ale pungii. La intersecția dintre cilindrul care reprezintă colonul normal și sfera care reprezintă punga a fost inclus un cordon rotund de rază \({r}_{f}=2\) mm pentru a elimina marginile ascuțite care ar putea să nu fie întâlnite în țesutul real și care ar putea modifica valorile tensiunii. S-a presupus că grosimea țesutului diverticulului este aceeași cu cea a colonului (adică \({R}_{o}-{R}_{i}\)). Au existat încercări de caracterizare a mecanicii colonului în diverticuloză23,30,40,41. Cu toate acestea, nu au fost stabilite proprietăți materiale în mod specific pentru țesutul diverticulului (adică caracterizarea pungii izolate a diverticulului). Proprietăți similare cu cele ale țesutului normal au fost astfel luate în considerare pentru secțiunea pungii.
Condiții limită
Condițiile limită au fost impuse pe baza studiului nostru anterior asupra colonului porcin15. Am observat că valorile presiunii luminale (presiunea pe peretele interior al țesutului cu presiunea exterioară egală cu zero) de peste 1,5 kPa au dus la deteriorarea permanentă a țesutului în condiții pasive. În cadrul simulărilor computaționale, această valoare a fost împinsă și mai mult pentru a vizualiza tendința de evoluție a valorilor de stres și a fost impusă o presiune luminală de până la 5 kPa pe suprafața internă a colonului, inclusiv pe pungă. Presiunea din peretele exterior a fost menținută egală cu 0. Deoarece lucrăm cu un țesut cu pereți subțiri (raportul dintre diametrul interior și grosime apropiat de 20), presiunea luminală ar putea fi astfel asimilată ca fiind delta dintre presiunea interioară (de exemplu, presiunea impusă în interiorul canalului țesutului de materiile fecale sau gazele intestinale) și presiunea exterioară (de exemplu, presiunea abdominală). În mod obișnuit, ar fi necesară o presiune mică pentru a deschide țesutul colonului care altfel s-ar prăbuși. Această presiune (~0,2 kPa) este, totuși, foarte mică în comparație cu intervalul luat în considerare aici, astfel încât este ignorată. Presiunea a fost aplicată în creșteri cvasi-statice. O întindere axială de aproximativ 10% a fost observată de obicei in vivo în colon în timpul studiului nostru anterior. Astfel, a fost impusă o întindere axială de 10% de-a lungul direcției \(Z\) înainte de presurizare. Am găsit valori foarte mici ale unghiului de deschidere în țesutul colonului de porc, iar tensiunile circumferențiale reziduale au fost astfel neglijate în model. Atât pentru modelele normale, cât și pentru cele bolnave, am folosit simetria în raport cu planul \(YZ\) și planul \(XY\) și am rezolvat doar un sfert din model pentru a obține o eficiență mai mare, înainte de a reconstrui întregul model în scopul vizualizării.
Mesurare
S-a făcut o secționare corespunzătoare a geometriei în fiecare caz pentru a ne asigura că plasarea este posibilă folosind elemente hexaedrice liniare. A fost selectată o formulare de element hibrid pentru a impune incompresibilitatea materialului (de exemplu, elementul C3D8H din biblioteca Abaqus). S-a efectuat un studiu de convergență a ochiurilor de plasă pentru modelul normal prin reducerea dimensiunii ochiurilor de plasă de referință pornind de la dimensiunea implicită sugerată automat de Abaqus (1,3 mm) până când variația valorii maxime a tensiunilor principale maxime între două ochiuri de plasă succesive a fost mai mică de 10 %. Dimensiunea convergentă a ochiurilor de plasă a fost apoi utilizată și în modelele bolnave, cu unele rafinări suplimentare în jurul pungii. Un total de 69608 elemente și 88350 noduri au fost utilizate pentru modelul normal la cea mai fină dimensiune a ochiurilor de plasă reținută pentru analiză (un sfert din geometrie). Pentru cele nouă modele bolnave au fost utilizate 69764 elemente și 80656 elemente (88605 și, respectiv, 102440 noduri) (un sfert din geometrie). O orientare locală a materialului a fost specificată folosind opțiunea de formulare discretă disponibilă în Abaqus pentru a defini în mod corespunzător orientarea fibrelor în fiecare element, așa cum este necesar în funcția de energie de deformare.
Analiza rezultatelor
Au fost luate în considerare un total de 10 modele (unul normal și nouă bolnave cu diferite dimensiuni ale pungilor). Fiecare model a fost simulat cu 5 seturi diferite de parametri de material, rezultând un total de 50 de cazuri simulate. Este necesară o măsură a nivelului de stres pentru a compara între modele. Având în vedere că nu stabilim un criteriu de eșec, ci pur și simplu analizăm variațiile în condiții diferite, orice măsură, cum ar fi tensiunea principală maximă sau tensiunea von Mises, ar fi adecvată. Rezultatele sunt prezentate aici doar în ceea ce privește solicitarea principală maximă pentru a evita informațiile redundante, deoarece observațiile generale cu solicitarea von Mises s-au dovedit a fi similare. Analiza și rezultatele cu tensiunea von Mises sunt incluse în datele asociate cu acest studiu, dacă cititorul dorește să le consulte. Pentru fiecare caz simulat, tensiunea principală maximă (la centroidul unui element al ochiului de plasă), denumită \({\sigma }_{MPS}\), a fost observată într-un interval de presiune de la 0 la 5 kPa, cu un increment de 0,25 kPa. A fost dezvoltat un script Python pentru a căuta în mod automat valoarea maximă a \({\sigma }_{MPS}\), desemnată ca \({\sigma }_{MPS,max}\) la diferite presiuni pentru fiecare caz simulat. Pentru fiecare model, valoarea medie obținută pentru fiecare caz simulat (adică pentru toate cele cinci seturi de parametri de material) a fost calculată și reținută pentru analiza ulterioară. Aceste valori sunt denumite \({\sigma }_{MPS,\,max}^{avg}\) și \({\sigma }_{MPS,\,max}^{avg,\,n}\) pentru modelele bolnave și, respectiv, normale.
Valorile \({\sigma }_{MPS,\,max}^{avg}\), au fost corelate, pentru presiuni între 1 și 5 kPa (în trepte de 0,25 kPa), cu diferiți parametri ai geometriei pungii (Fig. 1c) pentru a identifica relația dintre geometria pungii și tensiunea ridicată: lățimea gâtului pungii de-a lungul axei \(x\) \({D}_{x}\), lățimea gâtului pungii de-a lungul axei \(z\) \({D}_{z}\), suprafața de la baza gâtului pungii \({A}_{b}\) și suprafața pungii pe partea lumenului \({S}_{p}\). Coeficientul de corelație Pearson \({r}^{2}\) a fost calculat pentru fiecare parametru. S-a calculat, de asemenea, valoarea p cu două cozi pentru a evalua semnificația corelației.
Mai mult, am calculat presiunea \({P}_{d}\) în modelele bolnave, astfel încât pentru o presiune dată \({P}_{n}\) în colonul normal avem:
Următorul factor de presiune a fost ulterior calculat astfel:
Acest factor indică procentul cu care trebuie redusă presiunea \({P}_{n}\) în colonul bolnav astfel încât \({\sigma }_{MPS,max}\) să fie similară cu \({\sigma }_{MPS,\,max}^{n}\), adică.adică valoarea maximă a tensiunii principale maxime în modelul bolnav nu depășește valoarea corespunzătoare în modelul normal. Presiunea redusă este apoi desemnată ca \({P}_{d}\). Mai exact, pentru fiecare caz simulat, \({P}_{n}\) a fost redusă cu 0,25 kPa până când a fost satisfăcută relația din Ecuația 6, furnizând valoarea corespunzătoare a \({P}_{d}\) și a fost calculat factorul de presiune \(f\). Acest lucru s-a făcut pentru \({P}_{n}\) în intervalul 2 – 5 kPa (în trepte de 1 kPa). Valoarea medie obținută pentru modelele din toate cazurile simulate a fost calculată și reținută pentru analiza ulterioară.
În cele din urmă, am stabilit zona de influență a unei pungi, adică distanța din jurul unei pungi unde tensiunea este ridicată. Pentru a realiza acest lucru, am stabilit pentru fiecare caz bolnav o hartă a stresului principal maxim pe suprafața laterală a lumenului colonului (deoarece stresul este cel mai ridicat pe această suprafață) de-a lungul nodurilor de angrenare situate pe traseele longitudinale și circumferențiale care trec prin centrul pungii. Evoluția tensiunii principale maxime de-a lungul acestor trasee a fost observată în funcție de distanțele \({d}_{l}\) și \({d}_{c}\) măsurate în raport cu centrul sacului. Apoi, am calculat cele mai îndepărtate distanțe \({d}_{l}^{inf}\) de-a lungul traiectoriei longitudinale și \({d}_{c}^{inf}\) de-a lungul traiectoriei circumferențiale, după care tensiunea principală maximă scade în limita a 10% din valoarea corespunzătoare în modelul normal de-a lungul aceleiași traiectorii. Aceste distanțe au constituit măsura noastră a zonei de influență. Observând că acestea ating o valoare de platou după o anumită presiune, valorile lor medii pe toate cazurile pentru fiecare model la 5 kPa au fost corelate cu diferiți parametri ai geometriei pungii, în mod similar cu \({\sigma }_{MPS,\,max}^{avg}\).
Toate calculele au fost efectuate în limbajul de programare Python prin dezvoltarea de caiete Jupyter42. Datele au fost stocate și manipulate cu ajutorul bibliotecii Pandas din Python43. Operațiile matematice de bază au fost efectuate cu biblioteca NumPy44. Coeficientul de corelație \({r}^{2}\) și valoarea p cu două cozi au fost calculate cu ajutorul funcției „stats” a bibliotecii Scipy45. Toate graficele au fost realizate cu biblioteca Matplotlib46.
.