伽辽金法求解微分方程代码

伽辽金法（Galerkin Method）是一种求解微分方程的数值方法，下面是伽辽金法的求解微分方程的代码：

import numpy as np

def Galerkin_Method(a, b, n, m):
    h = (b - a) / n
    x = np.linspace(a, b, n+1)
    t = np.linspace(0, m, m+1)
    alpha = np.zeros(n+1)
    beta = np.zeros(n+1)
    u = np.zeros((n+1, m+1))
    
    # 初始化初始条件
    for i in range(n+1):
        alpha[i] = f(x[i])
        u[i, 0] = alpha[i]
    
    # 计算常数系数
    A = np.zeros((n+1, n+1))
    for i in range(n+1):
        for j in range(n+1):
            A[i, j] = h * np.sum(phi(i, x) * phi(j, x))
    
    # 计算右端项
    F = np.zeros(n+1)
    for i in range(n+1):
        F[i] = h * np.sum(f(x) * phi(i, x))
    
    # 求解矩阵方程
    beta = np.linalg.solve(A, F)
    
    # 计算数值解
    for j in range(1, m+1):
        for i in range(n+1):
            u[i, j] = alpha[i] + np.sum(beta * phi(i, x) * np.exp(-t[j]))
    
    return u

# 定义初始条件和方程
def f(x):
    return np.sin(x)

def phi(i, x):
    return np.sin((2 * i + 1) * x / 2)

def u_exact(x, t):
    return np.exp(-t) * np.sin(x)

# 调用函数并绘制图像
u = Galerkin_Method(0, np.pi, 100, 100)
x = np.linspace(0, np.pi, 100+1)
t = np.linspace(0, 1, 100+1)
for j in range(0, 101, 10):
    plt.plot(x, u[:, j], label='t={}'.format(round(t[j], 2)))
plt.plot(x, u_exact(x, 1), 'k--', label='Exact')
plt.legend()
plt.show()

其中，`a` 和 `b` 分别是区间的左右端点，`n` 是区间的划分数，`m` 是时间的划分数。`alpha` 和 `beta` 分别是常数系数，`u` 是数值解。`f(x)` 是初始条件，`phi(i, x)` 是基函数，`u_exact(x, t)` 是精确解。最后绘制出数值解的图像，以及精确解的图像。