Как изменить математический код ODE на Python - PullRequest
Новая
       26

Как изменить математический код ODE на Python

0 голосов
/ 04 июля 2019

Мне дан следующий код в Mathematica

solution = NDSolve[{\[CapitalOmega]\[CapitalLambda]'[x] == 
(1 - \[CapitalOmega]\[CapitalLambda][x]) \[CapitalOmega]\[CapitalLambda][x] ((1 - \[CapitalOmega]\[CapitalLambda][x])^((1 - \[Delta])/(4 -\2 \[Delta]))  \[CapitalOmega]\[CapitalLambda][x]^(1/(4 - 2 \[Delta]))), \[CapitalOmega]\[CapitalLambda][0] == \[CapitalOmega]\[CapitalLambda]f}, \[CapitalOmega]\\[CapitalLambda], {x, xi, xf}];

Я хочу написать это на Python, используя ODEINT но я действительно не понимаю способ написания этого, потому что он посылает `error1 комплексного числа. 

def OD_H(z, od, delt):
    dMdt = od * (1 - od) * ((1 - od)**((1 - delt)/(4-2*delt)) * od**(1/(2 *(2-delt))))
    return dMdt                                                                                                                      

def ant(z, od0, delt):
    z1 = 0
    od = odeint(OD_H, od0, [z1, z], args=(delt,))[-1]       
    return od     

for z in np.arange(0,3.1,0.1):
    print(ant(z, 0.7, 1.1))

ошибка

od = odeint(OD_H, od0, [z1, z], args=(delt,))[-1] File "C:\Python36-32\lib\site-packages\scipy\integrate\odepack.py", line 244, in odeint int(bool(tfirst))) TypeError: can't convert complex to float

1 Ответ

1 голос
/ 05 июля 2019

Из вызова odeint мы узнаем, что вы используете z в качестве переменной времени и od в качестве переменной пространства / состояния. Поскольку ваша производная функция имеет первое время, а по умолчанию в odeint указано первое состояние, вам нужно установить параметр tfirst=True. С этим изменением нет ошибок. Вы получаете те же значения с меньшими вычислениями, используя полную мощность odeint, предоставляя ей список всех желаемых точек выборки 

z = np.arange(0,3.1,0.1)
od = odeint(OD_H, od0, z, args=(delt,), tfirst=True)

где пары z[k], od[k] соответствуют парам аргумент-значение, которые вы вычисляете в цикле.

...