🎯要点

🎯二维热传导二阶偏微分方程 | 🎯调和函数和几何图曲率 | 🎯解潮汐波动方程 | 🎯解静止基态旋转球体流体运动函数 | 🎯水文空间插值 | 🎯流体流动模拟求解器 | 🎯随机算法解二维高压电容器电势 | 🎯解空心导电圆柱体交替电势 | 🎯稠密矩阵椭圆微分快速算法

📜拉普拉斯方程用例：Python火焰锋动力学和浅水表面波浪偏微分方程

📜泊松方程用例：Python低溫半导体电子束量子波算法计算

🍇Python椭圆微分-拉普拉斯方程

物理学中的许多问题与时间无关，但却具有丰富的物理意义：大质量物体产生的引力场、电荷分布的电势、拉伸膜的位移以及流体通过多孔介质的稳定流动……所有这些都可以用泊松方程建模：
$${\nabla }^{2}u=f$$
其中未知的 $u$ 和已知的 $f$ 是域 $\Omega$ 中的空间函数。为了找到解，我们需要边界条件。边值问题包括在给定上述信息的情况下找到 $u$。在数字上，我们可以使用松弛方法来做到这一点，该方法从对 $u$ 的初始猜测开始，然后迭代求解。

$f=0$（齐次情况）的特殊情况得出拉普拉斯方程：
$${\nabla }^{2}u=0$$
例如，稳定的二维热传导方程为：
$$\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}=0$$
其中 $T$​ 是已达到稳定状态的温度。拉普拉斯方程对系统在所提供的边界条件下的平衡状态进行建模。研究拉普拉斯方程解的学科称为势理论，解本身通常就是势场。从现在开始，我们用 $p$ 来表示我们的通用因变量，并再次写出拉普拉斯方程（二维）：
$$\frac{{\partial }^{2}p}{\partial x^2}+\frac{{\partial }^{2}p}{\partial y^2}=0$$
与扩散方程一样，我们用中心差离散化二阶导数
$$\frac{p_{i+1,j}-2p_{i,j}+p_{i-1,j}}{\Delta x^2}+\frac{p_{i,j+1}-2p_{i,j}+p_{i,j-1}}{\Delta y^2}=0$$
当 $\Delta x=\Delta y$ 时，我们最终得到以下等式：
$$p_{i+1,j}+p_{i-1,j}+p_{i,j+1}+p_{i,j-1}-4p_{i,j}=0$$
这告诉我们，网格点 $(i,j)$ 处的拉普拉斯微分算子可以使用该点处的 $p$ 值（因子 -4 ）和左右四个相邻点来离散计算，网格点 $(i,j)$ 上方和下方。

假设我们想在一块计算机芯片上模拟稳态热传递，该芯片一侧绝缘（零诺伊曼边界层），两侧保持固定温度（狄利克雷条件），一侧接触具有正弦温度分布的组件。我们需要求解拉普拉斯方程，其边界条件如下：
$$\begin{array}{c}p=0\text{ at }x=0\\ \frac{\partial p}{\partial x}=0\text{ at }x=L\\ p=0\text{ at }y=0\\ p=\sin \left(\frac{\frac{3}{2}\pi x}{L}\right)\text{ at }y=H.\end{array}$$
我们将 $L=1$ 和 $H=1$ 作为域在 $x$ 和 $y$ 方向上的大小。

from matplotlib import pyplot
import numpy
%matplotlib inline
from matplotlib import rcParams
rcParams['font.family'] = 'serif'
rcParams['font.size'] = 16from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cmdef plot_3D(x, y, p):fig = pyplot.figure(figsize=(11,7), dpi=100)ax = fig.gca(projection='3d')X,Y = numpy.meshgrid(x,y)surf = ax.plot_surface(X,Y,p[:], rstride=1, cstride=1, cmap=cm.viridis,linewidth=0, antialiased=False)ax.set_xlim(0,1)ax.set_ylim(0,1)ax.set_xlabel('$x$')ax.set_ylabel('$y$')ax.set_zlabel('$z$')ax.view_init(30,45)

具有上述边界条件的拉普拉斯方程有一个解析解，由下式给出
$$p(x,y)=\frac{\sinh \left(\frac{\frac{3}{2}\pi y}{L}\right)}{\sinh \left(\frac{\frac{3}{2}\pi H}{L}\right)}\sin \left(\frac{\frac{3}{2}\pi x}{L}\right)$$

def p_analytical(x, y):X, Y = numpy.meshgrid(x,y)p_an = numpy.sinh(1.5*numpy.pi*Y / x[-1]) /\(numpy.sinh(1.5*numpy.pi*y[-1]/x[-1]))*numpy.sin(1.5*numpy.pi*X/x[-1])return p_an

让我们尝试一下解析解，并用它来测试我们上面编写的plot_3D函数。

nx = 41
ny = 41x = numpy.linspace(0,1,nx)
y = numpy.linspace(0,1,ny)p_an = p_analytical(x,y)plot_3D(x,y,p_an)

👉参阅：计算思维 | 亚图跨际