I am given the following code within 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}];
I want to write it in Python, using ODEINT
but I really do not understand the way of writing it, because it sends `error1 of complex number.
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))
the error is
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
From the odeint
call we learn that you use z
as time variable and od
as space/state variable. As your derivatives function has the time first, and the default in odeint
is to have the state first, you need to set the option tfirst=True
. With that change there are no errors. You get the same values with less computation by using the full power of odeint
in giving it a list of all desired sampling points
z = np.arange(0,3.1,0.1)
od = odeint(OD_H, od0, z, args=(delt,), tfirst=True)
where then the pairs z[k], od[k]
correspond to the argument-value pairs that you compute in the loop.