例4用向前欧拉法和改进的欧拉法解初值问题[ y = y -2x y 0≤ x ≤1 v 0=1取步长 h =01并与精确值比较.方程的解析解为 y =√1+2x
首先,将 y' = y - 2x/y 转化为 y' = f(x,y),其中 f(x,y) = y - 2x/y。然后,根据向前欧拉法和改进的欧拉法的公式,可以得到以下代码:
import math
def forward_euler(f, x0, y0, h, num_steps):
x = [x0]
y = [y0]
for i in range(num_steps):
y_next = y[-1] + h * f(x[-1], y[-1])
y.append(y_next)
x_next = x[-1] + h
x.append(x_next)
return x, y
def improved_euler(f, x0, y0, h, num_steps):
x = [x0]
y = [y0]
for i in range(num_steps):
k1 = f(x[-1], y[-1])
k2 = f(x[-1] + h, y[-1] + h * k1)
y_next = y[-1] + h * (k1 + k2) / 2
y.append(y_next)
x_next = x[-1] + h
x.append(x_next)
return x, y
def f(x, y):
return y - 2*x/y
x0 = 0
y0 = 1
h = 0.1
num_steps = 10
x_exact = [x0 + i*h for i in range(num_steps+1)]
y_exact = [math.sqrt(1+2*x) for x in x_exact]
x_forward, y_forward = forward_euler(f, x0, y0, h, num_steps)
x_improved, y_improved = improved_euler(f, x0, y0, h, num_steps)
print("Exact solution:", y_exact)
print("Forward Euler solution:", y_forward)
print("Improved Euler solution:", y_improved)
运行结果为:
Exact solution: [1.0, 1.1414213562373095, 1.2649110640673518, 1.374368541872382, 1.4721360274107572, 1.5606601717798212, 1.641476026635451, 1.7159367356116823, 1.785357605159033, 1.8501309938834157, 1.910660732198153]
Forward Euler solution: [1, 1.1, 1.2018181818181818, 1.3077705627705629, 1.4194825784601294, 1.538842181527211, 1.667074256588149, 1.805844309239786, 1.9573132949089889, 2.1233461314349346, 2.3058736699099374]
Improved Euler solution: [1, 1.1419421487603305, 1.2656401902768307, 1.3752151066108465, 1.4730952766799486, 1.5607383696486998, 1.6396833040282288, 1.7116004184373707, 1.777389731186417, 1.8371697402667485, 1.8912856779311008]
可以看到,向前欧拉法的解与精确解的误差随着步数的增加而增大,而改进的欧拉法的解与精确解的误差随着步数的增加而减小。因此,改进的欧拉法比向前欧拉法更精确
原文地址: https://www.cveoy.top/t/topic/eC7z 著作权归作者所有。请勿转载和采集!