显式 Euler 方法不符合我的预期

Question

我在python中实现了以下显式欧拉方法：

def explicit_euler(df, x0, h, N):
    """Solves an ODE IVP using the Explicit Euler method.

    Keyword arguments:
    df  - The derivative of the system you wish to solve.
    x0 - The initial value of the system you wish to solve.
    h  - The step size.
    N  - The number off steps.
    """
    x = np.zeros(N)
    x[0] = x0

    for i in range(0, N-1):
        x[i+1] = x[i] + h * df(x[i])

    return x

关注 wikipedia I can plot the function and verify that I get the same plot: 上的文章。我相信这里我写的方法是正确的。

接下来我尝试用它来解决这个 page 上给出的最后一个系统，而不是那里显示的情节我得到了这个：

我不确定为什么我的情节与网页上显示的情节不符。当我用它来解决斜率不变的系统时，显式欧拉方法似乎工作正常，但对于振荡函数，它似乎根本无法模仿它。甚至没有显示链接网页上指示的预期误差增益。我不确定我实现的方法有什么问题。

这里是用于绘图和导数的代码：

def g(t):
    return -0.5 * np.exp(t * 0.5) * np.sin(5 * t) + 5 * np.exp(t * 0.5) 
    * np.cos(5 * t)

h = 0.001
x0 = 0
tn = 4
N = int(tn / h)

x = ee.explicit_euler(f, x0, h, N)
t = np.arange(0, tn, h)

fig = plt.figure()
plt.plot(t, x, label="Explicit Euler")
plt.plot(t, (np.exp(0.5 * t) * np.sin(5 * t)), label="Analytical 
solution")
#plt.plot(t, np.exp(0.5 * t), label="Analytical solution")
plt.xlabel('Timesteps t')
plt.ylabel('x(t)=e^(0.5*t) * sin(5*t)')
plt.legend()
plt.grid()
plt.show()

编辑：

此处要求的是我将此方法应用到的当前方程式：

y'-y=-0.5*e^(t/2)*sin(5t)+5e^(t/2)*cos(5t)

其中 y(0)=0。

不过，我想澄清一下，这种行为不仅仅发生在这个方程中，而是出现在所有方程中，其中斜率的符号或振荡行为发生变化。

编辑 2：好，谢谢。是的，下面的代码确实有效。但我还有一个问题。在指数函数的简单示例中，我定义了一个方法：

def f(x): 
    return x

对于系统 f'(x)=x。这给出了我的第一张图的输出，它看起来是正确的。然后我定义了另一个函数：

def k(x): 
    return cos(x)

对于系统 f'(x)=cos(x)，这不会给出预期的输出。但是当我将函数定义更改为

def k(t, x): 
    return cos(t)

我得到了预期的输出。如果我改变我的功能

def f(t, x): 
    return t

我的输出不正确。我是否总是在一个时间步实际评估函数，系统 x'=x 是否只是偶然在每个时间步的值只是 x 的值？

我已经明白欧拉法是为了得到下一个值，而使用上一个计算值的值。但是如果我运行为我的函数 k(x)=cos(x) 编写代码，我会得到如下图所示的输出，这一定是不正确的。这现在使用您提供的更新代码。

def k(t, x):
    return np.cos(x)

h  = 0.1         # Step size
x0 = (0, 0)        # Initial point of iteration
tn = 10        # Time step to iterate to
N  = int(tn / h)   # Number of steps

x = ee.explicit_euler(k, x0, h, N)
t = np.arange(0, tn, h)

Answer 1

问题是你错误的提出了函数g，你要解方程：

从我们观察到的地方：

y' = y -0.5*e^(t/2)*sin(5t)+5e^(t/2)*cos(5t)

那么我们定义函数f(t, y) = y -0.5*e^(t/2)*sin(5t)+5e^(t/2)*cos(5t)为：

def f(t, y):
    return y -0.5 * np.exp(t * 0.5) * np.sin(5 * t) + 5 * np.exp(t * 0.5) * np.cos(5 * t)

迭代初始点为f0=(t(0), y(0)):

f0 = (0, 0)

然后根据欧拉方程：

def explicit_euler(df, x0, h, N):
    """Solves an ODE IVP using the Explicit Euler method.

    Keyword arguments:
    df  - The derivative of the system you wish to solve.
    x0 - The initial value of the system you wish to solve.
    h  - The step size.
    N  - The number off steps.
    """
    x = np.zeros(N)
    t, x[0] = x0

    for i in range(0, N-1):
        x[i+1] = x[i] + h * df(t ,x[i])
        t += h

    return x

完整代码：

def explicit_euler(df, x0, h, N):
    """Solves an ODE IVP using the Explicit Euler method.

    Keyword arguments:
    df  - The derivative of the system you wish to solve.
    x0 - The initial value of the system you wish to solve.
    h  - The step size.
    N  - The number off steps.
    """
    x = np.zeros(N)
    t, x[0] = x0

    for i in range(0, N-1):
        x[i+1] = x[i] + h * df(t ,x[i])
        t += h

    return x

def df(t, y):
    return -0.5 * np.exp(t * 0.5) * np.sin(5 * t) + 5 * np.exp(t * 0.5) * np.cos(5 * t) + y

h = 0.001
f0 = (0, 0)
tn = 4
N = int(tn / h)

x = explicit_euler(df, f0, h, N)
t = np.arange(0, tn, h)

fig = plt.figure()
plt.plot(t, x, label="Explicit Euler")
plt.plot(t, (np.exp(0.5 * t) * np.sin(5 * t)), label="Analytical solution")
#plt.plot(t, np.exp(0.5 * t), label="Analytical solution")
plt.xlabel('Timesteps t')
plt.ylabel('x(t)=e^(0.5*t) * sin(5*t)')
plt.legend()
plt.grid()
plt.show()

截图：

Dump y' 右边是你应该放在df函数里的东西

我们将修改变量以保持变量的相同标准，并将y作为因变量，t作为自变量。

等式 2：在这种情况下，等式 f'(x)=cos(x) 将被重写为：

y'=cos(t)

然后：

def df(t, y):
    return np.cos(t)

总而言之，如果我们有以下形式的等式：

y' = f(t, y)

然后：

def df(t, y):
    return f(t, y)

显式 Euler 方法不符合我的预期

Explicit Euler method doesn't behave how I expect

python

numerical

numerical-methods

numerical-integration