Is it possible to solve this system of ODEs using Mathematica?

Question

I have a system of 6 1st order ODEs. Since, I don't have access to a Mathematica license yet, I have tried to solve the problem using Maxima but it failed. I have also tried Maple but it takes several hours then the software simply crashes.

I wonder if Mathematica could solve this system of ODEs:

diff(f(t),t) = -15.722 -89.755*f(t) + 1.192*x(t) -20.522*h(t) -101.456*f(t)^2 + 0.939*f(t)*x(t) + 41.345*f(t)*g(t) -70.537*f(t)*h(t) -0.271*x(t)*g(t) + 4.246*g(t)^2 + 8.767*g(t)*h(t) -3.25*h(t)^2,
diff(x(t),t) = -81429.279 + 114841.617*f(t) -549.97*x(t) + 149785.879*g(t) -599.101*y(t) -160497.335*h(t) -27.319*z(t) + 1721899.924*f(t)^2 -1117.648*f(t)*x(t) -386652.825*f(t)*g(t) -3087.313*f(t)*y(t) -2368.155*f(t)*h(t) + 425.026*f(t)*z(t) -3.547*x(t)^2 + 537.27*x(t)*g(t) + 1.582*x(t)*y(t) -805.526*x(t)*h(t) + 0.028*x(t)*z(t) -25833.096*g(t)^2 + 597.113*g(t)*y(t) + 97189.995*g(t)*h(t) -145.815*g(t)*z(t) -0.315*y(t)^2 -157.341*y(t)*h(t) -77306.535*h(t)^2 -39.105*h(t)*z(t) -0.37*z(t)^2,
diff(g(t),t) = 24.974 + 270.543*f(t) -0.304*x(t) -87.077*g(t) + 1.324*y(t) -1.43*h(t) + 682.7*f(t)^2 -1.172*f(t)*x(t) -448.244*f(t)*g(t) + 0.266*f(t)*y(t) + 7.795*f(t)*h(t) + 0.251*f(t)*z(t) + 54.566*g(t)^2 -0.171*g(t)*y(t) -15.812*g(t)*h(t) -10.325*h(t)^2 -0.18*h(t)*z(t),
diff(y(t),t) = 558255.344 + 3927732.918*f(t) -16311.34*x(t) -368818.217*g(t) -3562.042*y(t) -431711.16*h(t) + 1291.352*z(t) + 8775396.002*f(t)^2 -43186.863*f(t)*x(t) -3982744.391*f(t)*g(t) -14420.954*f(t)*y(t) + 306905.157*f(t)*h(t) + 1493.098*f(t)*z(t) + 22.638*x(t)^2 + 13587.668*x(t)*g(t) + 14.649*x(t)*y(t) -4819.379*x(t)*h(t) -8.782*x(t)*z(t) + 238574.183*g(t)^2 + 3562.853*g(t)*y(t) + 279799.995*g(t)*h(t) -2146.381*g(t)*z(t) -1.764*y(t)^2 -1174.223*y(t)*h(t) -0.609*y(t)*z(t) -452034.181*h(t)^2 -736.964*h(t)*z(t) + 1.509*z(t)^2,
diff(h(t),t) = -14.036 -226.65*f(t) -0.349*x(t) + 20.774*g(t) + 0.598*y(t) -8.65*h(t) + 1.276*z(t) -773.841*f(t)^2 -1.983*f(t)*x(t) + 47.953*f(t)*g(t) + 1.535*f(t)*y(t) -107.657*f(t)*h(t) + 1.423*f(t)*z(t) + 0.221*x(t)*g(t) -0.185*x(t)*h(t) + 13.429*g(t)^2 -0.387*g(t)*y(t) + 6.67*g(t)*h(t) -0.408*g(t)*z(t) + 0.317*y(t)*h(t) -4.393*h(t)^2,
diff(z(t),t) = -512379.34 -3748635.593*f(t) -32866.886*x(t) + 858047.636*g(t) + 5987.248*y(t) -759005.177*h(t) -321.996*z(t) -704671.07*f(t)^2 -92042.657*f(t)*x(t) + 2548473.477*f(t)*g(t) + 11364.094*f(t)*y(t) -490421.909*f(t)*h(t) -1747.524*f(t)*z(t) -2.934*x(t)^2 + 24426.966*x(t)*g(t) -2.318*x(t)*y(t) -17822.416*x(t)*h(t) + 18.989*x(t)*z(t) -335819.401*g(t)^2 -3789.785*g(t)*y(t) + 599080.396*g(t)*h(t) + 352.73*g(t)*z(t) + 0.966*y(t)^2 + 2579.875*y(t)*h(t) -1.916*y(t)*z(t) -183880.537*h(t)^2 -1103.295*h(t)*z(t) -1.478*z(t)^2

