Suppressing analytical evaluation in NDSolve

Question

I am trying to find a numerical solution to a set of coupled ODEs with NDSolve. Let us say $\boldsymbol{X}$ is a vector, and $\boldsymbol{F}$ is a non-linear vector function of $\boldsymbol{X}$: $$\dot{\boldsymbol{X}}=\boldsymbol{F}(\boldsymbol{X}),\qquad\dot{\boldsymbol{X}}(0)=\boldsymbol{X}_{0}$$.

In a nutshell, my problem is this: $\boldsymbol{F}$ is basically infeasible to evaluate analytically (if you pass it a symbolic list $\left\{ X_{1},X_{2},\ldots,X_{n}\right\} $ as an argument). However, if you pass it a list of numbers, $\boldsymbol{F}$ is no problem to evaluate. NDSolve seems to try to evaluate things analytically (I think), preventing it from calculating a result in feasible time. I'm looking for a way to tell NDSolve to stop trying any analytical evaluation of $\boldsymbol{F}$.

Here is an example just to illustrate the point. Let's say $\boldsymbol{F}$ involves some nasty matrix inversions: $$\boldsymbol{F}(\boldsymbol{X})=\left(\boldsymbol{1}+\left(\boldsymbol{1}+\boldsymbol{X}\otimes\boldsymbol{X}\right)^{-1}\right)^{-1}\boldsymbol{X}$$

ndim = 10;
Xvec[t_] = Table[x[k][t], {k, 1, ndim}];
F[xvec_] := 
 Nest[Inverse[IdentityMatrix[ndim] + #] &, TensorProduct[xvec, xvec], 
   2].xvec

In ndim = 10 dimensions, if you pass a symbolic argument just evaluating $\boldsymbol{F}$ takes longer than I have patience for:

    AbsoluteTiming[F[Xvec[t]];]
(* takes too long! *)

If instead you pass a numerical argument, it takes no time at all to evaluate:

AbsoluteTiming[F[RandomReal[{-5, 5}, ndim]]]

(* OUTPUT: {0.0000794, {3.72758, -2.68105, -1.71283, -3.62467, -3.39943, \
-1.28462, 2.47934, 0.216207, -0.944914, 1.36631}} *)

I think this problem carries over to NDSolve. The following calculation does not finish in feasible time for ndim = 10:

Xvec0 = RandomReal[{-5, 5}, ndim]; (* initial conditions *)
tend = 10;
sol = NDSolve[{Thread[Xvec'[t] == F[Xvec[t]]], 
   Thread[Xvec[0] == Xvec0]}, Xvec[t], {t, 0, tend}]

(* takes too long! *)

Is there a way I can instruct NDSolve to suppress analytical calculation and treat $\boldsymbol{F}$ with numerical arguments?

An unsuccessful solution attempt with NumericQ

I tried defining F with a pattern that checks for vector arguments (VectorQ) and that the argument is numerical (NumericQ):

ndim = 10;
Xvec[t_] = Table[x[k][t], {k, 1, ndim}];
F[xvec_?(VectorQ[#, NumericQ] &)] := 
 Module[{}, 
  Nest[Inverse[IdentityMatrix[ndim] + #] &, 
    TensorProduct[xvec, xvec], 2].xvec]

The behaviour is as expected: Now, symbolic arguments no longer evaluate

AbsoluteTiming[F[Xvec[t]]]

(* OUTPUT: {0.0000106, 
 F[{x[1][t], x[2][t], x[3][t], x[4][t], x[5][t], x[6][t], x[7][t], 
   x[8][t], x[9][t], x[10][t]}]} *)

... and numerical arguments evaluate quickly, just as previously.

However, NDSolve does not evaluate F, sadly:

Xvec0 = RandomReal[{-5, 5}, ndim]; (* initial conditions *)
tend = 10;
sol = NDSolve[{Thread[Xvec'[t] == F[Xvec[t]]], 
   Thread[Xvec[0] == Xvec0]}, Xvec[t], {t, 0, tend}]

(* OUTPUT: NDSolve[{{Derivative[1][x[1]][t] == F[{x[1][t], x[2][t], x[3][t], ... *)

I am hoping for some suggestions to make NDSolve calculate a solution in feasible time.

Answer 1

It is strange that NDSolve does not produce warning messages but the reason behind the failure is that F[{x[1][t], ...}] does not evaluate to an explicit List when supplied a symbolic argument. So, when you evaluate

Thread[Xvec'[t] == F[Xvec[t]]]

you end up with all x[k]'[t] equated to the same F[{x[1][t], ...}], each of which will later evaluate to a vector when supplied a numeric argument. This makes NDSolve confused because you have scalars in the initial conditions for x[k][0].

You don't have to use explicit vectors at all. NDSolve can treat an array as a single point:

ClearAll[ndim, Xvec0, F, x];
ndim = 10;
Xvec0 = RandomReal[{-5, 5}, ndim];
F[xvec_?(VectorQ[#, NumericQ] &)] := 
Module[{}, 
  Nest[Inverse[IdentityMatrix[ndim] + #] &, TensorProduct[xvec, xvec], 2] . xvec]
NDSolve[{x'[t] == F[x[t]], x[0] == Xvec0}, x, {t, 0, 1}]