2
$\begingroup$

I am using NDSolve on a vector function $\mathbf y'(x)=f[\mathbf y(x)]$ with initial condition $\mathbf y_0$, where the dimension of the vector should be user-defined.

Clear[y, y0]
y0 = {1, 2}
sol = NDSolve[{y'[x] + y[x] == 0, y[0] == y0}, {y}, {x, 0, 10}]

I'm attempting to use the resulting InterpolatingFunction as an input into a matrix ODE, $\mathbf z'(x) = g[\mathbf z(x)]$, with initial condition $\mathbf z_0$.

Clear[z,z0]
z0 = {{1., 2.},{3.,4.}}
failing = NDSolve[{z'[x] + foo[y[x] /. sol, {1, 1}] z[x] == 0, z[0] == z0}, {z}, {x, 0, 10}]

Where

foo[vectx_?(VectorQ[#] &), vecty_?(VectorQ[#] &)] := vectx.vecty

The function foo[x,y] is failing to compute the dot product correctly because the dimensions don't match up (Length[y[x] /. sol,] is 1), although VectorQ[y[x] /. sol] yields True. Naively, I expected that the result of NDSolve would be something like

$\mathbf y(x) = (\text{InterpolatingFunction}[y_1(x)], \text{InterpolatingFunction}[y_2(x)],...)$

So that the dot product

$\mathbf y \cdot (a,b,...)=a \, \text{InterpolatingFunction}[y_1(x)] + b \, \text{InterpolatingFunction}[y_2(x)]+...$

Rather, it's something more like {InterpolatingFunction[{{0.,10.}},<>][t]}, which is trivially a vector since it's enclosed by {...}.

How do I express my "dot product" logic above in Mathematica-ese? Or is there a more elegant way to get the code flow that I want?

$\endgroup$
4
  • 1
    $\begingroup$ You can change sol = NDSolve... to sol = First@NDSolve.. to get y[x]/.sol to have the form InterpolatingFunction[...] instead of {InterpolatingFunction[...]} $\endgroup$
    – ssch
    Jun 24, 2013 at 10:33
  • $\begingroup$ The (VectorQ[#] &) is not needed, VectorQ will do the job in this case. $\endgroup$
    – mmal
    Jun 24, 2013 at 12:20
  • $\begingroup$ Alternatively (VectorQ[#, NumericQ]&) is more precise. (Note that VectorQ[y[x] /. sol] is True, but VectorQ[y[0.4] /. sol] is False.) One could also add a condition, foo[..] /; Length[vectx] == Length[vecty] := .., but that's not really necessary. The major problem has been pointed out by ssch and partial81. $\endgroup$
    – Michael E2
    Jun 24, 2013 at 14:11
  • $\begingroup$ @ssch, this worked perfectly, thank you. $\endgroup$
    – cosmoguy
    Jun 24, 2013 at 21:41

2 Answers 2

1
$\begingroup$

I think the following could help:

Use in your code

failing = NDSolve[{z'[x] + 
foo[Evaluate[y[x] /. sol][[1]], {1, 1}] z[x] == 0, z[0] == z0}, {z}, {x, 0, 10}]

Please tell me if this works for you.

Btw. You can also skip the Evaluate and use only foo[y[x] /. sol[[1]], {1, 1}]. The result is the same.

$\endgroup$
2
  • $\begingroup$ This works, cheers! $\endgroup$
    – cosmoguy
    Jun 24, 2013 at 21:43
  • $\begingroup$ Great to hear! Btw. First[expr] (as ssch writes) and expr[[1]], are equivalent. So, chose the notation you like more ;-) $\endgroup$
    – partial81
    Jun 25, 2013 at 6:19
0
$\begingroup$

You asked if there's a more elegant way to achieve your output. I think the simplest approach is to combine your 2 ODEs into a single ODE, something like:

y0 = {1, 2};
z0 = {{1, 2}, {3, 4}};

{ysol, zsol} = NDSolveValue[
    {
    y'[x] + y[x] == 0, y[0] == y0,
    z'[x] + y[x].{1,1} z[x] == 0, z[0] == z0
    },
    {y, z},
    {x, 0, 10}
];

Note that I used NDSolveValue instead of NDSolve, and just used Dot instead of your auxiliary function foo.

Visualization:

Plot[zsol[x], {x, 0, 10}]

enter image description here

$\endgroup$

Your Answer

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

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