# ParametricNDSolve: Fails to differentiate an matrix differential equation

It seems like there is a problem with differentiating a ParametricNDSolve object when the underlying system contains an matrix differential equation.

Create a simple ParametricNDSolveValue object which returns x[2], depending on parameter k.

pfun1 = ParametricNDSolveValue[{x'[t] == -k*x[t], x[0] == 10}, x[2], {t, 0, 10}, k]
ParametricFunction[SequenceForm["<", ">"]]


Now we have the same system but we add a matrix differential equation independent of the state x[t]. Here P[t] should be zero for all t due to the inital condition.

pfun2 = ParametricNDSolveValue[{x'[t] == -k*x[t], x[0] == 10,
P'[t] == P[t], P[0] == {{0, 0}, {0, 0}}}, x[2], {t, 0, 10}, k]
ParametricFunction[SequenceForm["<", ">"]]


We can take out the values for k=0.4.

pfun1[0.4]
pfun2[0.4]
4.49329
4.49329


Now we differentiate the systems and take out the value for k=0.4.

pfun1'[0.4]
pfun2'[0.4]
-8.98658
Thread::tdlen: Objects of unequal length in {{0.},{0.},{0.},{0.},{0.}}+{{0.},{-10.}} cannot be combined. >>*
ExperimentalNumericalFunction::nlnum1: The function value {{{0.},{-10.}}+{{0.},{0.},{0.},{0.},{0.}}} is not a list of numbers with dimensions {5} when the arguments are {0.,{0.,0.,0.,0.,10.},{0.4},{{0.},{0.},{0.},{0.},{0.}}}.
`

The derivative of the ParametricNDSolveValue object seems to fail when there is a matrix differential equation in the system.

Is this a bug?

-
pretty bad; I filed this as a bug. –  user21 Sep 11 '13 at 7:39