论文部分内容阅读
【摘要】本文提出把一类椭圆区域上的偏微分方程变换成圆区域上的偏微分方程，然后在极坐标下作网格划分，使所有网格边界点与边界重合，并给出相应的差分方程和解法，从而避免了一般差分解法由边界条件转移问题所产生的误差.
【关键词】差分；方程；解
一、问题的提出
考虑方程：
λ22ux2 2uy2=f（x，y），（x，y）∈D，
u|D=α（x，y），
其中λ≠0，D=（x，y）x2y2 y2<1.（1）
对此类方程的一般差分解法是首先分别沿x，y轴方向取步长为h和k，作两族与坐标平行的直线：
xi=ih，i=0，±1，±2，….
yj=jk，j=0，±1，±2，….
两族直线的交点为网格点.
Dk={（xi，yj）|（xi，yj）∈D，并且（xi 1，yj），（xi-1，yj），（xi，yj 1），（xi，yj-1）∈D}为内点.
Dh={（xi，yj）|（xi，yj）∈D-Dh}为网格边界点.
显然Dh与D不重合.
在Dh上建立差分格式：
Lhuij=f（xi，yj），（xi，yj）∈Dh，
uij=g（xi，yj），（xi，yj）∈Dh.
（Lh为差分算子）于是出现了边界条件从D到Dh转移的问题，带来了误差并逐层传播，影响所得结果的准确性，所以考虑对这类问题做某种变换，同时采取某种特别的解域划分，使网格边界点与边界重合，并给出这种划分下相应的差分方程和解法.
二、建立差分方程
令x=λt，（1）成为
2ut2 2uy2=f1（t，y），（t，y）∈D1，
u|D1=α1（t，y），
其中D1={（t，y）|t2 y2<1}.（2）
令t=rcosθ，y=rsinθ，
即r=t2 y2，θ=arctanyt，（2）成为
2ur2cos2θ ur·sin2θr-2uθr·cosθ·sinθr 2ur2·sin2θ 2uθr·sinθ·cosθr ur·cos2θr-2urθ·sinθ·cosθr uθ·sinθ·cosθr2 2uθ2·sin2θr2 2urθ·sinθ·cosθr 2urθ·sinθ·cosθr-uθ·sinθ·cosθr2 2uθ2·cos2θr2-uθ·sinθ·cosθr2=1r·rrur 1r2·2uθ2=f（r，θ），（r，θ）∈D2，
u|D2=α2（1，θ）=α2（θ）.
其中D2={（r，θ）|0≤r<1，0≤θ≤2π}.（3）
當r=0时，方程（3）的系数无定义，因此，只有在r>0的情形下其解才有意义，为了定出有定义的解，需要补充在r=0处u有界.
沿径向r和角度θ方向取等步长Δr和Δθ作网格划分：
ri=i 12Δr，i=0，1，2，…，I.
θj=jΔθ，j=0，1，2，…，J-1.
Δr=1I 12，Δθ=2πJ.
于是方程（3）的差分方程为
1ri（Δr）2δr（rδruij） 1r2i（Δθ）2δ2θuij=f2（ri，θj），（4）
其中，δ和δθ分别表示r方向和θ方向的中心差分算子.
即，δruij=ui 12，j-ui-12，j；
δθuij=ui，j 12-ui，j-12.
利用循环边界条件u（r，θ）=u（r，θ 2π），具体写出方程（4）在j=0和j=J-1上的差分方程：
j=0，1ri（Δr2）[ri 12ui 1，0-（ri 12 ri-12）ui，0 ri-12ui-1，0] 1r2i（Δθ）2（ui，1-2ui，0 ui，j-1）=f2（ri，0），
i=1，2，…，I-1.（5）
j=J-1，1ri（Δr）2[ri 12ui 1，j-1-（ri 12 ri-12）ui，j-1 ri-12ui-1，j-1] 1r2i（Δθ）2（ui，0-2ui，j-1 ui，j-2）=f2（ri，θj-1），
i=1，2，…，I-1.（6）
当j=1，2，…，J-2时，差分方程为
1ri（Δr）2[ri 12ui 1，j-（ri 12 ri-12）ui，j ri-12ui-1，j] 1r2i（Δθ）2（ui，j 1-2ui，j ui，j-1）=f2（ri，θj），
i=1，2，…，I-1.（7）
令（7）中ri=Δr2，得（r0，θj）处的差分方程：
2（Δr）2（u1，j-u0，j） 4（Δr）2（Δθ）2（u0，j 1-2u0，j-u0，j-1）=f2（r0，θj），
j=0，1，2，…，J-1.（8）
至此建立了I×J个网格点上的差分方程.
在J个网络边界点上u（1，θj）=α2（θj），
j=0，1，2，…，J-1.（9）
至此得到了方程（4）的一个差分格式，并且使得网格边界点与边界点完全重合.
三、求解差分方程
把（5）（6）（7）（8）分别整理，并对每处（r0，θj）记为
αi，jui，j-αi 1，jui 1，j-αi-1，jui-1，j-αi，j 1ui，j 1-αi，j-1ui，j-1=-f2（ri，θj），
i=0，1，2，…，I-1；j=0，1，2，…，J-1.（10） 网格边界点上ui，j由（9）代入.
把D2的内点按如图所示的排列给予序号.
于是可得方程组（10）的矩阵形式
AU=F，（11）
其中A为方程组的系数矩阵，U为uij依序排列的矩阵，F为f2（ri，θj）依序排列的矩阵.
为了证明所得差分方程是可解的，需要用到下述定理.
极值定理：设差分算子
Lhuij=αijuij-∑blkulk，（xl，yk）∈G（xi，yj），
其中G（xi，yj）表示格式中除（xi，yi）之外相邻网格点集合，在Dh∪Dh上，算子的系数满足
aij>0，（xi，yj）∈Dh，
blk>0，（xi，yj）∈G（xi，yj），
dij=aij-∑blk≥0.
uij是定义在Dk上的函数，uij≠常数，如果Lhuij≤0，（xi，yj）∈Dh，则uij不可能在内点取正的最大值：
maxuij≤maxuij.
Dh Dh
如果Lhuij≥0，（xi，yj）∈Dh，則uij不可能在内点取得负的最小值：
minuij≤minuij.
Dh Dh
证明略.
可以验证差分方程（10）满足极值定理，由此可得差分格式的稳定性和收敛性.另外，由极值定理可知，满足差分方程（10）的uij除非在所有网格点上均取同一数值，不会在网格区域的内点取到最大值和最小值.由此立刻可以证明，相应于α2=0情形的差分方程，即在边界点取零值的差分方程，只有全为零的解，这说明相应于线性方程组（11）的齐次方程组只有零解，因此，线性代数方程组（11）必有唯一解，这说明了所得的差分方程恒有唯一解.
由此可以采用直接的方法来求解，但由于在实际计算中为了保证精确度，步长往往取得很小，网格点很多，所得的差分方程是一个高阶线性代数方程，而用直接求解法就会因为受计算机存储量的限制而产生困难.
由于每个内点上的值只依赖周围相邻四点上的值，相应地，线性方程组的系数矩阵有大量零元素，即是一个稀疏阵，采用迭代法求解就可以充分利用这个优点.为此这里给出如下办法求解这个差分方程.
首先任意给定在网格区域内点（xi，yi）上的数量{u（0）ij}及{u（1）i-1}（i=0，1，2，…，I-1）作为解的近似，把这组数值代入迭代式
u（k 1）ij=ωαi 1αiju（k）i 1，j αi-1αiju（k 1）i-1，j αi，j 1αiju（k）i，j 1
αi，j-1αiju（k 1）i，j-1 （1-ω）u（k）ij，（xi，yi）∈Dk，
u（k 1）ij=α2（θ1），（xi，yi）∈Dk
其中ω是一个常数，可在1<ω<2中选定.当k相当大时，{u（k）ij}就给出所要求的近似解.通常对充分大的k，当相邻两次迭代解{u（k-1）ij}，{u（k）ij}间的误差maxi，j|u（k）ij-u（k-1）ij|或1N∑i，j|u（k）ij-u（k-1）ij|，其中N为网格点的总数小于某个预先给定的适当小的控制数ε>0时，就可以结束迭代过程，完成差分方程的求解.
【关键词】差分；方程；解
一、问题的提出
考虑方程：
λ22ux2 2uy2=f（x，y），（x，y）∈D，
u|D=α（x，y），
其中λ≠0，D=（x，y）x2y2 y2<1.（1）
对此类方程的一般差分解法是首先分别沿x，y轴方向取步长为h和k，作两族与坐标平行的直线：
xi=ih，i=0，±1，±2，….
yj=jk，j=0，±1，±2，….
两族直线的交点为网格点.
Dk={（xi，yj）|（xi，yj）∈D，并且（xi 1，yj），（xi-1，yj），（xi，yj 1），（xi，yj-1）∈D}为内点.
Dh={（xi，yj）|（xi，yj）∈D-Dh}为网格边界点.
显然Dh与D不重合.
在Dh上建立差分格式：
Lhuij=f（xi，yj），（xi，yj）∈Dh，
uij=g（xi，yj），（xi，yj）∈Dh.
（Lh为差分算子）于是出现了边界条件从D到Dh转移的问题，带来了误差并逐层传播，影响所得结果的准确性，所以考虑对这类问题做某种变换，同时采取某种特别的解域划分，使网格边界点与边界重合，并给出这种划分下相应的差分方程和解法.
二、建立差分方程
令x=λt，（1）成为
2ut2 2uy2=f1（t，y），（t，y）∈D1，
u|D1=α1（t，y），
其中D1={（t，y）|t2 y2<1}.（2）
令t=rcosθ，y=rsinθ，
即r=t2 y2，θ=arctanyt，（2）成为
2ur2cos2θ ur·sin2θr-2uθr·cosθ·sinθr 2ur2·sin2θ 2uθr·sinθ·cosθr ur·cos2θr-2urθ·sinθ·cosθr uθ·sinθ·cosθr2 2uθ2·sin2θr2 2urθ·sinθ·cosθr 2urθ·sinθ·cosθr-uθ·sinθ·cosθr2 2uθ2·cos2θr2-uθ·sinθ·cosθr2=1r·rrur 1r2·2uθ2=f（r，θ），（r，θ）∈D2，
u|D2=α2（1，θ）=α2（θ）.
其中D2={（r，θ）|0≤r<1，0≤θ≤2π}.（3）
當r=0时，方程（3）的系数无定义，因此，只有在r>0的情形下其解才有意义，为了定出有定义的解，需要补充在r=0处u有界.
沿径向r和角度θ方向取等步长Δr和Δθ作网格划分：
ri=i 12Δr，i=0，1，2，…，I.
θj=jΔθ，j=0，1，2，…，J-1.
Δr=1I 12，Δθ=2πJ.
于是方程（3）的差分方程为
1ri（Δr）2δr（rδruij） 1r2i（Δθ）2δ2θuij=f2（ri，θj），（4）
其中，δ和δθ分别表示r方向和θ方向的中心差分算子.
即，δruij=ui 12，j-ui-12，j；
δθuij=ui，j 12-ui，j-12.
利用循环边界条件u（r，θ）=u（r，θ 2π），具体写出方程（4）在j=0和j=J-1上的差分方程：
j=0，1ri（Δr2）[ri 12ui 1，0-（ri 12 ri-12）ui，0 ri-12ui-1，0] 1r2i（Δθ）2（ui，1-2ui，0 ui，j-1）=f2（ri，0），
i=1，2，…，I-1.（5）
j=J-1，1ri（Δr）2[ri 12ui 1，j-1-（ri 12 ri-12）ui，j-1 ri-12ui-1，j-1] 1r2i（Δθ）2（ui，0-2ui，j-1 ui，j-2）=f2（ri，θj-1），
i=1，2，…，I-1.（6）
当j=1，2，…，J-2时，差分方程为
1ri（Δr）2[ri 12ui 1，j-（ri 12 ri-12）ui，j ri-12ui-1，j] 1r2i（Δθ）2（ui，j 1-2ui，j ui，j-1）=f2（ri，θj），
i=1，2，…，I-1.（7）
令（7）中ri=Δr2，得（r0，θj）处的差分方程：
2（Δr）2（u1，j-u0，j） 4（Δr）2（Δθ）2（u0，j 1-2u0，j-u0，j-1）=f2（r0，θj），
j=0，1，2，…，J-1.（8）
至此建立了I×J个网格点上的差分方程.
在J个网络边界点上u（1，θj）=α2（θj），
j=0，1，2，…，J-1.（9）
至此得到了方程（4）的一个差分格式，并且使得网格边界点与边界点完全重合.
三、求解差分方程
把（5）（6）（7）（8）分别整理，并对每处（r0，θj）记为
αi，jui，j-αi 1，jui 1，j-αi-1，jui-1，j-αi，j 1ui，j 1-αi，j-1ui，j-1=-f2（ri，θj），
i=0，1，2，…，I-1；j=0，1，2，…，J-1.（10） 网格边界点上ui，j由（9）代入.
把D2的内点按如图所示的排列给予序号.
于是可得方程组（10）的矩阵形式
AU=F，（11）
其中A为方程组的系数矩阵，U为uij依序排列的矩阵，F为f2（ri，θj）依序排列的矩阵.
为了证明所得差分方程是可解的，需要用到下述定理.
极值定理：设差分算子
Lhuij=αijuij-∑blkulk，（xl，yk）∈G（xi，yj），
其中G（xi，yj）表示格式中除（xi，yi）之外相邻网格点集合，在Dh∪Dh上，算子的系数满足
aij>0，（xi，yj）∈Dh，
blk>0，（xi，yj）∈G（xi，yj），
dij=aij-∑blk≥0.
uij是定义在Dk上的函数，uij≠常数，如果Lhuij≤0，（xi，yj）∈Dh，则uij不可能在内点取正的最大值：
maxuij≤maxuij.
Dh Dh
如果Lhuij≥0，（xi，yj）∈Dh，則uij不可能在内点取得负的最小值：
minuij≤minuij.
Dh Dh
证明略.
可以验证差分方程（10）满足极值定理，由此可得差分格式的稳定性和收敛性.另外，由极值定理可知，满足差分方程（10）的uij除非在所有网格点上均取同一数值，不会在网格区域的内点取到最大值和最小值.由此立刻可以证明，相应于α2=0情形的差分方程，即在边界点取零值的差分方程，只有全为零的解，这说明相应于线性方程组（11）的齐次方程组只有零解，因此，线性代数方程组（11）必有唯一解，这说明了所得的差分方程恒有唯一解.
由此可以采用直接的方法来求解，但由于在实际计算中为了保证精确度，步长往往取得很小，网格点很多，所得的差分方程是一个高阶线性代数方程，而用直接求解法就会因为受计算机存储量的限制而产生困难.
由于每个内点上的值只依赖周围相邻四点上的值，相应地，线性方程组的系数矩阵有大量零元素，即是一个稀疏阵，采用迭代法求解就可以充分利用这个优点.为此这里给出如下办法求解这个差分方程.
首先任意给定在网格区域内点（xi，yi）上的数量{u（0）ij}及{u（1）i-1}（i=0，1，2，…，I-1）作为解的近似，把这组数值代入迭代式
u（k 1）ij=ωαi 1αiju（k）i 1，j αi-1αiju（k 1）i-1，j αi，j 1αiju（k）i，j 1
αi，j-1αiju（k 1）i，j-1 （1-ω）u（k）ij，（xi，yi）∈Dk，
u（k 1）ij=α2（θ1），（xi，yi）∈Dk
其中ω是一个常数，可在1<ω<2中选定.当k相当大时，{u（k）ij}就给出所要求的近似解.通常对充分大的k，当相邻两次迭代解{u（k-1）ij}，{u（k）ij}间的误差maxi，j|u（k）ij-u（k-1）ij|或1N∑i，j|u（k）ij-u（k-1）ij|，其中N为网格点的总数小于某个预先给定的适当小的控制数ε>0时，就可以结束迭代过程，完成差分方程的求解.