数值抛物型方程的积分插值法
我们先从抛物型方程的基本形式开始。一维抛物型方程的标准形式是热传导方程:
\[\frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2} \]
其中 \(u\) 是未知函数(如温度),\(\alpha > 0\) 是扩散系数。数值求解此方程的关键在于对空间二阶导数进行离散。积分插值法的核心思想是:对守恒形式的方程在控制体上积分,并利用插值近似界面通量。
第一步,我们将方程改写为守恒形式。注意到 \(\frac{\partial^2 u}{\partial x^2} = \frac{\partial}{\partial x} \left( \frac{\partial u}{\partial x} \right)\),因此热方程可视为一种“通量”的守恒律:
\[\frac{\partial u}{\partial t} = \frac{\partial}{\partial x} \left( \alpha \frac{\partial u}{\partial x} \right) \]
现在,我们在空间计算域上划分控制体。对于一维问题,第 \(i\) 个控制体是区间 \([x_{i-1/2}, x_{i+1/2}]\),其中心点为 \(x_i\),步长为 \(\Delta x\)。
第二步,对控制体进行积分。将上述方程在第 \(i\) 个控制体上对空间积分,并利用散度定理:
\[\int_{x_{i-1/2}}^{x_{i+1/2}} \frac{\partial u}{\partial t} \, dx = \alpha \left( \frac{\partial u}{\partial x} \bigg|_{x_{i+1/2}} - \frac{\partial u}{\partial x} \bigg|_{x_{i-1/2}} \right) \]
左边项可近似为 \(\Delta x \frac{du_i}{dt}\),其中 \(u_i\) 是控制体中心 \(u\) 的平均值。于是我们得到半离散格式:
\[\frac{du_i}{dt} = \frac{\alpha}{\Delta x} \left( \left. \frac{\partial u}{\partial x} \right|_{i+1/2} - \left. \frac{\partial u}{\partial x} \right|_{i-1/2} \right) \]
现在,关键步骤是如何通过插值近似界面 \(x_{i+1/2}\) 处的导数 \(\frac{\partial u}{\partial x}\)。
第三步,采用插值法构造界面导数。积分插值法通常使用分段多项式插值。例如,采用线性插值:假设在控制体内 \(u\) 是线性分布的,即 \(u(x) = u_i + s_i (x - x_i)\),其中 \(s_i\) 是斜率。那么界面 \(x_{i+1/2}\) 处的导数需要由相邻节点值共同决定。一种常见做法是取相邻控制体的线性重构导数的平均,或直接使用中心差商:
\[\left. \frac{\partial u}{\partial x} \right|_{i+1/2} \approx \frac{u_{i+1} - u_i}{\Delta x} \]
代入半离散格式,得到:
\[\frac{du_i}{dt} = \alpha \frac{u_{i+1} - 2u_i + u_{i-1}}{\Delta x^2} \]
这正是熟知的中心差分离散。但积分插值法的优势在于可以灵活采用更高阶的插值多项式来提高精度。
第四步,讨论时间离散。上述半离散系统是一个常微分方程组,可采用各种时间积分方法,如显式欧拉法:
\[u_i^{n+1} = u_i^n + \frac{\alpha \Delta t}{\Delta x^2} (u_{i+1}^n - 2u_i^n + u_{i-1}^n) \]
稳定性要求 \(\frac{\alpha \Delta t}{\Delta x^2} \leq \frac{1}{2}\)。若采用隐式方法(如Crank-Nicolson格式),则无此限制,但需求解线性方程组。
第五步,分析方法的性质。积分插值法本质是一种有限体积法,天然保证守恒性。其精度取决于界面通量(此处是导数)的插值精度。线性插值得到二阶精度。若采用二次或三次插值(如利用更多相邻节点信息),则可达到更高阶精度。此外,该方法易于处理变系数问题 \(\alpha(x)\),只需在界面处取合适的平均(如调和平均)。
总结:积分插值法通过局部积分将方程转化为守恒形式,再利用插值多项式近似界面通量,结合时间离散,为抛物型方程提供了一种灵活、守恒且易于推广到高维和非均匀网格的数值解法。