I am trying to use diffeqpy to use DAESolver from my python code.
The following code works:
# the residual calculation function
def f(du, u, p, t):
resid1 = -0.04 * u[0] + 1e4 * u[1] * u[2] - du[0]
resid2 = +0.04 * u[0] - 3e7 * u[1] ** 2 - 1e4 * u[1] * u[2] - du[1]
resid3 = u[0] + u[1] + u[2] - 1.0
res = [resid1, resid2, resid3]
return res
# initial values
u0 = [1.0, 0.0, 0.0]
du0 = [-0.04, 0.04, 0.0]
tspan = (0.0, 1.0)
differential_vars = [True, True, False]
# crate the DAE problem and solve
prob = de.DAEProblem(f, du0, u0, tspan, differential_vars=differential_vars)
sol = de.solve(prob, de.IDA(), abstol=1e-8, reltol=1e-8)
The following fails:
# the residual calculation function
def f(du, u, p, t):
resid1 = -0.04 * u[0] + 1e4 * u[1] * u[2] - du[0]
resid2 = +0.04 * u[0] - 3e7 * u[1] ** 2 - 1e4 * u[1] * u[2] - du[1]
resid3 = u[0] + u[1] + u[2] - 1.0
res = np.array([resid1, resid2, resid3]) // <---- changed
return res
# initial values
u0 = np.array([1.0, 0.0, 0.0]) // <----- changed
du0 = np.array([-0.04, 0.04, 0.0]) // <----- changed
tspan = (0.0, 1.0)
differential_vars = [True, True, False]
# crate the DAE problem and solve
prob = de.DAEProblem(f, du0, u0, tspan, differential_vars=differential_vars)
sol = de.solve(prob, de.IDA(), abstol=1e-8, reltol=1e-8)
The following is the error: (I can share the complete error output if it helps):
juliacall.JuliaError: MethodError: Cannot `convert` an object of type Vector{Float64} to an object of type PyArray{Float64, 1, true, true, Float64}
The part calling the solver is a part of a larger system wherein we have generated the matrices for the system and wish to solve it. The residual function here (f) therefore needs to operate on numpy arrays.
I am new to Julia. Is this something that is possible at all? Alternatively, is it possible to modify the above in any way to make it work, like converting those numpy matrices to julia ones and doing the residual calculation on the julia side?