Julia diffeqpy issue - Julia version works, Python does not

diffeqpy_test.py

from diffeqpy import de
from julia import Main

jul_code = '''function(du,u,p,t)
du[1] = p[20].*u[19] + (-1)*(p[19].*u[1].*u[2])
du[2] = p[20].*u[19] + (-1)*(p[19].*u[1].*u[2])
du[3] = p[23].*u[21] + (-1)*(p[22].*u[20].*u[3])
du[4] = p[25].*u[22] + p[36].*u[38] + (-1)*(p[24].*u[20].*u[4]) + (-1)*(p[35].*u[4].*u[34])
du[5] = p[28].*u[24] + (-1)*(p[27].*u[23].*u[5])
du[6] = p[81].*u[57] + p[30].*u[25] + (-1)*(p[80].*u[6].*u[56]) + (-1)*(p[29].*u[23].*u[6])
du[7] = p[33].*u[29] + (-1)*(p[32].*u[27].*u[7])
du[8] = p[40].*u[30] + p[86].*u[58] + p[88].*u[54] + p[39].*u[30] + (-1)*(p[85].*u[56].*u[8]) + (-1)*(p[87].*u[52].*u[8]) + (-1)*(p[38].*u[27].*u[8])
du[9] = p[42].*u[31] + (-1)*(p[41].*u[27].*u[9])
du[10] = p[45].*u[26] + (-1)*(p[44].*u[23].*u[10])
du[11] = p[48].*u[32] + (-1)*(p[47].*u[11].*u[28])
du[12] = p[50].*u[33] + (-1)*(p[49].*u[12].*u[28])
du[13] = p[55].*u[40] + p[59].*u[42] + p[63].*u[44] + (-1)*(p[54].*u[13].*u[39]) + (-1)*(p[58].*u[13].*u[41]) + (-1)*(p[62].*u[13].*u[43])
du[14] = p[65].*u[45] + (-1)*(p[64].*u[14].*u[43])
du[15] = p[68].*u[47] + (-1)*(p[67].*u[15].*u[46])
du[16] = p[71].*u[48] + (-1)*(p[70].*u[16].*u[46])
du[17] = p[76].*u[53] + (-1)*(p[75].*u[17].*u[51])
du[18] = p[79].*u[56] + (-1)*(p[78].*u[18].*u[55])
du[19] = p[19].*u[1].*u[2] + (-1)*(p[21].*u[19]) + (-1)*(p[20].*u[19])
du[20] = p[21].*u[19] + p[26].*u[22] + p[23].*u[21] + p[25].*u[22] + (-1)*(p[22].*u[20].*u[3]) + (-1)*(p[24].*u[20].*u[4])
du[21] = p[22].*u[20].*u[3] + (-1)*(p[23].*u[21])
du[22] = p[24].*u[20].*u[4] + (-1)*(p[26].*u[22]) + (-1)*(p[25].*u[22])
du[23] = p[46].*u[26] + p[26].*u[22] + p[31].*u[25] + p[37].*u[38] + p[45].*u[26] + p[28].*u[24] + p[30].*u[25] + (-1)*(p[44].*u[23].*u[10]) + (-1)*(p[27].*u[23].*u[5]) + (-1)*(p[29].*u[23].*u[6])
du[24] = p[27].*u[23].*u[5] + (-1)*(p[28].*u[24])
du[25] = p[29].*u[23].*u[6] + (-1)*(p[31].*u[25]) + (-1)*(p[30].*u[25])
du[26] = p[44].*u[23].*u[10] + (-1)*(p[46].*u[26]) + (-1)*(p[45].*u[26])
du[27] = p[82].*u[57] + p[31].*u[25] + p[34].*u[29] + p[43].*u[31] + p[33].*u[29] + p[39].*u[30] + p[42].*u[31] + (-1)*(p[32].*u[27].*u[7]) + (-1)*(p[38].*u[27].*u[8]) + (-1)*(p[41].*u[27].*u[9])
du[28] = p[46].*u[26] + p[51].*u[33] + p[48].*u[32] + p[50].*u[33] + (-1)*(p[47].*u[11].*u[28]) + (-1)*(p[49].*u[12].*u[28])
du[29] = p[32].*u[27].*u[7] + (-1)*(p[34].*u[29]) + (-1)*(p[33].*u[29])
du[30] = p[38].*u[27].*u[8] + (-1)*(p[40].*u[30]) + (-1)*(p[39].*u[30])
du[31] = p[41].*u[27].*u[9] + (-1)*(p[43].*u[31]) + (-1)*(p[42].*u[31])
du[32] = p[47].*u[11].*u[28] + (-1)*(p[48].*u[32])
du[33] = p[49].*u[12].*u[28] + (-1)*(p[51].*u[33]) + (-1)*(p[50].*u[33])
du[34] = p[34].*u[29] + p[37].*u[38] + p[36].*u[38] + (-1)*(p[35].*u[4].*u[34])
du[35] = p[40].*u[30]
du[36] = p[43].*u[31]
du[37] = p[51].*u[33] + p[53].*u[39] + (-1)*(p[52].*u[37])
du[38] = p[35].*u[4].*u[34] + (-1)*(p[37].*u[38]) + (-1)*(p[36].*u[38])
du[39] = p[52].*u[37] + p[55].*u[40] + (-2)*(0.5*p[56].*u[39].^2) + (-1)*(p[53].*u[39]) + 2*(p[57].*u[41]) + (-1)*(p[54].*u[13].*u[39])
du[40] = p[54].*u[13].*u[39] + (-1)*(p[55].*u[40])
du[41] = 0.5*p[56].*u[39].^2 + p[59].*u[42] + (-2)*(0.5*p[60].*u[41].^2) + (-1)*(p[57].*u[41]) + 2*(p[61].*u[43]) + (-1)*(p[58].*u[13].*u[41])
du[42] = p[58].*u[13].*u[41] + (-1)*(p[59].*u[42])
du[43] = 0.5*p[60].*u[41].^2 + p[63].*u[44] + p[65].*u[45] + (-1)*(p[61].*u[43]) + (-1)*(p[62].*u[13].*u[43]) + (-1)*(p[64].*u[14].*u[43])
du[44] = p[62].*u[13].*u[43] + (-1)*(p[63].*u[44])
du[45] = p[64].*u[14].*u[43] + (-1)*(p[66].*u[45]) + (-1)*(p[65].*u[45])
du[46] = p[66].*u[45] + p[69].*u[47] + p[72].*u[48] + p[68].*u[47] + p[71].*u[48] + (-1)*(p[67].*u[15].*u[46]) + (-1)*(p[70].*u[16].*u[46])
du[47] = p[67].*u[15].*u[46] + (-1)*(p[69].*u[47]) + (-1)*(p[68].*u[47])
du[48] = p[70].*u[16].*u[46] + (-1)*(p[72].*u[48]) + (-1)*(p[71].*u[48])
du[49] = p[69].*u[47] + p[74].*u[51] + (-1)*(p[73].*u[49])
du[50] = p[72].*u[48] + p[84].*u[52] + (-1)*(p[83].*u[50])
du[51] = p[77].*u[53] + p[73].*u[49] + p[76].*u[53] + (-1)*(p[74].*u[51]) + (-1)*(p[75].*u[17].*u[51])
du[52] = p[83].*u[50] + p[88].*u[54] + (-1)*(p[84].*u[52]) + (-1)*(p[87].*u[52].*u[8])
du[53] = p[75].*u[17].*u[51] + (-1)*(p[77].*u[53]) + (-1)*(p[76].*u[53])
du[54] = p[87].*u[52].*u[8] + (-1)*(p[88].*u[54])
du[55] = p[77].*u[53] + p[79].*u[56] + (-1)*(p[78].*u[18].*u[55])
du[56] = p[82].*u[57] + p[78].*u[18].*u[55] + p[81].*u[57] + p[86].*u[58] + (-1)*(p[79].*u[56]) + (-1)*(p[80].*u[6].*u[56]) + (-1)*(p[85].*u[56].*u[8])
du[57] = p[80].*u[6].*u[56] + (-1)*(p[82].*u[57]) + (-1)*(p[81].*u[57])
du[58] = p[85].*u[56].*u[8] + (-1)*(p[86].*u[58])
nothing
end'''

