9
$\begingroup$

Consider the following test case:

mMax = 50;
Ffuncs[r_] = Table[F[n, m][r], {n, mMax}, {m, mMax}];
NDSolve[Flatten@{
    Thread[Flatten[D[Ffuncs[r], r]] == Flatten[Ffuncs[r].Ffuncs[r]]],
    Thread[Flatten@Ffuncs[0] == Table[0, {mMax*mMax}]]
    }, Flatten@Ffuncs[r], {r, 0, 2}];

When I set mMax to 49 or less, it works without any problems. But with 50 or higher, I get

NDSolve::ntdv: Cannot solve to find an explicit formula for the derivatives. Consider using the option Method->{"EquationSimplification"->"Residual"}.

That was Mathematica 11. Mathematica 9 gives a better message and a better suggestion, which actually works, unlike the suggestion above, which takes forever to complete the solution and eventually crashes the kernel. Here's Mathematica 9's message:

NDSolve::ndsdtc: The time constraint of 1.` seconds was exceeded trying to solve for derivatives, so the system will be treated as a system of differential-algebraic equations. You can use Method->{"EquationSimplification"->"Solve"} to have the system solved as ordinary differential equations. >>

But anyway, I've already gave NDSolve a system of equations explicitly solved for derivatives! Why does it still need to solve for derivatives?

Improve this question
$\endgroup$
0

2 Answers 2

Reset to default
7
$\begingroup$

In V10 and higher, you can increase the time limit on Solve with the system option "NDSolveOptions" -> "DefaultSolveTimeConstraint" -> time:

With[{opts = SystemOptions[]},
 Internal`WithLocalSettings[
  SetSystemOptions["NDSolveOptions" -> "DefaultSolveTimeConstraint" -> 10.`],
  sol = NDSolve[
     Flatten@{Thread[
        Flatten[D[Ffuncs[r], r]] == Flatten[Ffuncs[r].Ffuncs[r]]], 
       Thread[Flatten@Ffuncs[0] == Table[0, {mMax*mMax}]]}, 
     Flatten@Ffuncs[r], {r, 0, 2}];,
  SetSystemOptions[opts]
  ]]
(* Spurious message: SystemOptions::noset: ..SystemOption PreemptiveCheckUseThreads..  *)

sol[[1, 1]]
(*  F[1, 1][r] -> InterpolatingFunction[{{0., 2.}}, <>][r]  *)

I'm not sure what takes NDSolve[] so long, but it's much faster than Solve[]:

Solve[Thread[Flatten[D[Ffuncs[r], r]] == Flatten[Ffuncs[r].Ffuncs[r]]], 
   Flatten[D[Ffuncs[r], r]]]; // AbsoluteTiming
(*  {51.1422, Null}  *)
Improve this answer
$\endgroup$
8
  • 3
    $\begingroup$ But why does NDSolve even need to solve anything algebraically if the derivatives are already given to it? $\endgroup$
    – Ruslan
    Nov 18, 2016 at 18:55
  • 1
    $\begingroup$ @Ruslan NDSolve constructs a function that represents the "right-hand side" of the system ${\bf x}'(r) = f({\bf x}, r)$. It does that by solving the equations you pass it (I surmise). You don't actually pass derivatives to NDSolve. Each expression has Equal as its head. It has to inspect them somehow to see what's there. I don't understand why it should take so long, though. $\endgroup$
    – Michael E2
    Nov 18, 2016 at 20:24
  • $\begingroup$ Note that this doesn't work in Mathematica 9 (it continues to complain about 1 second being too small time). $\endgroup$
    – Ruslan
    Dec 8, 2016 at 6:08
  • $\begingroup$ @Ruslan Hmm, the code works for me without the time constraint option, which doesn't seem to have existed in V9: i.stack.imgur.com/4gvHh.png -- But maybe it's solving fast enough on my machine? $\endgroup$
    – Michael E2
    Dec 8, 2016 at 11:21
  • 1
    $\begingroup$ @Ruslan Yes, that was it. The workaround Method -> {"EquationSimplification" -> "Solve"} still works. $\endgroup$
    – Michael E2
    Dec 8, 2016 at 12:43
3
$\begingroup$

It helps to pose this in vector form, see: Vector form using NDSolve

This at least doesn't throw the solving for derivatives error. (tested up to 200..)

mMax = 60;
vf[vals : {_?NumberQ ..}] := 
 Module[{m = Partition[vals, mMax]}, Flatten[m.m]]

vsoln = NDSolveValue[{x'[r] == vf[x[r]], 
     x[0] == Table[RandomReal[.1], {mMax^2}]}, x[r], {r, 0, 1}];

or with the square matix as the unknowm:

mMax = 60;
vf[vals : {{_?NumberQ ..} ..}] := vals.vals
vsoln = NDSolveValue[{x'[r] == vf[x[r]], 
    x[0] == Table[RandomReal[{-1, 1}], {mMax}, {mMax}]}, 
    x[r], {r, 0, 1}];

I threw in a nonzero initial condition so you get a nonzero result. you should verify this is doing the right thing, but I think its correct.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.