Message Boards Message Boards

Does NDSolve has an upper limit requirement of the number of ODE equations?

GROUPS:

I founded a series of master equations(equation set) in Fe-Cu-Ni system, where I should use NDSolve function to solve these equations. The problem I encountered was that NDSolve raised an error when the number of equations in the set was over about 10000. The error said: enter image description here So, does NDSolve has an upper limit requirement of the number of ODE equations? Or how should I fix the problem? My code is in the attachment. Looking forward your help and much thanks!

POSTED BY: Lammond Wan
Answer
20 days ago

NDSolve tries to solve for the highest derivatives (or the derivatives of the equivalent first-order system). If solving fails, then you get the error. Try solving the system yourself. There may also be a time constraint on solving, but I'm not sure.

POSTED BY: Michael Rogers
Answer
20 days ago

As you can see in my code, the founded equations are all of first order, but its amount may reach to 200 000 if needed. The error occurs when the amount of ODEs gets greater than only 10 000. It really disappoints me. By the way, what do you mean by mentioning solving it myself? Is there other effective way to solve ODE set? I'm a beginner of MMA, sorry to ask you this question. Thanks!

POSTED BY: Lammond Wan
Answer
20 days ago

Try add option to NDSolve:

 Method->{"EquationSimplification"->"Residual"}

Did it help?

POSTED BY: Mariusz Iwaniuk
Answer
20 days ago

Wow, the error disappears! I will check the correctness of the answers later, but I need to thank you at first. Moreover, how does it works by setting the method of equation simplification?

POSTED BY: Lammond Wan
Answer
20 days ago

You can find the answer how it works probably here -> Differential-Algebraic Equations -> Introduction

Regards,MI

POSTED BY: Mariusz Iwaniuk
Answer
20 days ago

Thank you very much, it's very useful and helpful. I guess the key to the problem is that it will take a long long time to get a explicit formula when the number of equations gets so big, therefore residual form is a better way to solve the question. Thanks a lot.

POSTED BY: Lammond Wan
Answer
19 days ago

By solving, I mean something like

Solve[Sys, D[Through[Bc[t]], t]]

As for "As you can see in my code,..." actually, no I couldn't see. I didn't try. Do you really expect someone to closely analyze so much code and the 10000 equations you said it produced?

POSTED BY: Michael Rogers
Answer
19 days ago

Sorry to express it that way, my English is so poor that I can't describe what I mean to say correctly. No doubt I should fix my problem myself, and I'm really appreciate your help. Thank you very much!

POSTED BY: Lammond Wan
Answer
19 days ago

No problem, no worries.

I had an idea, which should work in principle. I'm pretty sure it fails because I ran out of memory. It actually almost completely freezes my laptop, and I have to force-quit Mathematica and several other applications to get back to where I can use it normally. That makes me think that the problem with your code is similar, except that your code fails semi-gracefully (no error reported, but does not freeze computer). Often running out of memory results in an error message, but I recall times when it did not.

I have only 16GB, but if you have more, the idea is to recast the system in terms of a vector variable cc comprising all your c[i, j] variables:

ca = CoefficientArrays[Sys[[;; Length@Bc, 2]], Through[Bc[t]]]
VectorQ[#["NonzeroValues"], NumericQ] & /@ ca  (* check there are no variables in the coefficients *)
(*  {True, True, True}  *)

ivs = Sys[[Length@Bc + 1 ;;, 2]];  (* initial values *)
VectorQ[ivs, NumericQ] (* check there are no variables in the values *)
(*  True  *)

MemoryConstrained[
 NDSolve[{
  cc'[t] == Fold[#1.cc[t] + #2 &, Reverse@ca], (* reconstruct polynomial RHS of the ode system (Horner's method) *)
  cc[0] == ivs}, 
  cc, {t, tMin, tMax; 0.001}(*Method->"StiffnessSwitching", MaxSteps -> Infinity*)],
  4*^9,   (* put your memory limit here *)
  "ran out of memory"]
(*  "ran out of memory"  *)

Should you be able to get a solution this way, the output value of the interpolation function will be a vector whose components correspond to the variables in Bc in the same order.

POSTED BY: Michael Rogers
Answer
17 days ago

Group Abstract Group Abstract