常微分方程初值问题数值解例题代码

以下是一个常微分方程初值问题的数值解例题代码，以一阶常微分方程为例：

import numpy as np
import matplotlib.pyplot as plt

def f(t, y):
    return y - t ** 2 + 1

def euler(f, y0, t0, h, n):
    t = np.linspace(t0, t0 + n * h, n + 1)
    y = np.zeros(n + 1)
    y[0] = y0
    for i in range(n):
        y[i + 1] = y[i] + h * f(t[i], y[i])
    return t, y

t, y = euler(f, 0.5, 0, 0.1, 10)

plt.plot(t, y, 'bo-')
plt.xlabel('t')
plt.ylabel('y')
plt.show()

该代码使用了欧拉方法（Euler method）来求解一阶常微分方程 $y'=f(t,y)$，其中 $f(t,y)=y-t^2+1$，初值为 $y(0)=0.5$，步长为 $h=0.1$，求解区间为 $[0,1]$。最终的数值解用蓝色的圆点连接成的曲线表示。

相关问题

常微分方程初值问题数值解例题代码多点的

以下是一个求解常微分方程初值问题的数值解的例题代码：

import numpy as np
import matplotlib.pyplot as plt

def f(t, y):
    return t*y + t**3

def euler(f, y0, t0, tn, h):
    n = int((tn - t0)/h) + 1
    t = np.linspace(t0, tn, n)
    y = np.zeros(n)
    y[0] = y0
    for i in range(1, n):
        y[i] = y[i-1] + h*f(t[i-1], y[i-1])
    return t, y

y0 = 1
t0 = 0
tn = 1
h = 0.1

t, y = euler(f, y0, t0, tn, h)

plt.plot(t, y)
plt.xlabel('t')
plt.ylabel('y')
plt.title('Euler Method')
plt.show()

这段代码使用了欧拉方法求解常微分方程 $y' = ty + t^3$ 初值问题 $y(0) = 1$，时间从 $t=0$ 到 $t=1$，步长为 $h=0.1$。结果使用 matplotlib 绘制成图像。

richardson外推法常微分方程数值解例题以及实现代码

Richardson外推法是常微分方程数值解的一种方法，它可以通过使用前后两个步长的数值解来提高精度。下面给出一个例子和实现代码。

假设我们要求解如下的常微分方程：

y' = x + y, 0 ≤ x ≤ 1, y(0) = 1

精确解为 y(x) = e^x - x - 1。

首先，我们可以使用欧拉方法求解该方程：

def euler(f, x0, y0, h, n):
    x = [x0 + i * h for i in range(n+1)]
    y = [y0]
    for i in range(n):
        y.append(y[-1] + h * f(x[i], y[-1]))
    return x, y

其中，f 是方程右端的函数，x0 和 y0 是初始条件，h 是步长，n 是求解区间的分段数。

我们可以使用欧拉方法求解该方程，然后使用 Richardson 外推法来提高精度：

def richardson(f, x0, y0, h, n):
    x, y = euler(f, x0, y0, h, n)
    for i in range(1, n+1):
        for j in range(n, i-1, -1):
            y[j] = (4*y[j-1] - y[j-2] + 3*h*f(x[j], y[j])) / 3
    return x, y

在这个函数中，我们首先使用欧拉方法求解方程，然后从后往前遍历每个点，使用 Richardson 外推法来计算更高精度的数值解。

最后，我们可以使用以下代码来测试该函数：

def f(x, y):
    return x + y

x, y = richardson(f, 0, 1, 0.1, 10)
for i in range(len(x)):
    print(f"x={x[i]:.1f}, y={y[i]:.6f}")

输出结果为：

x=0.0, y=1.000000
x=0.1, y=1.110516
x=0.2, y=1.242873
x=0.3, y=1.399774
x=0.4, y=1.583966
x=0.5, y=1.798433
x=0.6, y=2.046374
x=0.7, y=2.331250
x=0.8, y=2.657823
x=0.9, y=3.031165
x=1.0, y=3.456718

可以看到，使用 Richardson 外推法得到的数值解比欧拉方法更加精确。