一类二维有界域上稳态热传导方程侧边值问题的计算方法_刘唐伟.pdf

第 卷 第 期 年 月东 华 理 工 大 学 学 报（自 然 科 学 版）（）.收稿日期：基金项目：国家自然科学基金项目（，）作者简介：刘唐伟（），男，博士，教授，主要从事应用数学及地热学研究。：一类二维有界域上稳态热传导方程侧边值问题的计算方法刘唐伟，钟小雨，欧阳旺林，唐阿敏（东华理工大学 理学院，江西 南昌）摘 要：为了由二维区域部分边界地温观测数据推算区域内部地温场，建立一类具有非齐次边界条件稳态热传导方程侧边值问题的数学模型并进行数值求解。该数学模型是一类典型的不适定问题。利用齐次化原理，将问题中的非齐次边界条件齐次化。通过分离变量法将非齐次方程转化成第一类 积分方程。利用正则化方法求解不适定积分方程，得到未知边界条件，进一步求得泊松方程侧边值问题的数值解。依据所提出的数值方法，设计了三个数值算例，可由矩形域三条边界上的温度数据及其中一条边上的地温梯度数据，计算矩形区域上的地温场。本成果对地热资源勘探开发和岩石圈热结构研究中地温场的数值模拟具有参考意义。关键词：地温场；热传导方程；侧边值问题；积分方程；正则化方法中图分类号：；文献标志码：文章编号：（）刘唐伟，钟小雨，欧阳旺林，等，一类二维有界域上稳态热传导方程侧边值问题的计算方法 东华理工大学学报（自然科学版），（）：，（），（）：在地温场的研究中，经常遇到如何由低维观测数据推算高维温度场的问题。例如，已知两个观测点的垂向温度分布数据及两观测点间的地表温度和地温梯度数据，需要推算相应区域的二维地温场的温度分布。在一定条件下，该问题的数学模型为二维泊松方程侧边值模型。假设研究区为矩形域，具有分段光滑边界，（，），具体示意图见图。和 表示垂向边界，表示地表边界，控制方程及定解条件如式（），其中（，）、（）、（）、（）、（）均为已知函数，边界条件（，）（）未知，需求解未知温度函数（，），此类反问题可看作二维矩形域上第一类边界条件泊松方程侧边值问题。（，）（，），（，）（，）（），（，）（），（，）（），（，）（），|（）式中，。广义上来看，二维泊松方程侧边值问题是指在二维区域部分边界数据已知，而另外部分边界数据未知的问题，它包含二维有界域上的问题和无界域上的问题。泊松方程侧边值问题在无损探伤（，；，）、医学成像（，）、地质学（徐有缘等，）等领域都有涉及，是典型的不适定问题，某些边界上观测数据微小的扰动将会引起解发生巨大的改变。前人提出了各种方法去解决各种柯西问题的不适定性。如：正则化方法利用 小波求得半平面中的 方程柯西问题的正则化解（，）；利用拟可逆方法构造控制方程图 矩形域示意图 及边界条件，采用卡尔曼估计推导误差估计，再运用有限差分方法得到 方程侧边值问题的数值解（，）。对于二维和三维的 方程柯西问题，矩方法利用格林公式将 方程柯西问题转化为矩问题，通过求解矩问题来获得区域边界上的值，进而获得数值解（王泽文等，），该方法主要适用于方程源项为多项式函数的问题。基本解方法利用边界控制技术进行边界处理后，再结合基本解，用正则化方法获得数值结果（，；曹瑞华，）。变分正则化方法是利用格林函数求解 方程柯西问题，根据“曲线”准则选择最优参数，得到问题结果（，）。正则化方法通过构造正则解解决了在无限条状区域的带有非齐次 条件的 方程柯西问题，给出了近似解和精确解的误差估计，由偏差原理得到近似解的后验误差估计（曹笑笑等，）。人工神经网络方法利用多层网络作为近似，提出了一种非网格离散方法来解决柯西问题，该方法的优势在于更容易扩展到高维（，）。对具有混合边界条件的椭圆方程逆问题，最近有学者构造和验证相应的变分源条件，研究了吉洪诺夫正则化方法的收敛性，基于两种新的对数型稳定性，导出了逆问题求解范数的收敛性和收敛速度（，）。从已有文献可知，有较多学者对不同的泊松方程侧边值问题进行了研究，获得了不少理论研究成果，但二维有界域上侧边值问题的计算方法研究相对较少，简洁易行的数值算法并不多见。考虑数值求解矩形域上侧边值模型，表达式如式（），利用叠加原理和齐次化方法，将非齐次问题的边界条件齐次化，得到具有部分齐次边界的泊松方程侧边值问题；参考张宏武等（）方法，利用分离变量法求出齐次边界问题的通解，将所求问题转换成第一类 积分方程问题，再利用 正则化方法求解方程，进行解的存在唯一性分析和误差估计，并开展数值模拟实验。侧边值问题转化为积分方程 边界条件齐次化由泊松方程的边界条件齐次化原理，可设：（，）（）（）（）令（，）（，）（，），则可得反问题模型：（，），（，）（，），（，），（，）（）（，），（，）（）（，），|（）求解式（）可得（，）在边界 上的值（，），从而矩形域的四个边界条件均已知，通过正演计算得（，），即可得式（）的解（，）。泊松方程转化为积分方程下面对式（）进行分析，利用分离变量法，先设（，）（）为待定函数，则（，）（）（，）未知待定，令：（，）（）（）将式（）代入式（）中的稳态方程，在等式左右两边分别乘且对 积分得：（）()（）（）（），（）（）式中，（）（，），（）（，），（）（，）。根据二阶常系数线性微分方程理论，式（）的东 华 理 工 大 学 学 报（自 然 科 学 版）年齐次通解和特解都可以解出。根据（）可求得特解，设特解为?（），则：（）?（）（）代入式（）中边界条件得?（）?（）（）解式（）即可得，从而可得：（，）?（）（）由（，）（）（，），结合式（）知?（）?（）?（）?（）?（）?（）（）（，）（）继而可得（）（，）（）（）（）式中，为双曲正弦函数，（）（）（，）?（）?（）?（）?（）（）记式（）的左端（）（，）（）关于无穷级数，有如下引理。引理：级数 在，上收敛。证：因（）（，）属于空间，故（）（，）有上界，则 （）（）（）（）记式（）右端的无穷级数为，利用比式判别法对 进行收敛性判断。（）（）（）（）（）（）由，得，故 收敛，从而式（）的左端 是收敛的。进一步可得以下定理。定理：式（）的解在空间，上是唯一的。证：设式（）在空间，上有两个解（），（），则（）（，）（）（）（）（）（，）（）（）（）式（）、（）两式相减，得（）（）（）（）由函数系，在区间，上的正交完备性可知：第 期刘唐伟等：一类二维有界域上稳态热传导方程侧边值问题的计算方法（）（）（）（）从而（）（）（）易知（）（），故（）（），唯一性得证，定理成立。因此式（）可以转化成关于（）的第一类 积分方程，即（，）（）（，）（）（）式中，核函数（，）（）（）积分方程的正则化求解 积分方程的离散格式下面对式（）进行离散求解计算。记（，）（）（，）（）（，）。将式（）中的积分区间，划分为 等份，记 为步长，对变量，在区间，上离散取值，利用复化梯形数值求积公式：（，）（）（，）（）（，）（）（，）（）（）可得式（）的数值离散格式 （）式中，（，）（，）（，）（，）（，）（，）（，）（，）（，）|，（），（），（），（）（，），。由 节理论分析，结合张鹤琼等（）结论，知第一类 积分方程的解（）存在唯一。正则化求解与参数选取公式对式（）采用 正则化方法求解，设正则化参数为，正则化求解公式为：（），表示单位矩阵。利用后验参数选取方法中的“曲线”准则选取正则化参数。“曲线”准则是利用 尺度去刻画 和 的一条 形曲线，参考曲线的拐点（隅角）进行参数选择。参数 为 曲线在 尺度下的最大曲率（，）。令 ，表示 范数，具体表达式形如 。关于参数 的曲率函数定义为：（）（）（）（）大多数情况下，基于“曲线”准则，可取 和的乘积极小时的参数作为正则化参数值（，），即：（）（）式中，。为了便于计算，文中数值实验的正则化参数均采用式（）进行计算选取。根据后验参数选取理论，如果式（）右端扰动为，可得正则化解的误差估计为（张路寅等，）：（）（）有限差分格式和误差估计下面对式（）进行数值求解分析。当求得（，）在边界 上的数值解（）后，相应的控制方程及定解条件可化为：（，），（，）（，）（），（，）（），（，）（），（，）（），|（）利用有限差分法求解式（）。采用五点差分格式，将，分为 等份，分为 等份，沿 轴和 轴的步长分别为和。对矩形域进行网格剖分：，（）东 华 理 工 大 学 学 报（自 然 科 学 版）年在节点（，）处的网函数值（，）用，表示，（，）用，表示。在正则内点（，）处的离散格式为：，（）式中，。差分式（）的解存在唯一。在均匀网格下，即当 时，五点差分格式的截断误差为（李荣华等，）。（）（，）（）（）式中，是方程 的光滑解。数值实验设计了 个数值算例，其解析解表达式类型分别是三角函数、多项式函数、指数函数。采用加高斯白噪声的边界温度观测数据模拟计算，结合式（），扰动数据具体表达式为：（，）（）（）（）（）（）式中，（，）表示一个由正态分布随机数组成的（）向量。例设地表温度数据（，）（），垂向边界温度（，），（，），地表地温梯度（，），地表边界条件采用带噪声数据，所加噪声为高斯白噪声，精确解（，）（），。原问题齐次化边界条件后解的表达式为（，）（），正则化参数取 。例设地表温度数据（，），垂向边界温度（，），（，），地表地温梯度数据（，），地表边界条件采用带噪声数据，所加噪声为高斯白噪声，精确解（，），。原问题齐次化边界后的解表达式为（，），正则化参数取 。例设地表边界温度数据（，），垂向边界温度数据（，）（），和（，），地表地温梯度数据（，）（），地表边界条件采用带噪声数据，所加噪声为高斯白噪声，精确解（，）（），。原问题齐次化边界后的解的表达式为（，）（）（），正则化参数取 。图 是例 求解过程中正则化参数选取所用的 曲线。根据 节中的讨论和分析得到最优参数的位置，如图 中圆圈所示，进而得到正则化参数值，例 和例 中的正则化参数选取类似可得。具体数值解和精确解的对比如图 至图。为了更精确地分析和对比误差水平，计算了地温场的最大绝对误差和平均相对误差（表）。在计算过程中，有些离散点的平均相对误差会出现 的情况，这时计算所得为非数值，文中用 代替。由图 至图，以及表，可知正则化反演计算方法精度较高。图 例 参数选取的 曲线 表 （，）反演计算的误差 （，）数值算例最大绝对误差平均相对误差例 例 例 注：最大绝对误差为各离散点处计算值和精确值的差的绝对值的最大值；平均相对误差为各离散点处计算值和精确值的差除以精确值所得平均值。结论利用齐次化原理和分离变量法，求解了矩形区域上的第一类边界条件泊松方程侧边值问题，所提出的方法适用性较广，算法简洁且便于数值实现，数值结果显示此方法在求解精度和稳定性上效果第 期刘唐伟等：一类二维有界域上稳态热传导方程侧边值问题的计算方法图 例 的数值解和精确解对比图 图 例 的数值解和精确解对比图 良好。此计算方法可以用于由一维地温数据推算二维地温场，适用于均匀介质地层，对于非均匀介质区域地温场，可以考虑分层计算，将某一地层的温度场计算出来以后，获得相应边界温度和地温梯度，应用于下一层二维温度场的计算，依次算出多层介质地层各二维区域的地温场。参 考 文 献曹瑞华，多联通区域中的拉普拉斯方程柯西问题的一种数值计算方法 山西师范大学学报（自然科学版），（）：曹笑笑，毛东玲，程强，等，带有非齐次 条件的 方程 问题的一种傅里叶正则化方法 湖北大学学报（自然科学版），（）：李荣华，刘播，微分方程数值解法 版 北京：高等教育出版社：王泽文，刘唐伟，徐定华，方程 问题的一种数值解法 华东地质学院学报，（）：王泽文，钱海忠，徐定华，高维 方程 问题的一种数值解法 宁夏大学学报（自然科学版），（）：徐有缘，周健，缪俊发，方程若干问题的解及其在渗透破坏中的应用 上海地质（）：东 华 理 工 大 学 学 报（自 然 科 学 版）年图 例 的数值解和精确解对比图 张鹤琼，罗文陶，积分方程的解的存在性和唯一性 科技资讯（）：张宏武，朱睦正，方程 问题的修正 正则化方法 河西学院学报，（）：张路寅，张玉海，钱坤明，关于不适定问题的迭代 正则化方法 山东大学学报（理学版），（）：，：，（）：，：，（）：，（）：，（）：，（）：，（）：，：，（）：，：，：（），（）：第 期刘唐伟等：一类二维有界域上稳态热传导方程侧边值问题的计算方法 ，（，）：，：；东 华 理 工 大 学 学 报（自 然 科 学 版）年