2

Python3 を使用してオイラー法でこの微分方程式を解こうとしています。

ここに画像の説明を入力

Wolfram Alpha によると、これが正しい方程式のプロットです。

ここに画像の説明を入力

繰り返しになりますが、Wolfram Alpha によると、この場合、間隔の終わりまでにわかるように、従来のオイラー法は安定していないはずです。

ここに画像の説明を入力

ただし、私の実装では、オイラー法は安定した結果を提供しますが、これは奇妙なことです。何らかの理由で私の実装が間違っているのだろうか。それにもかかわらず、私はエラーを見つけることができません。

近似値と関数の分析出力を比較するいくつかの点とプロットを生成しました。青は対照群としての分析結果。赤で、私の実装の出力:

ここに画像の説明を入力

それが私のコードです:

import math
import numpy as np
from matplotlib import pyplot as plt
import pylab

def f(x):

    return (math.e)**(-10*x)

def euler(x):

    y_init = 1
    x_init = 0

    old_dy_dx = -10*y_init

    old_y = y_init 

    new_y = None

    new_dy_dx = None

    delta_x = 0.001

    limite = 0

    while x>limite:

        #for i in range(1,6):

        new_y = delta_x*old_dy_dx + old_y
        #print ("new_y", new_y)

        new_dy_dx = -10*new_y
        #print ("new dy_dx", new_dy_dx)

        old_y = new_y
        #print ("old_y", old_y)

        old_dy_dx = new_dy_dx
        #print ("old delta y_delta x", old_dy_dx)
        #print ("iterada",i)

        limite = limite +delta_x

    return new_y

t = np.linspace(-1,5, 80)

lista_outputs = []

for i in t:
    lista_outputs.append(euler(i))
    print (i)

# red dashes, blue squares and green triangles
plt.plot(t, f(t), 'b-', label='Output resultado analítico')
plt.plot(t , lista_outputs, 'ro', label="Output resultado numérico")
plt.title('Comparação Euler/Analítico - tolerância: 0.3')
pylab.legend(loc='upper left')
plt.show()

助けてくれてありがとう。

================================================== ==========

アップデート

@SourabhBhat の助けを借りて、実装が実際に正しいことを確認できました。それは確かに、不安定さを生み出していました。ステップサイズを大きくするだけでなく、ズームインしてそれが起こっていることを確認する必要がありました.

次の図は、それ自体を物語っています (ステップ サイズ 0.22)。

ここに画像の説明を入力

4

2 に答える 2