tspan = [0.0, 20000.0]
initials = [3.e+03, 2.e+02, 1.e+02, 2.e+04, 1.e+03, 1.e+04, 1.e+04, 1.e+05,
        1.e+06, 4.e+04, 2.e+04, 1.e+05, 2.e+04, 5.e+05, 5.e+05, 1.e+05,
        1.e+05, 1.e+05, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00]

params = [3.00000000e+03, 2.00000000e+02, 1.00000000e+02, 2.00000000e+04,
        1.00000000e+03, 1.00000000e+04, 1.00000000e+04, 1.00000000e+05,
        1.00000000e+06, 4.00000000e+04, 2.00000000e+04, 1.00000000e+05,
        2.00000000e+04, 5.00000000e+05, 5.00000000e+05, 1.00000000e+05,
        1.00000000e+05, 1.00000000e+05, 4.00000000e-07, 1.00000000e-03,
        1.00000000e-05, 1.00000000e-06, 1.00000000e-03, 1.00000000e-06,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-06, 1.00000000e-03,
        1.00000000e-07, 1.00000000e-03, 1.00000000e+00, 1.00000000e-06,
        1.00000000e-03, 1.00000000e+00, 3.00000000e-08, 1.00000000e-03,
        1.00000000e+00, 2.00000000e-06, 1.00000000e-03, 1.00000000e-01,
        1.00000000e-06, 1.00000000e-02, 1.00000000e+00, 1.00000000e-07,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-06, 1.00000000e-03,
        1.00000000e-07, 1.00000000e-03, 1.00000000e+00, 1.00000000e-02,
        1.00000000e-02, 1.42857143e-05, 1.00000000e-03, 2.85714286e-05,
        1.00000000e-03, 1.42857143e-05, 1.00000000e-03, 2.85714286e-05,
        1.00000000e-03, 1.42857143e-05, 1.00000000e-03, 1.42857143e-05,
        1.00000000e-03, 1.00000000e+00, 2.85714286e-05, 1.00000000e-03,
        1.00000000e+01, 2.85714286e-05, 1.00000000e-03, 1.00000000e+01,
        1.00000000e-02, 1.00000000e-02, 5.00000000e-07, 1.00000000e-03,
        1.00000000e+00, 5.00000000e-08, 1.00000000e-03, 5.00000000e-09,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-02, 1.00000000e-02,
        2.00000000e-06, 1.00000000e-03, 7.00000000e-06, 1.00000000e-03]

