无法通过 diffeqpy 求解微分方程

Question

最近了解到diffeqpy，想尝试一下。我有这个 Julia 代码，它有效并给出了预期的结果：

using DifferentialEquations
using Plots

function bcm4(resid, dx, x, params, t)
    # Algebraic equations
    resid[11] = params[3]*params[4]/16*x[9] - x[11]
    resid[12] = params[3]*params[4]/44*x[10] - x[12]
    resid[13] = x[11] + x[12] + params[11] - x[13]
    resid[14] = params[9]*(x[13] - params[15])*x[13]/params[15] - x[14]
    
    # Rate equations
    rate = [params[6]*x[5],
      params[8]*x[6],
      params[9]*x[7],
      params[7]*x[8],
      params[5]*(x[1] - 16*params[1]*x[11]),
      params[5]*(x[2] - 44*params[2]*x[12])]

    # Process equations
    process = [0.24819*rate[1] + 0.32208*rate[2] + 0.63928*rate[3] - rate[5],
      0.68087*rate[1] + 0.79543*rate[2] + 0.58172*rate[3]   + 0.13924*rate[4] - rate[6],
      -0.02065*rate[1] + 0.16892*rate[2] - 0.034418*rate[3],
      -0.045576*rate[1] - 0.45876*rate[2] - 0.41518*rate[3] - 0.022387*rate[4],
      -rate[1] + 0.05*rate[4],
      -rate[2] + 0.77257*rate[4],
      -rate[3] + 0.060578*rate[4],
      0.1125*rate[1] + 0.1723*rate[2] + 0.2286*rate[3] - rate[4],
      params[12]/params[13]*rate[5],
      params[12]/params[13]*rate[6]]
    
    # InOut
    in_out = [params[16]/params[12]*(params[17] - x[1]),
      params[16]/params[12]*(params[18] - x[2]),
      params[16]/params[12]*(params[19] - x[3]),
      params[16]/params[12]*(params[20] - x[4]),
      params[16]/params[12]*(params[21] - x[5]),
      params[16]/params[12]*(params[22] - x[6]),
      params[16]/params[12]*(params[23] - x[7]),
      params[16]/params[12]*(params[24] - x[8]),
      -x[9]*x[14]/params[13],
      -x[10]*x[14]/params[13]]
    
   # Differential Equations
   resid[1:10] = in_out + process - dx[1:10]
end

x0 = [0.01125, 0.5743, 1.6758, 761.37, 4.3936, 0.6606, 0.0830, 9.2012, 
    0.3455, 0.9169, 0.1, 0.1, 0.1, 1]
tspan = (0.0, 5.0)
params = [0.00116190, 0.0271467, 0.0831450, 273.150, 200, 0.50, 0.020, 0.50,
    50000, 0.500, 0.0556677,
    155, 20, 0, 1.03125,
    2.5, 0, 0, 0.5, 788.609, 134.22, 12.287, 1.9615, 0]
differentialVars = [true, true, true, true, true, true, true, true, true, true,
    false, false, false, false]
prob = DAEProblem(bcm4, x0/10, x0, tspan, params, differential_vars = differentialVars)
sol = solve(prob, saveat=1/12)
plot(sol, vars=[5])

