逆運算法在生物熱傳邊界估算之研究
蕭美枝
國立高雄應用科技大學 模具工程系
E-mail:[email protected]
摘 要
本文以循序運算法(sequential approach method)發展出解生物熱傳邊界直接解與逆向解,利用 Pennes
equation 描述生物組織熱傳導的模型和已知的生物組織熱物理參數,求解組織內部溫度場分佈數值解,再以溫
度分布數值解為已知,應用逆運算法由若干個量測點反求解生物組織的邊界溫度輸入或熱傳輸入;未來將可
應用於生物醫學上,提供非侵入式的熱診斷、熱治療、低溫手術與熱參數估算臨床相關資料。
關鍵辭:循序運算法、逆運算法、生物熱傳、Pennes equation
1. 前 言
在處理工程問題時,當邊界狀態無法直接量測時,以內部或外部幾點溫度反求邊界狀態者稱之為逆向熱
傳問題,藉由逆向熱傳分析可求解許多內部熱源、表面熱流、表面溫度,過去幾年已廣泛應用於許多的設計、
製造[1-11],現今應用生物醫學的生物熱傳有熱診斷[12]、熱治療[13-15]、低溫手術[16-17]與熱參數估算[18-21]
等,在生物醫學領域研究數量比起上述工程領域的研究相對的少,在表皮量取幾點溫度為已知,精確估算內
部組織溫度分佈是屬於直接數值解[22-29 ];若以估算之直接數值解為已知,反求生物組織的邊界輸入溫度或
熱傳的逆向數值解研究則更少[30-32],本文以 Pennes equation 為基礎描述熱傳導模式,利用循序運算法精確
的估算出生物組織內部溫度分佈直接解,再以估算之直接解為已知,可精確的求解未知邊界溫度或熱流的逆
向解,將來提供生物醫學在非侵入式的熱診斷、熱治療、低溫手術與熱參數估算臨床相關資料。
2. 數學模型
考慮一 1 維熱參數為常數的組織模型如圖 1 所示,為簡化模型假設組織為均質(homogeneous)、等方
(isotropic)健康的肌肉組織,在組織方向離表皮深度 L 無組織生成熱與空間熱,探討表皮(skin surface)處
施加溫度或熱流,在體內(body core)的溫度分布,解直接解與逆向解。
以 Pennes equation[33] 描述 1 維生物組織的熱傳導模式如下示所示:
L
Skin suface
Body core
x
1
圖
r m a b 2 2
Q
Q
)
T
T
(
c
x
T
k
t
T
c
+
−
+
+
∂
∂
=
∂
∂
ω
ρ
ρ
(1)
其中 ρ 、c
、 分別表示密度、定容與熱傳導係數,
kω
b表示血流灌注率,
是組織生成熱,
是空間
熱,
是動脈溫度假設 37oC。
m Q Qr a T邊界條件(B.C.)如下:
(1) 溫度輸入
1 T T =,
x
=
0
(2)
0 x T = ∂ ∂,
x
=
L
(3)
(2) 震盪 cos 型熱流輸入
) t ( f x T = ∂ ∂,
x=0(4)
0 x T = ∂ ∂,
x
=
L
(5)
因假設無組織生成熱與空間熱,所以
Qm與
Qr為零。令
θ
=T−Ta,式(1)整理如下式所示:
θ
ρ
ω
θ
α
θ
k
c
x
t
b 2 2−
∂
∂
=
∂
∂
(6)
其中
c kρ
α
=邊界條件式(2)~(5)分別整理如下式所示:
(1) 溫度輸入
1θ
θ
=,
x
=
0
(7)
0 x = ∂ ∂θ
,
x
=
L
(8)
(2) 震盪 cos 型熱流輸入
) ( f xθ
θ
= ∂ ∂,
x
=
0
(9)
0 x = ∂ ∂θ
,
x
=
L
(10)
2.1 數學方法與運算
2.1.1 直接解
以有限差分法將式(6)離散化,如下示所示:
k i b 2 k 1 i k i k 1 i k i 1 k i k c x 2 tθ
ρ
ω
Δ
θ
θ
θ
α
Δ
θ
θ
− + − = − − + −(11)
其中下標 為空間格點位置,上標
i
k
為時間格點位置,重新排列後得到向量表示式,如下:
} S { ] k [ } ]{ B [ ] k [ } {θ
k = −1θ
k−1 + −1(12)
其中
[K]、
[B]為溫度分佈向量之係數矩陣,
{S}為熱源向量之係數矩陣。
以循序運算法推導在時間格點
m
的溫度方程式:
}
S
]{
D
[
}
]{
C
[
}
{
θ
m=
θ
m−1+
(13)
其中
[C]=[k]−1[B],
[D]= [k]−1求解在組織內部的溫度分佈直接解。
2.1.2 逆向解
利用式(13)推導得
}
S
]{
D
][
C
[
}
S
]{
D
[
]
C
[
...
}
S
]{
D
[
]
C
[
}
S
]{
D
[
]
C
[
}
{
]
C
[
}
{
θ
m+r−2=
r−1θ
m−1+
r−1 m+
r−2 m+1+
+
2 m+r−1+
m+r(14)
} S ]{ D ][ C [ } S ]{ D [ ] C [ ... } S ]{ D [ ] C [ } S ]{ D [ ] C [ } { ] C [ } {θm+r−1 = r θm−1 + r m + r−1 m+1 + + 2 m+r−1 + m+r(15)
}
{
θ
m、
{
θ
m+1}
…
{
θ
m+r−1}
分別表示在時間
m、
…
的溫度分佈,得矩陣表示型式如下:
t
m 1 t +t
m+r−1{
m 1 r 1 r 2 1 r m 2 r m 1 m m 2 1 r r 2 r 1 r 2 2 1 r m 2 r m 1 m m]
C
[
]
C
[
]
C
[
]
C
[
}
S
{
}
S
{
}
S
{
}
S
{
]
D
[
]
C
[
]
C
[
]
C
[
]
C
[
0
]
C
[
]
C
[
]
C
[
]
C
[
]
C
[
]
C
[
0
0
]
C
[
}
{
}
{
}
{
}
{
− − − + − + + − − − − + − + +⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
+
⎪
⎪
⎪
⎪
⎭
⎪
⎪
⎪
⎪
⎬
⎫
⎪
⎪
⎪
⎪
⎩
⎪
⎪
⎪
⎪
⎨
⎧
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
=
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
θ
θ
θ
θ
θ
M
M
M
M
L
L
M
M
M
O
O
M
M
L
L
M
M
}
(16)
取出第 個時間格點,分量形式表示如下:
i
} S ]{ D ][ u [ } ]{ C ][ u [ } { i 1 m i m i =θ
− +θ
(17)
} S ]{ D ][ C ][ u [ } S ]{ D [ ] C ][ u [ ... } S ]{ D [ ] C ][ u [ } S ]{ D [ ] C ][ u [ } { ] C ][ u [ } {θm+r−1 = i r θm−1 + i r m + i r−1 m+1 + + i 2 m+r−1 + i m+r(18)
整理成下式:
m m , m i 0 , m i m i=
a
+
a
S
θ
(19)
m 1 m , 1 m i m m , 1 m i 0 , 1 m i 1 m ia
a
S
a
S
+ + + + +=
+
+
θ
(20)
...
...
S
a
S
a
a
im r 1,0 im r 1,m m im r 1,m 1 m 1 1 r m i=
+
+
++
+
+ − + − + − + − +θ
(21)
其中
aim,m=
=…=
=
(22)
1 m , 1 m i a + + aim+r−1,m+r−1 ei0 m , 1 m ia
+=
m 2,m 1=…=
=
(23)
i a + + aim+r−1,m+r−2e
1i m , 1 r m i a + −=
aim+r,m+1=…=
aim+r−1,m=
eir−1(24)
0 , m i a=
[ui][C]{θ
m−1}(25)
0 , 1 r m i a + −=
[ui][C]r{θ
m−1}(26)
為了提高逆運算過程的穩定度,Beck 等提出未來時間的論述(future time)[34],將未知參數在暫時區間內
假設為常數,因此得下式:
m
S
=
Sm+1=…=
Sm+r−1(27)
m 0 i 0 , m i m m , m i 0 , m i m i =a +a S =a +e S
θ
(28)
m 0 i 1 i 0 , 1 m i m 1 m , 1 m i m , 1 m i 0 , 1 m i 1 m i =a +(a +a )S =a +(e +e )S + + + + + +θ
(29)
m 0 i 2 r i 1 r i 0 , 1 r m i m 1 r m , 1 r m i 1 m , 1 r m i m , 1 r m i 0 , 1 r m i 1 r m i =a +(a +a +...a )S =a +(e +e +...+e )S − − − + − + − + + − + − + − + − +θ
(30)
定義
θ
im+k =aim+k,0 +EikSm =aim+k,0 +[Eik ]{Sm}(31)
其中
=
∑
,
= k 0 l l i k ie
E
k
=
0
,
1
,
2
,...
r
−
1
(32)
式(28)~(30)以向量矩陣式表示成下式:
S
Φ
ϑ =
(33)
⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ − − − − = − + = + − + = + + = 0 , 1 r m i 1 r m i 0 , 2 r m i 2 r m i 0 , 1 m i 1 m i 0 , m i m i a a a aθ
θ
θ
θ
ϑ
M M,
,
⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = − − 1 r i 2 r i 1 i 0 i E E E EM
M
Φ
S={Sm }(34)
逆向分析時僅需從溫度ϑ 中量取少數的分量來確定未知邊界,建立一個以
ϑ
meas、
Φ 的方程式,如下:
S
null null meas⎥
⎦
⎤
⎢
⎣
⎡
=
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
Φ
Φ
ϑ
ϑ
再藉由線性最小平方誤差法(linear least-squares error method),反求估算出
向量中的熱源強度,方法
如下:
S
meas T 1 T ) ( Sˆ =Φ
Φ
−Φ
ϑ
(35)
其中
T 1為轉換矩陣(reverse matrix)
) (Φ
Φ
−2.2 量測點的平均誤差百分比
量測點的平均誤差百分比定義如下:
)
100
%
Sˆ
Sˆ
Sˆ
(
N
1
N 1 i 2 exact meas exact×
∑
−
=
=η
,其中
N是步階數。
3. 結果與討論
本文以循序運算法分別求解 1 維生物組織未知邊界溫度或熱流輸入,其中組織參數參考[35-36]所示,
,
,
,
3 m / kg 1000=
ρ
c=
4000OCJ/kg k=
0.5OCW /mω
b=
0.0005l/s,組織深度
L=
3cm,模擬時間 6000
秒,每 20 秒為一時間步階。範例 1、2 分別為模擬未知邊界為溫度輸入與震盪 cos 型熱流輸入情況,求出直
接解後並反求出逆向解。
其中逆向解以直接解加入隨機亂數(
λ )成為量測溫度,如下式所示:
Ymeas=
Yexact+
λσ
,
λ 範圍如下:
576 . 2 576 . 2<
<
−
λ
,可靠度達 99%。
範例 1:邊界為溫度輸入
在表皮
x=
0處,輸入溫度
o,組織內部
1 =25θ
x=
3cm處,為絕熱邊界
0 x=
∂
∂θ
位於組織深度
L=
3cm處
取 1 點量測,模擬求得直接解,該處溫度分佈圖,如圖 2 所示。
圖 2
L=
3cm處,時間 6000 秒溫度分佈
若以所得直接解模擬反求邊界溫度之逆向解,在未來時間(
r=1)且無誤差(
σ
=
0)的理想狀態下,得估算
邊界溫度結果與實際邊界溫度輸入相符,如圖 3 所示。當未來時間
r= 並考慮量測誤差
8σ
=
0.1%,其平均誤
差百分比為 3.84790859671%與實際溫度輸入值非常接近,如圖 4 所示。
60.0000000000000 61.0000000000000 62.0000000000000 63.0000000000000 64.0000000000000 1 5 9 13 17 21 25 29 measure temperature exact temperature
圖 3
r=1、
σ
=0,1 點量測估算邊界溫度與實際邊界溫度
0.0000000000000 20.0000000000000 40.0000000000000 60.0000000000000 80.0000000000000 100.0000000000000 120.0000000000000 1 5 9 13 17 21 25 29 measure temperature exact temperature圖 4
r=8、
σ
=0.1%,1 點量測估算邊界溫度與實際邊界溫度
若以組織深度
L=3cm與
L=2.907cm處取 2 點量測,當未來時間
r=8考慮,提高測誤差達
σ
=1%、
% 10 =σ
, 其平均誤差百分比分別為 1.6044243901e-002%與 0.16044243901%,比 1 點量測所得實際溫度輸
入值更接近,如圖 5、圖 6 所示,相關數據如表 1 所示。
0 20 40 60 80 100 120 1 5 9 13 17 21 25 29 33 37 41 45 49 measure temperature exact temperature
圖 5
r=8、
σ
=1%,2 點量測估算邊界溫度與實際邊界溫度
0 20 40 60 80 100 120 1 5 9 13 17 21 25 29 33 37 41 45 49 measure temperature exact temperature圖 6
r=8、
σ
=10%,2 點量測估算邊界溫度與實際邊界溫度
表 1 邊界為溫度輸入估算值與實際值比較表
未來時間( )
r量測點數
誤差( σ %)
估算值與實際值平均誤差百分比(%)
1 1 0
1.2225791317e-006
8 1 0.1
3.84790859671
8 2 1
1.6044243901e-002
8 2 10
0.16044243901
範例 2:邊界為震盪 cos 型熱流輸入
在表皮
x=0處,輸入震盪 cos 型熱源
f(θ
)=900+500cos(0.01t),在組織內部
x =3cm處,為絕熱邊界
0 x = ∂ ∂θ
,位於組織深度
L=3cm處取 1 點量測,若以所得直接解模擬反求邊界熱流之逆向解,在未來時間(
r=1)
且無誤差(
σ
=0)的理想狀態下,得估算邊界熱流結果與實際邊界熱流輸入相符,如圖 7 所示。當未來時間
r=8並考慮量測誤差
σ
=0.001%, 其平均誤差百分比為 7.98790091308%與實際熱流輸入值非常接近,如圖 8 所
示。
0 200 400 600 800 1000 1200 1400 1600 1 4 7 10 13 16 19 22 25 28measure heat flux exact heat flux
圖 7
r=1、
σ
=0,1 點量測估算邊界熱流與實際邊界熱流
0 200 400 600 800 1000 1200 1400 1600 1800 1 4 7 10 13 16 19 22 25 28measure heat flux exact heat flux
若以組織深度
L=3cm與
L=2.907cm處取 2 點量測,當未來時間
r=8考慮量提高測誤差達
σ
=1%、
% 10 =σ
,其平均誤差百分比為 0.15499177375%與 1.28128974394%,明顯比 1 點量測所得實際熱流輸入值更
接近,如圖 9、圖 10 所示,相關數據如表 2 所示。
0 200 400 600 800 1000 1200 1400 1600 1 4 7 10 13 16 19 22 25 28measure heat flux exact heat flux
圖 9
r=8、
σ
=1%,2 點量測估算邊界熱流與實際邊界熱流
-2000 -1000 0 1000 2000 3000 4000 5000 1 4 7 10 13 16 19 22 25 28measure heat flux exact heat flux
表 2 邊界為溫度輸入估算值與實際值比較表
未來時間( )
r量測點數
誤差( σ %)
估算值與實際值平均誤差百分比(%)
1 1 0
1.4662778507e-005
8 1 0.001
7.98790091308
8 2 1
0.15499177375
8 2 10
1.28128974394
4. 結 論
在生物醫學上,當組織溫度達 41℃~44℃[37-38 ]以上,將破壞組織細胞的生長達到組織細胞的壞死,因
此模擬邊界溫度與熱流輸入時,均須以實際需求考量,考慮若以未知邊界為熱治療時,在組織內部鄰近
處應達到溫度約 41℃~44℃,方能達到腫瘤細胞壞死,但又必須避免溫度過高而傷害過多的體內正
常細胞組織。
cm 3 x=以上述二範例模擬結果顯示,固定未來時間
r=8,1 個量測點的估算效果已經很好,但 2 個量測點估算
結果還比 1 個量測點好;2 點量測時,其誤差值
σ 可取高達 10 倍至 100 倍於 1 點量測,且平均誤差比也比 1
點量測小。模擬結果也顯示,邊界為熱流輸入其敏感度較溫度輸入為大,故當邊界為熱流輸入誤差值須比溫
度輸入小。
由上述二範例可確認本文所提出逆運算法可以精確反求未知邊界溫度或熱流輸入,未來可提供生物醫學
在非侵入式的熱診斷、熱治療、低溫手術與熱參數估算臨床相關資料。
參考文獻
[1] Ching-yu Yang,“A linear inverse model for the temperature-dependent thermal conductivity determination in one-dimensional problems,”Applied Mathermatical Modelling, vol.22, pp.1-9, 1998.
[2] Ching-yu Yang,“ Solving the two-dimensional inverse heat source problem through the linear least-squares error method,” Heat and Mass transfer, vol.41, No.2, pp. 393-398, 1998.
[3] Ching-yu Yang,,“A sequential method to estimate the strength of the heat source based on symbolic computation,” Heat and Mass transfer, vol.41, No.14, pp.2245-2252.
[4] Ching-yu Yang,,“The determination of two heat sources in an inverse heat conduction problem,” Heat and Mass transfer, vol.42, pp. 345-356, 1999.
[5] Ching-yu Yang,, “Estimation of the temperature-dependent thermal conductivity in inverse heat conduction problems,” Applied Mathematics Modelling, vol. 23, pp. 469-478, 1999.
[6] Ching-yu Yang,, “Determination of the temperature dependent thermophysical properties from temperature responses measured at medium's boundaries,” Heat and Mass transfer, Vol.43, pp. 1261-1270, 2000.
[7] Ching-yu Yang,, “The determination of two moving heat sources in two-dimensional inverse heat problem,” Applied Mathematics Modelling, vol. 30, pp. 278-292, 2006.
[8] Shidfar, A. and Pourgholi , R., “Numerical approximation of solution of an inverse heat conducton problem based on Legendre polynomials,” Applied Mathermatics and Computeration, vol. 175, pp.1366-1374, 2006.
[9] Shih-Yu Shen ,“A numerical study of inverse heat conduction problems,” Computers & Mathematics with applications, vol. 38, pp. 173-188, 1999.
[10] Park, H. M. and Yoon, T. Y., “Solutions of the inverse radiation problem using a conjugate gradient method,” Heat and Mass transfer, vol. 43, pp. 1767-1776, 2000.
[11] Prud'homme, M. and Hung Nguyen, T., “On the iterative regularization of inverse heat conduction problems by conjugate gradient method,” Heat and Mass transfer, vol. 25, No. 7, pp.999-1008, 1998.
[12] Liu, J. and Xu, L. X., “Boundary information based diagnostics on the thermalstates of biological bodies,” Int. J. Heat Mass Transfer, vol. 43, pp. 2827-2839, 2000.
[13] Nelson, D. A., Nelson, M. T., Walters, T. J. and Mason, P. A., “Skin heating effects of millimeter-wave irradiation—thermal modeling results,” IEEE Trans. Microwave Theory Techniques, vol. 48, pp.2111–2120, 2000.
[14] Deng, Z. S. and Liu, J., “Analytical study on bioheat transfer problems with spatial or transient heating on skin surface or inside biological bodies,” J. Biomed. Engrg., vol. 124, pp. 638–649, 2002.
[15] Hu, L., Gupta, A., Gore, J. P. and Xu, L. X., “Effect of forced convection on the skin thermal expression of breast cancer,” J. Biomed. Engrg., vol. 126, pp.204–211, 2004.
[16] Rabin, Y. and Shitzer, A., “Numerical solution of the multidimensional freezing problem during cryosurgery,” J. Biomed. Engrg., vol. 120, pp.32–37, 1998.
[17] Zhang, J., Sandison, G. A., Murthy, J. Y. and Xu, L. X., “Numerical simulation for heat transfer in prostate cancer cryosurgery,” J. Biomed. Engrg., vol. 127 , pp.279–294, 2005.
[18] Chato, J. C., “Measurement of thermal properties of growing tumors,” Ann. N.Y. Acad. Sci., vol. 355, pp.67–85, 1980. [19] Chen, M. M., Holmes, K. R. and Rupinskas, V., “Pulse-decay method for measuring the thermal conductivity of living
tissues,” J. Biomed. Engrg., vol.103, pp.253–260, 1981.
[20] Valvano, J. W., Allen, J. T., and Bowman, H. F., “The simultaneous measurement of thermal conductivity, thermal diffusivity, and perfusion in small volumes of tissue,” J. Biomed. Engrg., vol.106, pp.192–197, 1984.
[21] Liu, J. and Xu, L. X., “Estimation of blood perfusion using phase shift in temperature response to sinusoidal heating at the skin surface, “ IEEE Trans. Microwave Theory Techniques, vol. 46, pp.1037–1043, 1999.
[22] Liu, J., Xu, L. X., “Boundary information based diagnostics on the thermal states of biological bodies,” Int. J. Heat Mass Transfer, vol. 43, pp.2827–2839, 2000.
[23] Babbs, C. F., Fearnot, N. E., Marchosky, J. A., Moran, C. J., Jones, J. T. and Plantenga, T. D., “Theoretical basis for controlling minimal tumor temperature during interstitial conductive heat therapy,” IEEE Trans. Bio-Med. Eng., vol.37, pp.662–672, 1990.
[24] Young, L. A. and Boehm, R. F., “Finite di5erence heat transfer analysis of a percutaneous transluminal microwave angioplasty system,” ASME J. Biomech. Eng., vol. 115, pp.441–446, 1993.
[25] Haider, S. A., Cetas, T. C. and Roemer, R. B., “Temperature distribution in tissues from a regular array of hot source implants,” IEEE Trans. Bio-Med. Eng., vol. 40, pp.408–417 , 1993.
[26] Meaney, P. M., Clarke, R. L., Ter Haar, G. R. and Rivens, I. H., “3-D Dnite-element model for computation of temperature proDles and regions of thermal damage during focused ultrasound surgery exposures,” Ultrasound Med. Biol., vol. 24, pp.1489–1499, 1998.
[27] Chatterjee, I., Adams, R. E., “Finite element thermal modeling of the human body under hyperthermia treatment for cancer,” Int. J. Comput. Appl. Technol., vol. 7, pp. 151–159, 1994.
[28] Chan, C. L., “Boundary element method analysis for the bioheat transfer equation,” ASME J. Biomech. Eng., vol. 114, pp.358–365, 1992.
[29] Ren, Z. P., Liu, J. and Wang, C. C., “Boundary element method (BEM) for solving normal or inverse bio-heat transfer problem of biological bodies with complex shape,” J. Therm. Sci., vol. 4, pp.117–124, 1995.
[30] Cheng-Hung Huang, Shao-Pei Wang, “A three-dimensional inverse heat conduction problem in estimating surface heat flux by conjugate gradient method,” Heat and mass transfer, Vol. 42, pp. 3387-3403, 1999.
[31] Zhong-Shan Deng, Jing Liu, “Mathematical modeling of temperature mapping over skin surface and its implementation in thermal disease diagnostics,” Computers in biology and Medicine, vol. 34, pp. 495-521, 2004.
[32] Paruch, M. and Majchrzak, E., “Identification of tumor region parameters using evolutionary algorithm and multiple reciprocity boundary element method,” Engineering Applications of Artificial Intelligence, vol. 20, pp.647-655, 2007. [33] Pennes, H. H., “Analysis of tissue and arterial blood temperatures in the resting human forearm,” J. Appl. Physiol., vol.1,
pp.93–122., 1948.
[34] Beck, J. V., Blackwell, B. and Clair, C. R., st., “Inverse Heat Conduction: Ill-posed problem,”Wiley, New York, 1985. [35] Wensheng Shen , Jun Zhang, Fuqian Yang, “Modeling and Numerical Simulation of Bioheat Transfer and Biomechanics in
Soft Tissue,” Mathematical and Computer Modelling, vol. 41, pp. 1251-1265, 2005.
[36] Liu J, and Xu LX., “Boundary information based diagnostics on the thermal states of biological bodies.” Int J Heat Mass Transfer, Vol. 43, pp.2827–39, .2000.
[37] Rajan Dua, and Suman Chakraborty, “A novel modeling and simulation technique of photo-thermal interactions between lasers and living biological tissues undergoing multiple changes in phase,” Computers in Biology and Medicine, vol. 35, pp.447-462, 2005.
[38] Cedric Thiebaut, and Denis Lemonnier, “Three-dimensional modeling and optimisation of thermal fields incluced in a human body during hyperthermia,” International Journal of Thermal Sciences, vol.41, pp. 500-508, 2002.