熱方程Robin系數(shù)反演的優(yōu)化方法

王兵賢,杜海清

(淮陰師范學(xué)院 數(shù)學(xué)科學(xué)學(xué)院，江蘇 淮安 223300)

考慮熱傳導(dǎo)問題

(1)

這是一個(gè)帶混合邊界條件的熱傳導(dǎo)系統(tǒng)，其中σ(t)為Robin系數(shù).當(dāng)σ(t)，φ(x)已知且滿足一定光滑性條件時(shí)，問題(1)的解存在而且唯一[1].

當(dāng)給定附加條件

u(x0,t)=h(t),t>0，

(2)

其中：x0∈(0,π)，或者給定h(t)的測(cè)量值hδ(t)滿足

‖hδ(t)-h(t)‖<δ,t>0，

(3)

反演邊界上Robin系數(shù)為σ(t)，這是一個(gè)不適定問題[2]，也是一個(gè)非常有意義的研究課題.Robin系數(shù)度量了邊界材料對(duì)內(nèi)部擴(kuò)散的阻尼作用，當(dāng)該系數(shù)與傳導(dǎo)介質(zhì)邊界上的位置有關(guān)時(shí)，即跟空間變量有關(guān)[2].Liu等[3]研究二維區(qū)域的熱傳導(dǎo)系統(tǒng)，在由邊界上的非局部測(cè)量數(shù)據(jù)反演邊界的Robin系數(shù)時(shí)，發(fā)現(xiàn)該系數(shù)與局部的空氣流動(dòng)以及外部條件有關(guān)，即跟時(shí)間變量有關(guān).對(duì)于Robin系數(shù)反演的研究成果見文獻(xiàn)[4-9].

定義容許集

M:={σ(t)∈C[0,T],0<σ1≤σ(t)≤σ2}，

(4)

則反問題(1)～(3)有唯一性結(jié)果[6]，即引理1.

引理1假設(shè)φ(x)>0且φ(x)不恒等于0，σi(t)(i=1,2)∈M，對(duì)于u[σi,φ](x,t)∈C2,1(Q)，如果在[0,T]上，u[σ1,φ](x0,t)=u[σ2,φ](x0,t)，則σ1(t)=σ2(t)∈C[0,T].

1 優(yōu)化問題

基于反問題(1)～(3)解的唯一性，將反問題轉(zhuǎn)化為優(yōu)化問題.文獻(xiàn)[10]對(duì)于近似測(cè)量數(shù)據(jù)hδ(t)，定義正則化泛函

(5)

其中：α>0為正則化參數(shù)，Θ(σ)為正則化罰項(xiàng).

(6)

下面的定理說明在Θ(σ)意義下，目標(biāo)泛函(5)極小元的存在性.

(7)

下面假設(shè){σn:n∈N}?M為Jα(σ*)的極小化序列，即

因此，有

定理得證.

2 目標(biāo)泛函的梯度

定理2目標(biāo)泛函Jα[σ]是Frechét可微的，而且其Frechét導(dǎo)數(shù)為

其中：P(x,t)滿足問題

(8)

(9)

對(duì)于極限

(10)

(11)

(12)

在問題(12)中方程的兩側(cè)同時(shí)乘以函數(shù)P(x,t)，在區(qū)域[x0,π]×[0,T]上積分，并運(yùn)用分部積分化簡(jiǎn)得

由問題(8)，(11)，(12)中的初邊值條件得Jα[σ]是Frechét可微的，且其導(dǎo)數(shù)為

其中：P(x,t)滿足初邊值問題(8)，即定理得證.

σJ=P(π,t)u(π,t)-P(x0,t)u(x0,t)+α(σ(t)-σ″(t)).

(13)

3 迭代算法

對(duì)于反問題(1)，(2)，設(shè)σ0為初始猜測(cè)函數(shù)，由泛函Jα[σ]的Frechét可微性，可以構(gòu)造如下迭代格式

σk+1=σk+Dkhk，k=0,1,2,…，

(14)

其中：hk搜索因子，這里通過Wolfe線性搜索方法[12]得到；Dk為搜索方向，初始搜索方向取負(fù)梯度方向，第二次以后的方向用文獻(xiàn)[13]中混合的下降算法確定.有

D0=-g0=-σJ|σ0，σ1=σ0+D0h0，

Dk=-gk+ζk-1Dk-1，σk+1=σk+Dkhk,k=1,2,3,….

具體的迭代算法如下：

步驟1 給定σ0，先驗(yàn)選取正則化參數(shù)α，并給定任意小的正數(shù)ε和最大迭代次數(shù)Nmax；

步驟2 數(shù)值求解問題(1)，得到u[σ](xi,tj)，并代入h(t)數(shù)值求解問題(8)，得到P[σk](xi,tj)，計(jì)算出Dk，由式(14)得到σk+1，再計(jì)算出J(σk)，J(σk+1)；

步驟3 判斷J(σk)>J(σk+1)是否成立：

①是，則判斷J(σk)<ε或者k>=Nmax是否成立，如果是，則得到σk+1≈σ*,并輸出，如果否，則回到步驟2；②否，則重新選取式(14)的hk；

步驟4 輸出σ*的近似結(jié)果.

4 數(shù)值模擬

hδ(t)=u(x0,t)[1+δ·randn(x)]，

(15)

時(shí)的反演結(jié)果;時(shí)的反演結(jié)果.圖1 反演結(jié)果

需要說明的是，在計(jì)算目標(biāo)泛函梯度的時(shí)候，其中第二個(gè)罰項(xiàng)中存在σ(t)的二階導(dǎo)數(shù)，給數(shù)值計(jì)算帶來了很大困難，當(dāng)選取試驗(yàn)中σ(t)是光滑函數(shù)時(shí)，用離散近似計(jì)算的結(jié)果相對(duì)很好，即先將數(shù)據(jù)用插值方法計(jì)算，然后再計(jì)算二階導(dǎo)數(shù).如果σ(t)是分段光滑函數(shù)，在間斷點(diǎn)附近的近似結(jié)果則不是很理想，增加正則化罰項(xiàng)以后，雖然迭代次數(shù)增多了，但是在間斷點(diǎn)的重建效果將會(huì)有所改善.