I appreciate any help.

EDIT: The initial conditions:

t0 = 2.083195363,
f(t0)=2.7604518140e-02,
x(t0)=6.39998111386,
g(t0)=6.07690217883e-01,
y(t0)=1.149214543557e+02,
h(t0)=-1.1758467477,
z(t0)=-7.980181509e+01

Only with NDSolve ,but you need a initial or boundary conditions. — Mariusz Iwaniuk, Sep 21 '21 at 09:26
@user64494: Sorry, I have never used Mathematica before, I don't know the syntax yet. — Iyach tharwa nambarek, Sep 21 '21 at 10:04
@user64494 OK, I have seen that, what next? I don't have Mathematica, and I don't see anything in this website guidline that prevents such questions. — Iyach tharwa nambarek, Sep 21 '21 at 10:23
Please read What topics can I ask about here?. The guidelines say, "some proof of minimal Mathematica knowledge: be familiar with the language, its syntax, and the ways you can look for help." — creidhne, Sep 21 '21 at 10:31
@creidhne: How Am I supposed to know Mathematica if I don't have a license? — Iyach tharwa nambarek, Sep 21 '21 at 10:33
@creidhne: Any help on how can I convert that to Mathematica syntax? — Iyach tharwa nambarek, Sep 21 '21 at 10:51
You tried Maple and Maxima without success. Possibly your ode might be incorrect? — Ulrich Neumann, Sep 21 '21 at 10:52
@UlrichNeumann: What does it mean for an ODE to be incorrect? — Iyach tharwa nambarek, Sep 21 '21 at 10:53
What I mean is: Did you check your ode's , perhaps a little mistake causes the problems — Ulrich Neumann, Sep 21 '21 at 10:59
@UlrichNeumann: I have generate these ODEs using Sparse Identification of Nonlinear Dynamics method, when solved numerically there it looks ok. — Iyach tharwa nambarek, Sep 21 '21 at 11:03
@Signiycello That means you are looking for a symbolical solution? By the way, I tried to solve your odes numerically starting with t==2.0832 but Mathematica gives an error "NDSolveValue::ndsz: At t == 2.086949080697118`, step size is effectively zero; singularity or stiff system suspected." — Ulrich Neumann, Sep 21 '21 at 11:10
@UlrichNeumann: Could you please show the syntax so I can try it my self. I will try with other initial values: t0=2.132505169999999950e+00, f(t0)=3.011115631999999875e-02, x(t0)=9.535997718542907009e+00, g(t0)=9.042972289999999802e-01, y(t0)=1.004099163458574395e+02, h(t0)=-1.278683797100000019e+00, z(t0)=-1.960911326308590219e+01 — Iyach tharwa nambarek, Sep 21 '21 at 11:17
Also, nowadays you can use the free Wolfram engine: https://mathematica.stackexchange.com/q/198839/1871 But as mentioned by Ulrich, there's probably something wrong with your ODE system, you'd better double-check if you've made any simple mistake. — xzczd, Sep 21 '21 at 11:32
@xzczd: Thank you! I have tried wolframalpha but it didn't accept it due to the size limit. I hope the Wolfram Kernel within jupyter will not have the same limitations. — Iyach tharwa nambarek, Sep 21 '21 at 11:37

score 2 · Answer 1 · answered Sep 21 '21 at 11:24