jul_fn = Main.eval(jul_code)

prob = de.ODEProblem(jul_fn, initials,
                     tspan,
                     params)
result = de.solve(prob)

julia_equivalent.jl

using DifferentialEquations

function earm(du,u,p,t)
du[1] = p[20].*u[19] + (-1)*(p[19].*u[1].*u[2])
du[2] = p[20].*u[19] + (-1)*(p[19].*u[1].*u[2])
du[3] = p[23].*u[21] + (-1)*(p[22].*u[20].*u[3])
du[4] = p[25].*u[22] + p[36].*u[38] + (-1)*(p[24].*u[20].*u[4]) + (-1)*(p[35].*u[4].*u[34])
du[5] = p[28].*u[24] + (-1)*(p[27].*u[23].*u[5])
du[6] = p[81].*u[57] + p[30].*u[25] + (-1)*(p[80].*u[6].*u[56]) + (-1)*(p[29].*u[23].*u[6])
du[7] = p[33].*u[29] + (-1)*(p[32].*u[27].*u[7])
du[8] = p[40].*u[30] + p[86].*u[58] + p[88].*u[54] + p[39].*u[30] + (-1)*(p[85].*u[56].*u[8]) + (-1)*(p[87].*u[52].*u[8]) + (-1)*(p[38].*u[27].*u[8])
du[9] = p[42].*u[31] + (-1)*(p[41].*u[27].*u[9])
du[10] = p[45].*u[26] + (-1)*(p[44].*u[23].*u[10])
du[11] = p[48].*u[32] + (-1)*(p[47].*u[11].*u[28])
du[12] = p[50].*u[33] + (-1)*(p[49].*u[12].*u[28])
du[13] = p[55].*u[40] + p[59].*u[42] + p[63].*u[44] + (-1)*(p[54].*u[13].*u[39]) + (-1)*(p[58].*u[13].*u[41]) + (-1)*(p[62].*u[13].*u[43])
du[14] = p[65].*u[45] + (-1)*(p[64].*u[14].*u[43])
du[15] = p[68].*u[47] + (-1)*(p[67].*u[15].*u[46])
du[16] = p[71].*u[48] + (-1)*(p[70].*u[16].*u[46])
du[17] = p[76].*u[53] + (-1)*(p[75].*u[17].*u[51])
du[18] = p[79].*u[56] + (-1)*(p[78].*u[18].*u[55])
du[19] = p[19].*u[1].*u[2] + (-1)*(p[21].*u[19]) + (-1)*(p[20].*u[19])
du[20] = p[21].*u[19] + p[26].*u[22] + p[23].*u[21] + p[25].*u[22] + (-1)*(p[22].*u[20].*u[3]) + (-1)*(p[24].*u[20].*u[4])
du[21] = p[22].*u[20].*u[3] + (-1)*(p[23].*u[21])
du[22] = p[24].*u[20].*u[4] + (-1)*(p[26].*u[22]) + (-1)*(p[25].*u[22])
du[23] = p[46].*u[26] + p[26].*u[22] + p[31].*u[25] + p[37].*u[38] + p[45].*u[26] + p[28].*u[24] + p[30].*u[25] + (-1)*(p[44].*u[23].*u[10]) + (-1)*(p[27].*u[23].*u[5]) + (-1)*(p[29].*u[23].*u[6])
du[24] = p[27].*u[23].*u[5] + (-1)*(p[28].*u[24])
du[25] = p[29].*u[23].*u[6] + (-1)*(p[31].*u[25]) + (-1)*(p[30].*u[25])
du[26] = p[44].*u[23].*u[10] + (-1)*(p[46].*u[26]) + (-1)*(p[45].*u[26])
du[27] = p[82].*u[57] + p[31].*u[25] + p[34].*u[29] + p[43].*u[31] + p[33].*u[29] + p[39].*u[30] + p[42].*u[31] + (-1)*(p[32].*u[27].*u[7]) + (-1)*(p[38].*u[27].*u[8]) + (-1)*(p[41].*u[27].*u[9])
du[28] = p[46].*u[26] + p[51].*u[33] + p[48].*u[32] + p[50].*u[33] + (-1)*(p[47].*u[11].*u[28]) + (-1)*(p[49].*u[12].*u[28])
du[29] = p[32].*u[27].*u[7] + (-1)*(p[34].*u[29]) + (-1)*(p[33].*u[29])
du[30] = p[38].*u[27].*u[8] + (-1)*(p[40].*u[30]) + (-1)*(p[39].*u[30])
du[31] = p[41].*u[27].*u[9] + (-1)*(p[43].*u[31]) + (-1)*(p[42].*u[31])
du[32] = p[47].*u[11].*u[28] + (-1)*(p[48].*u[32])
du[33] = p[49].*u[12].*u[28] + (-1)*(p[51].*u[33]) + (-1)*(p[50].*u[33])
du[34] = p[34].*u[29] + p[37].*u[38] + p[36].*u[38] + (-1)*(p[35].*u[4].*u[34])
du[35] = p[40].*u[30]
du[36] = p[43].*u[31]
du[37] = p[51].*u[33] + p[53].*u[39] + (-1)*(p[52].*u[37])
du[38] = p[35].*u[4].*u[34] + (-1)*(p[37].*u[38]) + (-1)*(p[36].*u[38])
du[39] = p[52].*u[37] + p[55].*u[40] + (-2)*(0.5*p[56].*u[39].^2) + (-1)*(p[53].*u[39]) + 2*(p[57].*u[41]) + (-1)*(p[54].*u[13].*u[39])
du[40] = p[54].*u[13].*u[39] + (-1)*(p[55].*u[40])
du[41] = 0.5*p[56].*u[39].^2 + p[59].*u[42] + (-2)*(0.5*p[60].*u[41].^2) + (-1)*(p[57].*u[41]) + 2*(p[61].*u[43]) + (-1)*(p[58].*u[13].*u[41])
du[42] = p[58].*u[13].*u[41] + (-1)*(p[59].*u[42])
du[43] = 0.5*p[60].*u[41].^2 + p[63].*u[44] + p[65].*u[45] + (-1)*(p[61].*u[43]) + (-1)*(p[62].*u[13].*u[43]) + (-1)*(p[64].*u[14].*u[43])
du[44] = p[62].*u[13].*u[43] + (-1)*(p[63].*u[44])
du[45] = p[64].*u[14].*u[43] + (-1)*(p[66].*u[45]) + (-1)*(p[65].*u[45])
du[46] = p[66].*u[45] + p[69].*u[47] + p[72].*u[48] + p[68].*u[47] + p[71].*u[48] + (-1)*(p[67].*u[15].*u[46]) + (-1)*(p[70].*u[16].*u[46])
du[47] = p[67].*u[15].*u[46] + (-1)*(p[69].*u[47]) + (-1)*(p[68].*u[47])
du[48] = p[70].*u[16].*u[46] + (-1)*(p[72].*u[48]) + (-1)*(p[71].*u[48])
du[49] = p[69].*u[47] + p[74].*u[51] + (-1)*(p[73].*u[49])
du[50] = p[72].*u[48] + p[84].*u[52] + (-1)*(p[83].*u[50])
du[51] = p[77].*u[53] + p[73].*u[49] + p[76].*u[53] + (-1)*(p[74].*u[51]) + (-1)*(p[75].*u[17].*u[51])
du[52] = p[83].*u[50] + p[88].*u[54] + (-1)*(p[84].*u[52]) + (-1)*(p[87].*u[52].*u[8])
du[53] = p[75].*u[17].*u[51] + (-1)*(p[77].*u[53]) + (-1)*(p[76].*u[53])
du[54] = p[87].*u[52].*u[8] + (-1)*(p[88].*u[54])
du[55] = p[77].*u[53] + p[79].*u[56] + (-1)*(p[78].*u[18].*u[55])
du[56] = p[82].*u[57] + p[78].*u[18].*u[55] + p[81].*u[57] + p[86].*u[58] + (-1)*(p[79].*u[56]) + (-1)*(p[80].*u[6].*u[56]) + (-1)*(p[85].*u[56].*u[8])
du[57] = p[80].*u[6].*u[56] + (-1)*(p[82].*u[57]) + (-1)*(p[81].*u[57])
du[58] = p[85].*u[56].*u[8] + (-1)*(p[86].*u[58])
nothing
end

