欧拉近似

Question

我正在尝试为微分方程的欧拉逼近编写一个程序（我在一本书中遇到了一个问题告诉我）。我设法写了一些适用于所有东西的东西（我至少测试过) 但由于某种原因它在 step = 0,1 时失败了...... 这是代码：

from math import pow as p

y0: int = 3

h = input("Step size:")
h =float(h)
h1 = h
x_n=0
while x_n < (1):
    yn = y0 + h*(3*(p(x_n ,2))*(2-y0))
    x_n=h1
    h1=h1+h
    y0=yn
else:
    print(yn)

# h(1),f(1)=3
# h(0.1),f(1)=2,39279
# h(0.01),f(1)=2.37011
# h(0.001),f(1)=2.36810

它适用于 h=1、h=0.5、h=0.2、h=0.05、h=0.01、h=0.001 等等...出于某种原因，它仅在 h=0.1 时失败。

所讨论的微分方程是y'=3(x^2)(2-y)

感谢您的帮助！

Answer 1

你的代码没有问题。（您也可以使用 ** 而不是 pow；例如 pow(2,3) = 2**3 = 8

Answer 2

（如果变量名有一个统一的结构会更好看）

如果你也打印出最后一步的时间

for h in [1.0,0.5,0.2,0.1,0.05]:
    x_n, y_0 = 0.0, 3.0
    while x_n < (1):
        y_n = y_0 + h*(3*(x_n**2)*(2-y_0))
        x_n = x_n + h
        y_0 = y_n
    else:
        print(f"h={h}: f({x_n:12.10f}) = {y_n:12.10f}")

然后你在结果中找到

h= 1.00000: f(1.0000000000) = 3.0000000000
h= 0.50000: f(1.0000000000) = 2.6250000000
h= 0.20000: f(1.0000000000) = 2.4261034230
h= 0.10000: f(1.1000000000) = 2.2749557398
h= 0.05000: f(1.0000000000) = 2.3795619396

确实在 h=0.1 处，迭代仅在 x=1.1 处停止。这是由于变量的浮点格式，接近 x=1 的点 x_n 可以落在下面不久，因此执行另一个步骤。防止这种情况的一种方法是调整最后一步（另一种方法是使用“接近”比较）

for h in [1.0,0.5,0.2,0.1,0.05]:
    x_n, y_0, x_f = 0.0, 3.0, 1.0
    y_n = y_0
    while x_n < x_f:
        if x_n + 1.001*h > x_f: h = x_f - x_n
        y_n = y_n + h*(3*(x_n**2)*(2-y_n))
        x_n = x_n + h
    else:
        print(f"h={h:8.5f}: f({x_n:12.10f}) = {y_n:12.10f}")

给出更一致的结果

h= 1.00000: f(1.0000000000) = 3.0000000000
h= 0.50000: f(1.0000000000) = 2.6250000000
h= 0.20000: f(1.0000000000) = 2.4261034230
h= 0.10000: f(1.0000000000) = 2.3927939140
h= 0.05000: f(1.0000000000) = 2.3795619396

为了探究数字的浮点表示的影响，打印出 long 浮点数的真实含义，

print(", ".join("%20.17f"%x for x in np.cumsum(10*[0.1])))

给予

0.10000000000000001,  0.20000000000000001,  0.30000000000000004,  
0.40000000000000002,  0.50000000000000000,  0.59999999999999998,  
0.69999999999999996,  0.79999999999999993,  0.89999999999999991,  
0.99999999999999989

这表明浮点运算后最后一个尾数位的舍入可能会产生相当违反直觉的结果

欧拉近似

Euler's Approximation

python

math

differential-equations

calculus