python求解一个偏微分方程例子

以下是一个求解二维泊松方程的 Python 代码示例：

import numpy as np
import matplotlib.pyplot as plt

# 定义问题域
Lx = 1.0
Ly = 1.0
nx = 101
ny = 101
dx = Lx / (nx - 1)
dy = Ly / (ny - 1)

# 定义边界条件
u_top = 0.0
u_bottom = 0.0
u_left = 0.0
u_right = 1.0

# 定义初始解
u = np.zeros((ny, nx))

# 设置边界条件
u[0, :] = u_top
u[-1, :] = u_bottom
u[:, 0] = u_left
u[:, -1] = u_right

# 定义迭代函数
def iterate(u, dx, dy):
    un = u.copy()
    r = 0.25 * (dx ** 2 + dy ** 2)
    for j in range(1, ny - 1):
        for i in range(1, nx - 1):
            u[j, i] = (un[j, i + 1] + un[j, i - 1] + un[j + 1, i] + un[j - 1, i]) / 4.0 + r * f[i, j]
    return u

# 定义右端项
f = np.zeros((ny, nx))
f[:, :] = -4.0 * np.pi ** 2 * np.sin(np.pi * np.arange(ny)[:, np.newaxis] / (ny - 1)) ** 2 * np.sin(np.pi * np.arange(nx) / (nx - 1)) ** 2

# 迭代求解
n_iterations = 10000
for n in range(n_iterations):
    u = iterate(u, dx, dy)

# 绘制解
X, Y = np.meshgrid(np.linspace(0, Lx, nx), np.linspace(0, Ly, ny))
plt.contourf(X, Y, u, cmap='coolwarm')
plt.colorbar()
plt.show()

在此示例中，我们使用迭代法求解二维泊松方程，其中边界条件为 $u(x,0)=u(x,1)=0$ 和 $u(0,y)=0, u(1,y)=1$。我们使用了 $101 \times 101$ 的网格，并迭代求解 $10000$ 次。最终将求得的解绘制出来。