现在，如果我尝试通过 diffeqpy 在 Python 中解决它，它不起作用。

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
BCM4, solved in Julia via diffeqpy, coded in Python
@author: winkmal
"""
import numpy
import matplotlib.pyplot as plt
from diffeqpy import de

bcm4_Julia = de.eval("""
function bcm4(resid, dx, x, params, t)
#Past code from above
end""")

x0 = numpy.array([0.01125, 0.5743, 1.6758, 761.37, 4.3936, 0.6606, 0.0830, 9.2012, \
                  0.3455, 0.9169, 0.1, 0.1, 0.1, 1])
tspan = (0.0, 5.0)
params = numpy.array([0.00116190, 0.0271467, 0.0831450, 273.150, 200, 0.50, 0.020, 0.50, \
    50000, 0.500, 0.0556677, \
    155, 20, 0, 1.03125, \
    2.5, 0, 0, 0.5, 788.609, 134.22, 12.287, 1.9615, 0])
differentialVars = [True, True, True, True, True, True, True, True, True, True, \
                    False, False, False, False]
prob = de.DAEProblem(bcm4_Julia, x0/10, x0, tspan, params, differential_vars = differentialVars)
sol = de.solve(prob, saveat=1/12)

我收到以下错误消息：

    prob = de.DAEProblem(bcm4_Julia, x0/10, x0, tspan, params, differential_vars = differentialVars)
ValueError: Julia exception: ArgumentError: reducing over an empty collection is not allowed

你能告诉我这个错误是什么意思，以及如何解决吗？

Answer 1

您需要从 Main 开始评估。工作代码：

import numpy
import matplotlib.pyplot as plt
from diffeqpy import de
from julia import Main

bcm4 = Main.eval("""
function bcm4(resid, dx, x, params, t)
    # Algebraic equations
    resid[11] = params[3]*params[4]/16*x[9] - x[11]
    resid[12] = params[3]*params[4]/44*x[10] - x[12]
    resid[13] = x[11] + x[12] + params[11] - x[13]
    resid[14] = params[9]*(x[13] - params[15])*x[13]/params[15] - x[14]

    # Rate equations
    rate = [params[6]*x[5],
      params[8]*x[6],
      params[9]*x[7],
      params[7]*x[8],
      params[5]*(x[1] - 16*params[1]*x[11]),
      params[5]*(x[2] - 44*params[2]*x[12])]

    # Process equations
    process = [0.24819*rate[1] + 0.32208*rate[2] + 0.63928*rate[3] - rate[5],
      0.68087*rate[1] + 0.79543*rate[2] + 0.58172*rate[3]   + 0.13924*rate[4] - rate[6],
      -0.02065*rate[1] + 0.16892*rate[2] - 0.034418*rate[3],
      -0.045576*rate[1] - 0.45876*rate[2] - 0.41518*rate[3] - 0.022387*rate[4],
      -rate[1] + 0.05*rate[4],
      -rate[2] + 0.77257*rate[4],
      -rate[3] + 0.060578*rate[4],
      0.1125*rate[1] + 0.1723*rate[2] + 0.2286*rate[3] - rate[4],
      params[12]/params[13]*rate[5],
      params[12]/params[13]*rate[6]]

    # InOut
    in_out = [params[16]/params[12]*(params[17] - x[1]),
      params[16]/params[12]*(params[18] - x[2]),
      params[16]/params[12]*(params[19] - x[3]),
      params[16]/params[12]*(params[20] - x[4]),
      params[16]/params[12]*(params[21] - x[5]),
      params[16]/params[12]*(params[22] - x[6]),
      params[16]/params[12]*(params[23] - x[7]),
      params[16]/params[12]*(params[24] - x[8]),
      -x[9]*x[14]/params[13],
      -x[10]*x[14]/params[13]]

   # Differential Equations
   resid[1:10] = in_out + process - dx[1:10]
end""")

x0 = numpy.array([0.01125, 0.5743, 1.6758, 761.37, 4.3936, 0.6606, 0.0830, 9.2012, \
                  0.3455, 0.9169, 0.1, 0.1, 0.1, 1])
tspan = (0.0, 5.0)
params = numpy.array([0.00116190, 0.0271467, 0.0831450, 273.150, 200, 0.50, 0.020, 0.50, \
    50000, 0.500, 0.0556677, \
    155, 20, 0, 1.03125, \
    2.5, 0, 0, 0.5, 788.609, 134.22, 12.287, 1.9615, 0])
differentialVars = [True, True, True, True, True, True, True, True, True, True, \
                    False, False, False, False]
prob = de.DAEProblem(bcm4, x0/10, x0, tspan, params, differential_vars = differentialVars)
sol = de.solve(prob, saveat=1/12)

无法通过 diffeqpy 求解微分方程

Trouble getting Differential Equation to solve via diffeqpy

python-3.x

differential-equations

julia

differentialequations.jl