Реализация метода Эйлера в Python дает стабильный результат, но он должен быть нестабильным - PullRequest
0 голосов
/ 28 сентября 2018

Я пытаюсь решить это дифференциальное уравнение методом Эйлера, используя Python3:

enter image description here

Согласно Wolfram Alpha, это сюжет правильногоуравнение.

enter image description here

Опять-таки, согласно Вольфраму Альфе, в этом случае классический метод Эйлера должен не быть стабильным, как вы можетесм. к концу интервала:

enter image description here

Однако в моей реализации метод Эйлера обеспечивает стабильный результат, что странно .Интересно, что моя реализация по какой-то причине неверна.Тем не менее, я не могу найти ошибку.

Я сгенерировал несколько точек и график, сравнивающий мое приближение и аналитический вывод функции.В синем цвете аналитический результат в качестве контрольной группы.Красным цветом вывод моей реализации:

enter image description here

Вот мой код:

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):

enter image description here

1 Ответ

0 голосов
/ 28 сентября 2018

В зависимости от размера временного шага интеграция Эйлера может быть стабильной или нестабильной, поскольку это явный метод.Вы выбрали довольно маленький временной шаг.Если вы увеличите его, вы начнете видеть колебания, как показано на рисунке ниже.

Oscillations

Вот небольшая программа тестирования, которую я написал (попробуйте медленноувеличение переменной steps [20,30,40,50 ....]):

import numpy as np
import matplotlib.pyplot as plt

steps = 20


def exact_solution(t):
    return np.exp(-10.0 * t)


def numerical_solution(y0, dt, num_steps):
    y = np.zeros(num_steps + 1)
    y[0] = y0
    for step in range(num_steps):
        y[step + 1] = y[step] - 10.0 * y[step] * dt

    return y


if __name__ == "__main__":
    t0 = 0
    time = np.linspace(t0, 5, steps + 1)
    num_sol = numerical_solution(exact_solution(t0), time[1] - time[0], steps)
    exact_sol = exact_solution(time)

    plt.plot(time, num_sol, ".--", label="numerical")
    plt.plot(time, exact_sol, label="exact")
    plt.legend(loc="best")
    plt.show()
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...