# How to use DSolve with vectors without decomposing into vector components?

Given a basic, two-dimensional ballistic trajectory problem, I can solve the equations of motion using DSolve (or NDSolve) by decomposing the vector equations of motion into equations with the scalar components, as follows:

DSolve[
Flatten@{
Thread[{rx''[t], ry''[t]} == {rx[t], -9.8 ry[t]}],
Thread[{vx[0], vy[0]} == {40 Cos[30 °],
40 Sin[30 °]}]}
},
{rx[t], ry[t], vx[t], vy[t]},
t
];


But which general approaches are possible to do this just with "vectors"? (I'm not interested in having anyone solve this particular example per se.)

For example, something like this would be desirable:

DSolve[
{
Thread[m r''[t] == {0, -9.8} r[t]],
Thread[v[0] == {40 Cos[30 °], 40 Sin[30 °]}]}
},
{r[t], v[t]},
t
]


Surely this must be possible in some way? I understand that "vectors" are just Lists, and that Mathematica doesn't know that I want r[t], r'[t] or v[t] to be "vectors". However, using {r[t][[1]],r[t][[2]]} hasn't worked for me yet, either. I would also prefer to implement the method myself "manually" rather than rely on packages which I would have to see as a black box.

• I've never figured out how to do this with DSolve. One can write a function to replace vector functions with component/coordinate functions and then reassemble the component solutions into a vector solution. I don't know what's the matter with such an approach. You'd think users would have asked for it and WRI considered it. Mar 22, 2015 at 13:38

The idea (as in my comment) is to replace a vector-valued dependent variable by coordinate functions K[i] for i = 1, 2, ...; then call DSolve; and finally reassemble the solution coordinates into a vector-valued solutions.

ClearAll[reassembleSols, vectorize];
reassembleSols[sols_, fns_, var_, dim_] := Function[sol0, MapThread[
Function[{f, s},
If[FreeQ[f, var],
f -> Function @@ {var, s},
f -> s]],
{fns, Partition[sol0, dim]}]
] /@ (Table[K[i][var], {i, 1 + Length[fns]*dim}] /. sols);

SetAttributes[vectorize, HoldAll];
vectorize[dim_, DSolve[ode_, f_, var_, rest___]] :=
With[{fns = Flatten[{f}]},
Block @@
With[{depvars = InternalProcessEquationsFindDependentVariables[ode, var]},
Hold[
depvars,
Module[{rules, sols},
Do[(#[\[FormalT]_] = Table[K[i][\[FormalT]], {i, 1 + j*dim - dim, j*dim}]) &@
depvars[[j]], {j, Length[depvars]}];
rules =
Table[depvars[[j]] -> Table[K[i], {i, 1 + j*dim - dim, j*dim}],
{j, Length[depvars]}];
sols = DSolve[ode, Flatten[fns /. rules], var, rest];
Clear @@ depvars;
If[FreeQ[sols, DSolve],
reassembleSols[sols, fns, var, dim],
sols]
]]]


OP's example.

vectorize[2,
DSolve[{m r''[t] == {{0, 0}, {0, -10}}. r[t],
v[t] == r'[t], {r[0] == {0, 0},
v[0] == {40 Cos[30 °], 40 Sin[30 °]}}}, {r[t],
v[t]}, t]
]
(*
{{r[t] -> {20 Sqrt[3] t,
-I Sqrt[10] *
E^(-((I Sqrt[10] t)/Sqrt[m])) (-1 + E^((2 I Sqrt[10] t)/Sqrt[m])) Sqrt[m]},
v[t] -> {20 Sqrt[3],
10 E^(-((I Sqrt[10] t)/Sqrt[m])) (1 + E^((2 I Sqrt[10] t)/Sqrt[m]))}}}
*)


Note: The variables ({r, v} in the OP's example) are defined to be functions of the form

r[t_] = {K[1][t], K[2][t]}


inside the Block. These definitions cause the ode to evaluate to a componentwise system of equations when DSolve is called. When defined as functions this way, the derivatives of r are automatically computed properly in the system passed to DSolve. If we were to apply rules, this would not happen. Afterwards, these definitions are cleared in order to reassemble the solutions for r[t] etc. so that r[t] stays r[t] and does not evaluate to a vector of K's.

• I should have pointed out that in the OP's ODE, there was {0, -9.8} r[t], which I changed to a Dot product. The problem is that Times threads prematurely, before r[t] is replaced by the coordinates. Dot, however, won't evaluate until r[t] becomes a vector, which is the desired behavior. Mar 22, 2015 at 15:18