Analizirane su tri moguće metode za određivanje parametara Voigtova modela arterijske stijenke na temelju mjerenja vremenske promjene dvije veličine u poprečnom presjeku arterije. Prva metoda zahtijeva mjerenje tlaka i površine poprečnog presjeka, druga tlaka i prostorne derivacije protoka, a treća protoka i prostorne derivacije površine poprečnog presjeka. Metode su uspoređene na način da su se umjesto mjerenja koristili rezultati numeričke simulacije u koje je uvedena slučajna pogreška s ciljem oponašanja realnih mjerenja. Zaključeno je da je prva metoda najmanje osjetljiva na pogrešku mjerenja, a najlošija je treća

Korade, Ivan

Krizmanić, Severino

Virag, Zdravko

FAMENA Repository

Šesti susret Hrvatskoga društva za mehaniku Rijeka, 29-30. svibnja 2014.     PROCJENA MOGUĆNOSTI MJERENJA LOKALNIH PARAMETARA VISKOELASTIČNE ARTERIJSKE STIJENKE  Korade, I., Virag, Z. & Krizmanić, S.  Sažetak: Analizirane su tri moguće metode za određivanje parametara Voigtova modela arterijske stijenke na temelju mjerenja vremenske promjene dvije veličine u poprečnom presjeku arterije. Prva metoda zahtijeva mjerenje tlaka i površine poprečnog presjeka, druga tlaka i prostorne derivacije protoka, a treća protoka i prostorne derivacije površine poprečnog presjeka. Metode su uspoređene na način da su se umjesto mjerenja koristili rezultati numeričke simulacije u koje je uvedena slučajna pogreška s ciljem oponašanja realnih mjerenja. Zaključeno je da je prva metoda najmanje osjetljiva na pogrešku mjerenja, a najlošija je treća.  Ključne riječi: Hemodinamika, Arterijsko stablo, Metoda mjerenja parametara modela stijenke  1 UVOD  Jedna od prepreka široj primjeni jednodimenzijskih modela strujanja krvi u arterijskom stablu je nepoznavanje lokalnih svojstava stijenke arterija. Poznato je da velike arterije pokazuju svojstvo viskoelastičnosti. U slučaju Voigtova modela arterijska stijenka je karakterizirana s popustljivošću (elastičnošću) arterijske stijenke, i viskoznošću (otporom) stijenke. U nedostatku mjernih metoda za neinvazivno mjerenja viskoznosti stijenke najčešće se pretpostavlja čisto elastična priroda stijenke. Ovdje će se analizirati mogućnost mjerenja parametara viskoelastične stijenke. U radu će se koristiti rezultati simulacije strujanja krvi u arterijskom stablu dobiveni rješavanjem jednodimenzijskog matematičkog modela temeljenog na Voigtovoj konstitutivnoj relaciji. Izračunate vrijednosti protoka, tlaka i površine poprečnog presjeka cijevi, te vrijednosti njihovih prostornih derivacija, će se uzimati za mjerene vrijednosti. Analizirat će se iz koje kombinacije mjerenih vrijednosti će se moći najtočnije rekonstruirati traženi parametri Voigtova modela. Uvođenjem slučajne pogreške u „mjerene“ rezultate, ocijenit će se osjetljivost metoda za određivanje parametara arterijske stijenke.  2 SIMULACIJA STRUJANJA KRVI U VELIKIM ARTERIJAMA   2.1 Matematički model jednodimenzijskog strujanja krvi Problem jednodimenzijskog strujanja krvi u velikim arterijama je opisan sljedećim sustavom jednadžbi [1]: 0A Qt x∂ ∂+ =∂ ∂,      (1) ( )QvQ A pfQt x xρ∂∂ ∂+ + = −∂ ∂ ∂,     (2) 0 01( )Ap p A AC tη∂− = − +∂,     (3) gdje su: x i t prostorna i vremenska koordinata, p je tlak, Q je protok, A je površina poprečnog presjeka arterije (A=D2π/4), v je prosječna brzina po presjeku (v=Q/A), C je površinska popustljivost arterijske stijenke (C=dA/dp), η je otpor stijenke, ρ je gustoća fluida, f je faktor trenja, ovdje definiran Hagen-Pouisseuilevim zakonom f=32µ/(ρD2), a µ je viskoznost fluida. 0A  označuje konstantnu vrijednost površine A pri tlaku 0p p= .  2.2 Metoda karakteristika za rješavanje jednodimenzijskog modela  Matematički model, opisan sustavom jednadžbi od (1) do (3), efikasno se rješava numeričkom metodom karakteristika. Za potrebe te metode arterijsko stablo se podijeli na određeni broj elemenata duljine ∆x, na čijim se rubovima nalaze čvorovi. Slika 1. prikazuje dva susjedna elementa (označena sa j i k) koja su omeđena čvorovima  I, J i K. Vrijednosti tlaka su definirane u čvorovima, površina A poprečnog presjeka elementa se definira u sredini elementa (i smatra se konstantnom duž elementa duljine ∆x), a protoci Q se definiraju na oba ruba svakog elementa. Na taj način se uz svaki element pamte četiri varijable. Na primjer, nepoznanice koje se odnose na element j su: tlak Jp , protok na lijevom rubu jLQ ,  i desnom jDQ  te površina jA , kao što je prikazano na slici 1. Detaljan opis metode je dan u [2].   Sl.1: Prikaz dijela diskretiziranog arterijskog stabla s rasporedom varijabli.  2.3 Simulacija mjerenih rezultata Rezultate simulacije dobivene opisanom metodom će se koristiti umjesto mjerenih podataka, temeljem kojih će se određivati parametri C i η arterijske stijenke. S obzirom da u stvarnim mjerenjima uvijek postoji pogreška u dobivene rezultate simulacije ćemo uvesti slučajnu pogrešku (primjenom generatora slučajnih brojeva) u obliku: ( )max 2 1Mp p p r= + ∆ − ,      (4) gdje je Mp simulirani tlak s unesenom pogreškom, p je tlak dobiven iz simulacije, maxp∆  maksimalna pogreška tlaka, a r slučajni broj u intervalu nula do jedan. Analogno se uvode pogreške i u rezultate za Q i A. S obzirom da su u numeričkoj metodi p i Q definirani na rubu elementa, a površina na sredini elementa, p i Q u sredini elementa će se definirati kao aritmetička sredina vrijednosti s ruba, a parcijalne derivacije na sredini elementa (npr. elemenata j, prema slici 1), definiraju se izrazima:  j J I K k jA  kA  Jp  Ip KpjDQjLQ kLQ kDQ  TJQ  jD jLJ I     i     M MM MM M Q Qp pp Qx x x x−−∂ ∂= =∂ ∆ ∂ ∆,    (5)  3 METODE MJERENJA SVOJSTAVA ARTERIJSKE STIJENKE  3.1 Metode određivanja svojstava arterijske stijenke Cilj je definirati postupak određivanja popustljivosti C i otpora η arterijske stijenke, koji se direktno pojavljuju u jednadžbi (3). Jasno je da se iz simultanog mjerenja vremenske promjene tlaka i vremenske promjene površine poprečnog presjeka arterije mogu odrediti ova dva parametra. Iz jednadžbe (3) je jasno da za slučaj η =0, tlak i površina moraju biti u fazi, te bi se C mogao odrediti iz mjerenih signala tlaka i promjene površine kao nagib pravca u p-A ravnini (npr. primjenom metode najmanjih kvadrata). U slučaju kada postoji otpor stijenke, signali tlaka i vremenske promjene površine su povezani linearnom diferencijalnom jednadžbom, te je problem pogodno promatrati u frekvencijskoj domeni. U tom pristupu se vremenski signal tlaka i površine razviju u Fourierov red: ( ) ( ) ( )0 00sin cosNp pn nnp t S n t C n tω ω= = + ∑ ,   (6) gdje koeficijenti Fourierova reda čine fazor tlaka u kompleksnom prostoru ˆ p pn n np S iC= +  na kružnoj frekvenciji 0nω ω= , gdje n označuje redni broj harmonika, a 0 2 /Tω π=  osnovnu kružnu frekvenciju definiranu periodom T srčanog takta. Analogno se može definirati i fazor ˆ A An n nA S iC= +  iz mjerenog signala vremenske promjene površine A. Temeljem jednadžbe (3) se može definirati 1. metoda za određivanje C i η, na osnovu mjerenja p i A (ovdje ne ulazimo u raspravu kako bi se to mjerilo) u jednom presjeku, koja glasi: 0ˆ1,    gdje je ,  a 1,  2, 3, ....ˆpi n nC Aωη ω ω+ = = = ,   (7) Teorijski gledano vrijednosti C i η bi trebale biti jedne te iste na svim frekvencijama, što naravno u mjerenjima neće biti slučaj. Može se ili definirati neka srednja vrijednost iz vrijednosti dobivenih na različitim frekvencijama, ili vrijednost dobivenu na osnovnoj frekvenciji proglasiti za najtočniju vrijednost, budući da se na višim frekvencijama može očekivati sve veći utjecaj slučajne pogreške. 2. metoda se može izvesti iz vremenski derivirane jednadžbe (3) u koju se uvrsti jednadžba (1), te se dobije: 21p Q Qt C x t xη∂ ∂ ∂= − −∂ ∂ ∂ ∂,     (8) ili u frekvencijskoj domeni 0ˆ1   gdje je ,  a 1,  2, 3, ....ˆi pi n nC Qxωωη ω ω+ = − = =∂∂,  (9) Jasno je da u ovoj metodi treba mjeriti tlak i prostornu derivaciju protoka (protok u dvije točke ili razliku protoka u dvije bliske točke).  3. metoda se dobije iz jednadžbe (3), derivirane po prostornoj koordinati, nakon što se u nju uvrsti pojednostavljena jednadžba (2) (sa zanemarenim trenjem i konvekcijskim članom) koja glasi: 0 0AQ pt xρ∂ ∂+ =∂ ∂,      (10) Gledano u frekvencijskoj domeni, dobije se: 00ˆ1   gdje je ,  a 1,  2, 3, ....ˆi Qi n nC A Axρ ωωη ω ω+ = − = =∂∂,  (11) Jasno je da su se pri izvodu prve dvije metode koristile originalne jednadžbe jednodimenzijskog modela, dok se kod ove metode koristila pojednostavljena jednadžba (2).  4 USPOREDBA METODA  Primjenom opisane metode simulirano je strujanje fluida u račvastoj cijevi, koja je shematski prikazana na slici 2.   Sl. 2: Shematski prikaz elemenata račvaste cijevi  Cijev je podijeljena na 38 elemenata duljine ∆x=6 cm. Prvih deset elemenata je jednolike površine A0=5 cm2 (na pretlaku p0=0), a poprečni presjek granajućih cijevi je A1=3 cm2. Na ulazu u cijev (čvor 1), propisan je protok, oblika polusinusnog vala trajanja T/3 (T=0,3 s) amplitude 50 ml/s. Na izlazima iz cijevi (u čvorovima 25 i 39), zadan je otpor R25=1 mmHg·s/ml i R39 =3 mmHg·s/ml. U proračunu su parametri C i η stijenke cijevi zadani indirektno konstantnom brzinom širenja vala cw=10 m/s, pri čemu je C=A/(ρ·cw2) i vremenskom konstantom τ=C·η =0.05 s. Treba primijetiti da je uz konstantni cw parametar C vremenski promjenjiv jer se mijenja površina A, a onda je uz konstantni τ promjenjiv i parametar η. Međutim, s obzirom da su vremenske promjene površine A male, promjene parametara će se zanemariti. Period integracije ovog problema je 10T, a za simulaciju mjerenja se uzimaju rezultati na petom elementu i to iz zadnjeg perioda integracije, kada se postiglo ustaljeno oscilirajuće strujanje u cijevi. 27 ... 38 39 1 2 3 ... 10 11 12 13 ... 24 25 26 27 ... 38 39 2 3 ... 10 11 12 13 ... 24 25 26 nC(ml/mmHg)pogreskaC(%)1 2 3 4 5 6-0.2-0.15-0.1-0.0500.050.10.15-1000100200300400500CeCerCMETODA 2nC(ml/mmHg)pogreskaC(%)1 2 3 4 5 600.020.040.060.080.10.12050100150200CeCerCMETODA 1nη(ml/mmHg)pogreskaη(%)1 2 3 4 5 600.10.20.30.40.50.6-20020406080100ηeηerηs dodanom pogreskom 5%nη(ml/mmHg)pogreskaη(%)1 2 3 4 5 600.10.20.30.40.50.6-20020406080100ηeηerηs dodanom pogreskom 5%nη(ml/mmHg)pogreskaη(%)1 2 3 4 5 600.10.20.30.40.50.6-20020406080100ηeηerηs dodanom pogreskom 5%nC(ml/mmHg)pogreskaC(%)1 2 3 4 5 6-0.2-0.100.1-1000100200300400500CeCerCMETODA 3 Sl. 3. Vrijednosti parametara C i η u funkciji rednog broja harmonika, rekonstruiranih s tri metode, pri čemu je u simuliranom signalu uvedena pogreška 5 %.  Tablica 1. prikazuje vrijednosti zadanih (točnih) parametara C i η izračunatih pri srednjoj površini poprečnog presjeka A0. U istoj tablici su dane rekonstruirane vrijednosti tih parametara primjenom triju metoda, te njihovo odstupanje od točnih vrijednosti i to za dva slučaja: bez uvođenja slučajne pogreške u rezultat simulacije i s pogreškom od 5 %. Slika 3. prikazuje rezultate za slučaj uvedene pogreške od 5 %. Jasno je da u slučaju da u signalu nema pogreške sve tri metode vrlo dobro rekonstruiraju vrijednosti parametara C i η (vrijednosti relativne pogreške za η su za prva tri harmonika ispod 1 %, jedino treća metoda pokazuje malo veće odstupanje parametra C na višim harmonicima). Kad se u „mjereni signal“ uvede pogreška od 5%, očito je da sve metode na višim harmonicima rekonstuiraju parametre stijenke sa znatnijim odstupanjem. Pa ipak prve dvije metode daju dobre rezultate na prvom harmoniku.    Tablica 1. Vrijednosti rekonstruiranih parametara stijenke, temeljem simulacije i pogreška u odnosu na zadane (točne) vrijednosti pri A0, te odstupanja u odnosu na točne vrijednosti, erC=100(Ce- C)/Ce, er η =100(ηe- C)/ ηe  n Točno Bez uvođenja pogreške Uvedena pogreška od 5 % 100 Ce ηe 100 C erC η erη 100 C erC η erη Metoda 1 1 9.91 0.504 9.95 -0.366 0.502 0.559 10.03 -1.15 0.501 0.626 2 9.93 -0.203 0.504 0.081 10.35 -4.38 0.500 0.779 3 9.90 0.160 0.508 0.687 8.88 10.5 0.514 -1.87 4 9.78 1.325 0.515 -2.00 14.1 -42.5 0.613 -21.5 5 10.1 -1.81 0.517 -2.52 822 -8189 0.478 5.29 6 10.1 -1.84 0.529 -4.78 7.09 28.5 0.422 16.4 Metoda 2 1 9.91 0.504 9.93 -0.199 0.503 0.383 9.97 -0.544 0.550 -9.10 2 9.98 -0.693 0.503 0.360 8.74 11.8 0.518 -2.61 3 10.1 -1.45 0.503 0.323 4.06 59.1 0.420 16.7 4 10.2 -2.65 0.503 0.234 9.82 0.95 0.213 57.9 5 10.2 -3.18 0.502 0.470 -7.08 171 -0.180 135.7 6 10.4 -4.98 0.503 0.385 -5.06 151 0.133 73.6 Metoda 3 1 9.91 0.504 9.84 0.709 0.503 0.204 9.28 6.35 0.457 9.49 2 9.68 2.32 0.504 0.140 6.18 37.7 0.538 -6.75 3 9.35 5.67 0.504 0.071 24.5 -146 0.363 28.1 4 8.43 15.0 0.505 -0.176 3.25 67.1 -0.067 113.2 5 10.4 -4.89 0.500 0.880 -7.29 173 -0.010 119.8 6 9.04 8.82 0.501 0.598 -6.13 162 0.037 92.7  5 ZAKLJUČAK  U radu su uspoređene tri metode za određivanje parametara stijenke u Voigtovom modelu. Za mjerenja su uzeti rezultati numeričke simulacije u koje su unesene slučajne pogreške. Metoda koja se temelji na simultanom mjerenju tlaka i površine poprečnog presjeka se pokazala najmanje osjetljivom na uvedenu pogrešku, s tim da se za parametre uzimaju vrijednosti iz prvog harmonika.  Literatura [1] Kitawaki, T., “Numerical Simulation Model with Viscoelasticity of Arterial Wall”,  Chapter 9, http://dx.doi.org/10.5772/49976, (2012). [2] Korade, I., Virag, Z., Šavar, M., “Numerical Simulation Of One-Dimensional Flow In Elastic And Viscoelastic Branching Tube” prihvaćen rad za 11th World Congress on Computational Mechanics (WCCM XI), Barcelona, Spain, 20.-25. July 2014.  Autori: Ivan Korade, Zdravko Virag, Severino Krizmanić. Sveučilište u Zagrebu Fakultet strojarstva i brodogradnje, Zavod za energetska postrojenja, energetiku i ekologiju, Katedra za mehaniku fluida, Ivana Lučića 5, 10000 Zagreb, tel. 01/6168-222, e-mail: ivan.korade@fsb.hr, zdravko.virag@fsb.hr, severino.krizmanic@fsb.hr. 

Procjena mogućnosti mjerenja lokalnih parametara viskoelastične arterijske stijenke

Abstract

Similar works

Full text

Available Versions

FAMENA Repository