热方程Robin系数反演的优化方法

王兵贤,杜海清

(淮阴师范学院 数学科学学院，江苏 淮安 223300)

考虑热传导问题

(1)

这是一个带混合边界条件的热传导系统，其中σ(t)为Robin系数.当σ(t)，φ(x)已知且满足一定光滑性条件时，问题(1)的解存在而且唯一[1].

当给定附加条件

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

(2)

其中：x0∈(0,π)，或者给定h(t)的测量值hδ(t)满足

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

(3)

反演边界上Robin系数为σ(t)，这是一个不适定问题[2]，也是一个非常有意义的研究课题.Robin系数度量了边界材料对内部扩散的阻尼作用，当该系数与传导介质边界上的位置有关时，即跟空间变量有关[2].Liu等[3]研究二维区域的热传导系统，在由边界上的非局部测量数据反演边界的Robin系数时，发现该系数与局部的空气流动以及外部条件有关，即跟时间变量有关.对于Robin系数反演的研究成果见文献[4-9].

定义容许集

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

(4)

则反问题(1)～(3)有唯一性结果[6]，即引理1.

引理1假设φ(x)>0且φ(x)不恒等于0，σi(t)(i=1,2)∈M，对于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 优化问题

基于反问题(1)～(3)解的唯一性，将反问题转化为优化问题.文献[10]对于近似测量数据hδ(t)，定义正则化泛函

(5)

其中：α>0为正则化参数，Θ(σ)为正则化罚项.

(6)

下面的定理说明在Θ(σ)意义下，目标泛函(5)极小元的存在性.

(7)

下面假设{σn:n∈N}⊂M为Jα(σ*)的极小化序列，即

因此，有

定理得证.

2 目标泛函的梯度

定理2目标泛函Jα[σ]是Frechét可微的，而且其Frechét导数为

其中：P(x,t)满足问题

(8)

(9)

对于极限

(10)

(11)

(12)

在问题(12)中方程的两侧同时乘以函数P(x,t)，在区域[x0,π]×[0,T]上积分，并运用分部积分化简得

由问题(8)，(11)，(12)中的初边值条件得Jα[σ]是Frechét可微的，且其导数为

其中：P(x,t)满足初边值问题(8)，即定理得证.

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

(13)

3 迭代算法

对于反问题(1)，(2)，设σ0为初始猜测函数，由泛函Jα[σ]的Frechét可微性，可以构造如下迭代格式

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

(14)

其中：hk搜索因子，这里通过Wolfe线性搜索方法[12]得到；Dk为搜索方向，初始搜索方向取负梯度方向，第二次以后的方向用文献[13]中混合的下降算法确定.有

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

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

具体的迭代算法如下：

步骤1 给定σ0，先验选取正则化参数α，并给定任意小的正数ε和最大迭代次数Nmax；

步骤2 数值求解问题(1)，得到u[σ](xi,tj)，并代入h(t)数值求解问题(8)，得到P[σk](xi,tj)，计算出Dk，由式(14)得到σk+1，再计算出J(σk)，J(σk+1)；

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

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

步骤4 输出σ*的近似结果.

4 数值模拟

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

(15)

时的反演结果;时的反演结果.图1 反演结果

需要说明的是，在计算目标泛函梯度的时候，其中第二个罚项中存在σ(t)的二阶导数，给数值计算带来了很大困难，当选取试验中σ(t)是光滑函数时，用离散近似计算的结果相对很好，即先将数据用插值方法计算，然后再计算二阶导数.如果σ(t)是分段光滑函数，在间断点附近的近似结果则不是很理想，增加正则化罚项以后，虽然迭代次数增多了，但是在间断点的重建效果将会有所改善.