首先,将 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 著作权归作者所有。请勿转载和采集!

免费AI点我,无需注册和登录