Matrix issue with NDSolve

Question

I am solving a differential matrix equation with NDSolve.

A = NDSolveValue[{y'[t] + (J\[Transpose]).y[t] + y[t].(J\[Transpose]) - 
 2*y[t].(J\[Transpose]).y[t] == 0, y[0] == \[Gamma]i}, {y[T]}, {t, 0, T}]

Here, initial value is a matrix. (This is the time to admit that I am very new to Mathematica.) First strange thing is that A appears to be not a matrix of the same size as initial was. It has three dimensions, first one is equal to 1, others are the same as for initial matrix.

Secondly, when I have a function g[y[t]] in the equation, it returns an error. The code is:

A = NDSolveValue[{y'[t] + (J\[Transpose] + g[y[t]]\[Transpose]).y[t] +
  y[t].(J\[Transpose] + g[y[t]]\[Transpose]) - 2*y[t].(J\[Transpose] + 
  g[y[t]]\[Transpose]).y[t] == 0, y[0] == \[Gamma]i}, {y[T]}, {t, 0, T}]

It appears that y[t] not to be a matrix the same size as A:

Part::take: Cannot take positions 5 through 8 in Norm[Symbol[]]. >>

The question is: how can I send the actual solution y[t] in a moment t in function g[t]?

(Possibly this question duplicates NDSolve with vector function though I am not sure.)

Answer 1

In the absence of specific definitions for the matrices in the question, I simply made up my own. Here is a way you can make it work:

J = {{1, 2}, {2, 1}};

T = 1;

γi = {{1, 0}, {0, 0}};

g = Identity;

Clear[y, Y]

y = Function[t, ##] &[Array[Y[##][t] &, Dimensions[J]]];

eqn1 = Thread /@ 
   Thread[y'[t] + (J\[Transpose] + g[y[t]]\[Transpose]).y[t] + 
      y[t].(J\[Transpose] + g[y[t]]\[Transpose]) - 
      2*y[t].(J\[Transpose] + g[y[t]]\[Transpose]).y[t] == 0];

eqn2 = Thread /@ Thread[y[T] == γi];

A = NDSolveValue[Flatten[{eqn1, eqn2}], Flatten[y[t]], {t, 0, T}]

The first little issue was that your second argument in NDSolveValue was wrapped in an extra {...} which is why you saw the first dimension being 1.

The second issue is that although NDSolve and NDSolveValue try to identify the dimensions of the unknowns from the structure of the equations (including initial conditions), this fails in your second example because the unknowns are too deep inside other functions.

Therefore, it seems to be necessary to revert to methods that make the structure of the equations more explicit. This is how it used to be in version 8 anyway, because the ability to recognize y as vectors or matrices is a relatively new feature, which apparently isn't always foolproof.

The method I chose here is to define the differential equation as eqn1, and the initial condition as eqn2, but before putting them into NDSolveValue I make sure that y in these equations is replaced by a function of t that explicitly has the same dimension as J. The components are labeled by a dummy variable Y[i,j]. The argument t is not included in the definition of y because otherwise the derivatives aren't done correctly. The two Threads are there to convert a matrix equation into a matrix of equations. Furthermore, I then flatten all the equations inside NDSolveValue so that they are correctly counted.