Bab III Model Difusi Oksigen di Jaringan dengan Laju Konsumsi Konstan
Pada bab ini, akan dibahas penyebaran oksigen di pembuluh kapiler dan jaringan, dimana sel-sel di jaringan diasumsikan mengkonsumsi oksigen secara konstan. Berdasarkan fungsi Michaelis-Menten, konsumsi oksigen konstan ini terjadi ketika nilai konsentrasi di jaringan cukup besar, atau sel-sel di jaringan telah mencapai keadaan jenuh (saturated ).
III.1
Solusi Keadaan Tunak
Darah yang kaya akan oksigen mengalir dari jantung melalui arteri kemudian masuk ke kapiler. Karena telah diasumsikan bahwa darah merupakan cairan yang homogen, maka konsentrasi oksigen di kapiler dalam arah radial adalah konstan. Sehingga untuk daerah kapiler, perbedaan konsentrasi hanya terjadi dalam arah aksial. Perbedaan konsentrasi ini disebabkan oleh adanya aliran darah dari ujung awal kapiler. Dimisalkan konsentrasi oksigen di kapiler pada posisi z˜ adalah c˜k (˜ z ), dan konsentrasi di ujung awal kapiler (˜ z = 0) adalah cin .
Jaringan per satuan volume mengkonsumsi oksigen dengan laju konstan, misalkan sebesar g0 . Perhatikan Gambar 3.1. Kita tinjau posisi z˜ = z ∗ dan segmen jaringan dengan panjang ∆˜ z . Dalam keadaan tunak, konsumsi oksigen di segmen tersebut adalah sebesar π(b2 − a2 )∆˜ z g0 . Besarnya konsumsi oksigen di jaringan diimbangi oleh penyediaan oksigen dari kapiler sebanyak huiπa2 [˜ ck (z ∗ ) − c˜k (z ∗ + ∆˜ z )], dimana hui adalah rata-rata kecepatan darah. Oleh karena itu, di kapiler berlaku hubungan: π(b2 − a2 )∆˜ z g0 = huiπa2 [˜ ck (z ∗ ) − c˜k (z ∗ + ∆˜ z )] . π(b2 − a2 )g0 = −huiπa2
[˜ ck (z ∗ + △˜ z ) − c˜k (z ∗ )] . △˜ z
(3.1) (3.2)
12 Δz
z*
z=0
z=l
Gambar 3.1: Silinder Kapiler-Jaringan. Untuk △˜ z → 0, persamaan konsentrasi untuk daerah kapiler adalah: π(b2 − a2 )g0 = −huiπa2
∂˜ ck . ∂ z˜
(3.3)
Untuk daerah jaringan, karena b ≪ l, kita dapat mengabaikan difusi dalam
arah aksial. Dalam kasus konsumsi oksigen konstan dan t˜ → ∞, persamaan difusi di jaringan (2.6) adalah: 0 = Dj
c ∂ 2 c˜ 1 ∂˜ − g0 . + ∂˜ r 2 r˜ ∂˜ r
(3.4)
Dinding kapiler, yang menjadi batas antara kapiler dan jaringan diasumsikan tidak memiliki ketahanan perpindahan massa. Akibatnya, nilai konsentrasi oksigen pada r˜ = a sama dengan nilai konsentrasi oksigen di kapiler c˜k . Oleh karena itu, syarat batas untuk persamaan (3.3) dan (3.4) adalah: c˜(a, z˜) = c˜k (˜ z ), ∂˜ c (b, z˜) = 0, ∂˜ r c˜k (0) = cin .
Masalah di atas dapat dituliskan dalam bentuk persamaan tidak berdimensi, dengan penskalaan: c˜ = cin c,
r˜ = ar,
z˜ = az.
Sehingga diperoleh: ∂ 2 c 1 ∂c + = M, ∂r 2 r ∂r
(3.5)
dan π(b2 − a2 )g0 = −huiπacin
dck , dz
(3.6)
13 dimana M =
g0 a2 , c a Dj
dengan batas:
c(1, z) = ck (z), ∂c b , z = 0, ∂r a ck (0) = 1.
Solusi masalah syarat batas tersebut, untuk daerah kapiler adalah:
ck (z) = 1 − N dimana N =
g0 a . huicin
b2 − 1 z, a2
(3.7)
Untuk daerah jaringan diperoleh solusi (lihat Lampiran
A): c(r, z) = 1 − M
r2 − 1 b2 ln r − 2a2 4
−N
b2 − 1 z. a2
(3.8)
Gambar 3.2 menunjukkan penyebaran konsentrasi dalam arah radial, sedangkan Gambar 3.3 dalam arah aksial. Pada kedua gambar tersebut, digunakan nilai M = 2.39 × 10−4 , N = 3.19 × 10−4 , dan
b a
= 11.
z=0 z=10 z=20
1.2
1
c
0.8
0.6
0.4
0.2
0
0
2
4
6
8
10
r
Gambar 3.2: Grafik Konsentrasi Oksigen dalam Arah Radial untuk Laju Konsumsi Konstan.
14
1 kapiler r=4 r=11
0.9 0.8 0.7
c
0.6 0.5 0.4 0.3 0.2 0.1 0
0
5
10 z
15
20
Gambar 3.3: Grafik Konsentrasi Oksigen dalam Arah Aksial untuk Laju Konsumsi Konstan. Berdasarkan Gambar 3.3, dapat dinyatakan bahwa nilai konsentrasi oksigen di kapiler turun sebanding dengan jarak terhadap inlet. Karena ketersediaan oksigen semakin menurun, maka untuk daerah jaringan pun semakin jauh dari inlet, nilai konsentrasi oksigen semakin kecil.
Dalam arah radial, semakin jauh jarak dari dinding kapiler, konsentrasi oksigen di jaringan bernilai semakin kecil. Karena setiap daerah jaringan mengkonsumsi oksigen dengan laju konstan, maka selisih antara konsentrasi di dinding kapiler dan di dinding luar jaringan menunjukkan nilai yang sama untuk semua z.
III.2
Solusi Keadaan Tidak Tunak
Untuk keadaan tidak tunak ini, hanya akan ditinjau untuk daerah jaringan. Dengan mengabaikan arah aksial, persamaan difusi untuk daerah jaringan adalah: ∂˜ c = Dj ∂ t˜
c ∂ 2 c˜ 1 ∂˜ − g0 , + ∂˜ r 2 r˜ ∂˜ r
a ≤ r˜ ≤ b.
(3.9)
Ketika t˜ = 0, konsentrasi oksigen di dinding kapiler adalah ca , dan jaringan dalam keadaan setimbang. Sehingga untuk kondisi awal di daerah jaringan
15 adalah: g0 c˜(˜ r , 0) = ca − Dj
b2 r˜ r˜2 − a2 ln − 2 a 4
.
(3.10)
Kemudian di dalam kapiler mengalir darah yang kaya akan oksigen. Misalkan konsentrasi oksigen di dinding kapiler adalah ci , dimana ci > ca . Oleh karena itu, dalam selang waktu tertentu (0 ≤ t˜ ≤ ε) konsentrasi oksigen di dinding kapiler mengalami kenaikan dari ca menuju ci , diasumsikan kenaikannya linier. Selanjutnya (t˜ > ε) konsentrasi oksigen di dinding kapiler adalah konstan, yaitu sebesar ci . Dalam model matematika, keadaaan di dinding kapiler dapat dituliskan: c(a, t˜) =
ci −ca ˜ t + ca , ε
c, i
jika 0 ≤ t˜ < ε; jika t˜ ≥ ε.
Pada r˜ = b, tidak terdapat aliran konsentrasi yang menembus dinding jaringan, sehingga: ∂˜ c (b) = 0. ∂˜ r
(3.11)
Dengan penskalaan: 2
c˜ = ca c,
r˜ = ar,
a t˜ = t , Dj
(3.12)
masalah nilai awal dan syarat batas (3.9)-(3.11) menjadi tidak berdimensi, yaitu: −
∂c ∂ 2 c 1 ∂c + + = M. ∂t ∂r 2 r ∂r
c(1, t) = f (t) =
b c1 −1 t δ
b c1 ,
∂c b ( , t) = 0. ∂r a
g0 a2 , c a Dj
b c1 =
ci , ca
dan δ =
+ 1, jika 0 ≤ t < δ; jika t ≥ δ. (3.14)
c(r, 0) = 1 − M dimana M =
(3.13)
b2 ln r r 2 − 1 , − 2a2 4
(3.15)
εDj . a2
Grafik dari c(1, t), untuk parameter δ = 1.6, b c1 = 1.25 ditunjukkan oleh Gambar 3.4.
16
1.35
1.3
1.25
f(t)
1.2
1.15
1.1
1.05
1
0
2
4
6
8
10
t
Gambar 3.4: Syarat Batas pada r=1.
III.2.1
Solusi Analitik
Untuk menyelesaikan masalah (3.13)-(3.15), kita tulis c(r, t) dalam bentuk: c(r, t) = u(r) + w(r, t).
(3.16)
Persamaan (3.16) memenuhi: ∂u ∂w ∂c = + , ∂r ∂r ∂r ∂c ∂w = , ∂t ∂t sehingga persamaan (3.13) menjadi ∂ 2 u 1 ∂u ∂w ∂ 2 w 1 ∂w + − + 2 + = M. ∂r 2 r ∂r ∂t ∂r r ∂r
(3.17)
Oleh karena itu untuk menyelesaikan masalah nilai awal dan syarat batas (3.13)-(3.15), adalah dengan mencari solusi u(r) dan w(r, t) dari masalah: 1. ∂ 2 u 1 ∂u + = M, ∂r 2 r ∂r
(3.18)
dengan batas:
∂u ∂r
u(1) = 1, b , t = 0. a
(3.19) (3.20)
17 2. −
∂w ∂ 2 w 1 ∂w + 2 + = 0, ∂t ∂r r ∂r
(3.21)
dengan syarat batas dan nilai awal: w(1, t) = c(1, t) − u(1) = f (t) − 1, ∂w b ,t = 0, ∂r a w(r, 0) = c(r, 0) − u(r).
Berdasarkan solusi keadaan tunak, solusi untuk u(r) adalah: 2 r2 − 1 b . ln r − u(r) = 1 − M 2a2 4
(3.22) (3.23) (3.24)
(3.25)
Karena u(r) = c(r, 0) maka persamaan (3.24) menjadi w(r, 0) = 0.
(3.26)
Untuk mencari solusi w(r, t), interval t kita bagi menjadi dua bagian, yaitu 0 ≤ t ≤ δ dan t > δ. Untuk 0 ≤ t ≤ δ, persamaan (3.21) beserta nilai awal dan syarat batasnya, akan diselesaikan dengan teorema Duhamel. Untuk itu, perlu dicari terlebih dahulu s˜(r, t), dimana s˜(r, t) adalah solusi dari persamaan (3.21)(3.24) dengan syarat batas s˜(1, t) = 1. Menggunakan metode transformasi Laplace, yang secara rinci terdapat pada lampiran E, diperoleh: s˜(r, t) = 1 − π
∞ X n=1
exp
−α2n t
2 b G(r, αn ) αn J1 (αn ) , a F (αn )
dimana J1 merupakan fungsi Bessel orde 1, αn merupakan akar dari persamaan: b b −Y1 (αn )J0 (αn ) + J1 (αn )Y0 (αn ) = 0, a a
(3.27)
18 dan G(r, αn ) = Y0 (rαn )J0 (αn ) − J0 (rαn )Y0 (αn ), b F (αn ) = (αn J0 (αn ))2 − (αn J1 (αn ))2 . a Grafik untuk s˜(r, 0) adalah: sHr, 0L 1.5
1.0
0.5
0.0
r 2
3
4
5
6
-0.5
Gambar 3.5: Grafik s˜(r, 0) untuk 40 suku pertama. Maka berdasarkan teorema Duhamel, untuk 0 ≤ t ≤ δ, berlaku: w(r, t) =
Z
τ =t
τ =0
s˜(r, t − τ )
b c1 − 1 δ
dτ.
" # 2 ∞ X b G(r, α ) b c1 − 1 2 n t−π (1 − exp−αn t ) J1 (αn ) . w(r, t) = δ a F (α ) n n=1 (3.28) Untuk t > δ, kita punya masalah nilai awal dan syarat batas: −
∂w ∂ 2 w 1 ∂w + 2 + = 0. ∂t ∂r r ∂r
w(1, t) = c(1, t) − u(1) = b c1 . ∂w b ,t = 0. ∂r a Z τ =δ b c1 − 1 dτ. w(r, 0) = s(r, δ − τ ) δ τ =0
(3.29)
(3.30) (3.31) (3.32)
Menggunakan metode pemisahan variabel, diperoleh: w(r, t) =
∞ X n=1
2
An exp−αn t G(r, αn ) + b c1 ,
(3.33)
19 dimana An =
=
R
b a
1
R
(w(r, 0) − b c1 )rG(r, αn )dr Rb a rG(r, αn )dr 1 b a
1
(w(r, 0) − b c1 )rG(r, αn )dr (παn J1 (αn ab ))2 . 2 J0 (αn )2 − J1 (αn ab )2
Sehingga c(r, t) = 1 − M
b2 ln r r 2 − 1 − 2a2 4
+ w(r, t),
(3.34)
dengan w(r, t) merupakan persamaan (3.28) dan (3.33).
Untuk memverifikasi solusi analitik yang telah diperoleh, masalah (3.13)-(3.15) akan diselesaikan dengan metode numerik.
III.2.2
Metode Numerik
Untuk menyelesaikan masalah (3.13)-(3.15) digunakan metode beda maju. Domain R = {(r, t) : 1 ≤ r ≤ ab , 0 ≤ t ≤ te } dipartisi menjadi n × m buah persegi panjang, dengan sisi ∆r = h dan ∆t = k. Nilai c(r, t) diaproksimasi pada titik-titik ri = 1 + ih dan tj = jk, dimana i = 0, 1, ..., n dan j = 0, 1, ..., m. Untuk penyederhanaan c(ri , tj ) dinotasikan dengan ci,j .
Aproksimasi untuk cr (r, t), crr (r, t) dan ct (r, t) adalah:
1 ∂c = [c(r + h, t) − c(r − h, t)] + O(h2 ). ∂r 2h 1 ∂2c = 2 [c(r + h, t) − 2c(r, t) + c(r − h, t)] + O(h2 ). 2 ∂r h ∂c 1 = [c(r, t + k) − c(r, t)] + O(k). ∂t k
(3.35) (3.36) (3.37)
Dengan mengabaikan O(k) dan O(h2 ), persamaan (3.35)-(3.37) disubstitusi ke persamaan (3.13), diperoleh: −
ci,j+1 − ci,j ci+1,j − 2ci,j + ci,j+1 ci+1,,j − 2ci,j + ci,j+1 = M. + + k h2 2hri,j
(3.38)
20 Berdasarkan persamaan (3.38), nilai ci,j+1 dapat ditentukan dari nilai-nilai ci,j , ci−1,j dan ci+1,j yang telah diketahui. Sehingga untuk i = 1, ..., n − 1 dan j = 0, ..., m: ci,j+1 = ci,j + ˚ A [ci+1,j − 2ci,j + ci−1,j ] + dimana ˚ A=
k ,ß h2
ß [ci+1,j − ci−1,j ] − kM. ri,j
(3.39)
k . 2h
=
Nilai ci,0 , dengan i = 0, 1, ..., n, dapat diperoleh dari nilai awal yaitu: ci,0 = c(ri , 0).
(3.40)
Sedangkan c0,j , dengan j = 0, 1, ..., m diperoleh dari nilai batas pada r = 1, yaitu: c0,j = f (tj ).
(3.41)
Karena c′ ( ab ) = 0, maka untuk i = n, dan j = 0, ..., m berlaku persamaan: ci,j+1 = ci,j + 2˚ A [ci−1,j − ci,j ] − kM.
(3.42)
Persamaan (3.39) dan (3.42) akan stabil jika dan hanya jika 0 ≤
≤ 21 .
k h2
Gambar 3.6 menunjukkan perbandingan solusi analitik (garis penuh) dan numerik (garis ∗) untuk parameter
b a
t=0
t=10
numerik analitik
numerik analitik
1.3
1.3
1.2
1.2
1.1
1.1
c
c
= 11, M = 1.2 × 10−3, δ = 1.6, b c1 = 1.25.
1
1
0.9
0.9
0.8
1
2
3
4
5
6 r
7
8
9
10
11
0.8
1
2
3
4
5
6 r
7
8
9
10
11
Gambar 3.6: Perbandingan Solusi Analitik dan Numerik untuk t = 0 dan t = 10. Selanjutnya, Gambar 3.7 menunjukkan proses penyebaran konsentrasi oksigen di jaringan.
21
1.3 t=0 t=0.5 t=1.6 t=20 t=200
1.25 1.2 1.15
c
1.1 1.05 1 0.95 0.9 0.85 0.8
1
2
3
4
5
6 r
7
8
9
10
11
Gambar 3.7: Penyebaran Konsentrasi Oksigen di Jaringan untuk Laju Konsumsi Konstan. Pada selang waktu dimana konsentrasi oksigen di dinding kapiler meningkat, daerah yang dekat dengan dinding kapiler mengikuti peningkatan tersebut. Selanjutnya untuk waktu yang cukup lama, sebaran konsentrasi oksigen mendekati keadaan tunak.
