Есть ли способ изменить мой ODE внутри ODEINT в Python? - PullRequest
Новая
       98

Есть ли способ изменить мой ODE внутри ODEINT в Python?

0 голосов
/ 15 апреля 2020

Я пытаюсь смоделировать колебательную систему с небольшой разницей: я хочу, чтобы она использовала конкретное c уравнение движения (ОДУ), когда тело движется вверх, и другое уравнение движения, когда тело двигаясь вниз. Чтобы решить эти уравнения, я использую ODEINT от Scypi.

Например, давайте рассмотрим систему масс-пружин classi c. Я пытался заставить систему работать с уравнением движения для внешних возбуждений на теле, когда оно движется вверх, и с простым уравнением, когда оно движется вниз. 

def function (x,t):

    F0 = 10.00
    w = 1.00
    m = 2.00
    c = 1.00
    k = 20.00
    s = x[0]
    dsdt = x[1]

    if x[1] >= 0:
        d2sdt2 = (F0*np.sin(w*t)-c*dsdt-k*s)/m
    else:
        d2sdt2 = (-c*dsdt-k*s)/m

    result = [dsdt,d2sdt2]
    return result

initial = [3.00,0.00]
t = np.linspace(0.00,10.00,101)
y = odeint(function, initial, t)

Полученные результаты показывают, что на тело воздействует только второе уравнение движения ( Полученные результаты ). Я ожидал более хаотического c паттерна движения, когда тело движется вверх, из-за внешней силы.

Есть ли лучший способ реализовать это?

1 Ответ

0 голосов
/ 15 апреля 2020

Простое добавление некоторых параметров для увеличения плотности внутренних шагов и плотности вывода 

t = np.linspace(0.00,10.00,301)
y = odeint(function, initial, t, hmax=0.1, atol = 1e-8, rtol=1e-10)

без использования пунктирных линий дает график enter image description here

где изломы во второй половине отчетливо видны.

Таким образом, ваш код верен, но вы должны принять во внимание, что решатели ODE построены для гладких правых частей, начиная с математических методов и заканчивая и особенно для предикторов для оптимального размера шага. Так как здесь есть локально непредсказуемые негладкие изменения в правой части, нужно сказать, чтобы решатель не использовал очень большие размеры шагов, так как это может перепрыгнуть через изменение фазы модели и обратно.

Как Кроме того, для гладких графиков выходная плотность 100 точек по всей горизонтальной оси обычно слишком мала. Используйте от 200 до 600 точек в зависимости от окончательного размера печати в печатном виде.

...