tspan = (0.0, 20000.0)
initials = [3.e+03, 2.e+02, 1.e+02, 2.e+04, 1.e+03, 1.e+04, 1.e+04, 1.e+05,
        1.e+06, 4.e+04, 2.e+04, 1.e+05, 2.e+04, 5.e+05, 5.e+05, 1.e+05,
        1.e+05, 1.e+05, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00, 0.e+00,
        0.e+00, 0.e+00]
params = [3.00000000e+03, 2.00000000e+02, 1.00000000e+02, 2.00000000e+04,
        1.00000000e+03, 1.00000000e+04, 1.00000000e+04, 1.00000000e+05,
        1.00000000e+06, 4.00000000e+04, 2.00000000e+04, 1.00000000e+05,
        2.00000000e+04, 5.00000000e+05, 5.00000000e+05, 1.00000000e+05,
        1.00000000e+05, 1.00000000e+05, 4.00000000e-07, 1.00000000e-03,
        1.00000000e-05, 1.00000000e-06, 1.00000000e-03, 1.00000000e-06,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-06, 1.00000000e-03,
        1.00000000e-07, 1.00000000e-03, 1.00000000e+00, 1.00000000e-06,
        1.00000000e-03, 1.00000000e+00, 3.00000000e-08, 1.00000000e-03,
        1.00000000e+00, 2.00000000e-06, 1.00000000e-03, 1.00000000e-01,
        1.00000000e-06, 1.00000000e-02, 1.00000000e+00, 1.00000000e-07,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-06, 1.00000000e-03,
        1.00000000e-07, 1.00000000e-03, 1.00000000e+00, 1.00000000e-02,
        1.00000000e-02, 1.42857143e-05, 1.00000000e-03, 2.85714286e-05,
        1.00000000e-03, 1.42857143e-05, 1.00000000e-03, 2.85714286e-05,
        1.00000000e-03, 1.42857143e-05, 1.00000000e-03, 1.42857143e-05,
        1.00000000e-03, 1.00000000e+00, 2.85714286e-05, 1.00000000e-03,
        1.00000000e+01, 2.85714286e-05, 1.00000000e-03, 1.00000000e+01,
        1.00000000e-02, 1.00000000e-02, 5.00000000e-07, 1.00000000e-03,
        1.00000000e+00, 5.00000000e-08, 1.00000000e-03, 5.00000000e-09,
        1.00000000e-03, 1.00000000e+00, 1.00000000e-02, 1.00000000e-02,
        2.00000000e-06, 1.00000000e-03, 7.00000000e-06, 1.00000000e-03]
prob = ODEProblem(earm, initials, tspan, params)
sol = solve(prob)