As a special service ;-) here your problem in Mathematica syntax:

 odes = {f'[t] == -15.722 - 89.755*f[t] + 1.192*x[t] - 20.522*h[t] - 
     101.456*f[t]^2 + 0.939*f[t]*x[t] + 41.345*f[t]*g[t] - 
     70.537*f[t]*h[t] - 0.271*x[t]*g[t] + 4.246*g[t]^2 + 
     8.767*g[t]*h[t] - 3.25*h[t]^2, 
   x'[t] == -81429.279 + 114841.617*f[t] - 549.97*x[t] + 
     149785.879*g[t] - 599.101*y[t] - 160497.335*h[t] - 27.319*z[t] + 
     1721899.924*f[t]^2 - 1117.648*f[t]*x[t] - 386652.825*f[t]*g[t] - 
     3087.313*f[t]*y[t] - 2368.155*f[t]*h[t] + 425.026*f[t]*z[t] - 
     3.547*x[t]^2 + 537.27*x[t]*g[t] + 1.582*x[t]*y[t] - 
     805.526*x[t]*h[t] + 0.028*x[t]*z[t] - 25833.096*g[t]^2 + 
     597.113*g[t]*y[t] + 97189.995*g[t]*h[t] - 145.815*g[t]*z[t] - 
     0.315*y[t]^2 - 157.341*y[t]*h[t] - 77306.535*h[t]^2 - 
     39.105*h[t]*z[t] - 0.37*z[t]^2, 
   g'[t] == 
    24.974 + 270.543*f[t] - 0.304*x[t] - 87.077*g[t] + 1.324*y[t] - 
     1.43*h[t] + 682.7*f[t]^2 - 1.172*f[t]*x[t] - 448.244*f[t]*g[t] + 
     0.266*f[t]*y[t] + 7.795*f[t]*h[t] + 0.251*f[t]*z[t] + 
     54.566*g[t]^2 - 0.171*g[t]*y[t] - 15.812*g[t]*h[t] - 
     10.325*h[t]^2 - 0.18*h[t]*z[t], 
   y'[t] == 
    558255.344 + 3927732.918*f[t] - 16311.34*x[t] - 368818.217*g[t] - 
     3562.042*y[t] - 431711.16*h[t] + 1291.352*z[t] + 
     8775396.002*f[t]^2 - 43186.863*f[t]*x[t] - 
     3982744.391*f[t]*g[t] - 14420.954*f[t]*y[t] + 
     306905.157*f[t]*h[t] + 1493.098*f[t]*z[t] + 22.638*x[t]^2 + 
     13587.668*x[t]*g[t] + 14.649*x[t]*y[t] - 4819.379*x[t]*h[t] - 
     8.782*x[t]*z[t] + 238574.183*g[t]^2 + 3562.853*g[t]*y[t] + 
     279799.995*g[t]*h[t] - 2146.381*g[t]*z[t] - 1.764*y[t]^2 - 
     1174.223*y[t]*h[t] - 0.609*y[t]*z[t] - 452034.181*h[t]^2 - 
     736.964*h[t]*z[t] + 1.509*z[t]^2, 
   h'[t] == -14.036 - 226.65*f[t] - 0.349*x[t] + 20.774*g[t] + 
     0.598*y[t] - 8.65*h[t] + 1.276*z[t] - 773.841*f[t]^2 - 
     1.983*f[t]*x[t] + 47.953*f[t]*g[t] + 1.535*f[t]*y[t] - 
     107.657*f[t]*h[t] + 1.423*f[t]*z[t] + 0.221*x[t]*g[t] - 
     0.185*x[t]*h[t] + 13.429*g[t]^2 - 0.387*g[t]*y[t] + 
     6.67*g[t]*h[t] - 0.408*g[t]*z[t] + 0.317*y[t]*h[t] - 
     4.393*h[t]^2, 
   z'[t] == -512379.34 - 3748635.593*f[t] - 32866.886*x[t] + 
     858047.636*g[t] + 5987.248*y[t] - 759005.177*h[t] - 
     321.996*z[t] - 704671.07*f[t]^2 - 92042.657*f[t]*x[t] + 
     2548473.477*f[t]*g[t] + 11364.094*f[t]*y[t] - 
     490421.909*f[t]*h[t] - 1747.524*f[t]*z[t] - 2.934*x[t]^2 + 
     24426.966*x[t]*g[t] - 2.318*x[t]*y[t] - 17822.416*x[t]*h[t] + 
     18.989*x[t]*z[t] - 335819.401*g[t]^2 - 3789.785*g[t]*y[t] + 
     599080.396*g[t]*h[t] + 352.73*g[t]*z[t] + 0.966*y[t]^2 + 
     2579.875*y[t]*h[t] - 1.916*y[t]*z[t] - 183880.537*h[t]^2 - 
     1103.295*h[t]*z[t] - 1.478*z[t]^2} // Rationalize[#, 0] &

The initial conditions

t0 = 2.083195363
ics = {f[t0] == 2.7604518140 10^-02, x[t0] == 6.39998111386, 
   g[t0] == 6.07690217883 10^-01, y[t0] == 1.149214543557 10^+02, 
   h[t0] == -1.1758467477, z[t0] == -7.980181509 10^+01} // 
  Rationalize[#, 0] &

Numerical solution (attempt)

NDSolveValue[{odes, ics}, {f, x, g, y, h, z}, {t, t0, 10},Method -> {"EquationSimplification" -> "Residual"}]

returns error message "NDSolveValue::ndsz: At t == 2.086949080697118`, step size is effectively zero; singularity or stiff system suspected."

Hope it helps!

Thanks for you post. I got the same error when using t0=2.13... Apparently, Mathematica can't solve it as well. — Iyach tharwa nambarek, Sep 21 '21 at 11:31
These are the plots: https://i.imgur.com/wB84oz0.png and these are the equations I have obtained: https://i.imgur.com/PKKVwIy.png — Iyach tharwa nambarek, Sep 21 '21 at 11:42
Thanks, but my interest concerns the numerical methods you used. — Ulrich Neumann, Sep 21 '21 at 11:44

Is it possible to solve this system of ODEs using Mathematica?

1 Answers1