# Working with a system of differential equations that cannot be solved explicitly

I have to work a lot with three functions $\;o_1(t), o_2(t), o_3(t)\;$ that are solutions to a certain system of differential equations:

Halp = {
D[o1[t], t] == o1[t]*o2[t] + o1[t]*o3[t] - o2[t]*o3[t],
D[o2[t], t] == o1[t]*o2[t] + o2[t]*o3[t] - o1[t]*o3[t],
D[o3[t], t] == o1[t]*o3[t] + o2[t]*o3[t] - o1[t]*o2[t]}
sol = DSolve[Halp, {o1[t], o2[t], o3[t]}, {t}]


The system cannot be solved explicitly with Mathematica, but I do not need the solution. What I would like to do is to compute higher order derivatives assuming that $o_1(t)$, $o_2(t)$, and $o_3(t)$ are solutions to the above system.

How could I tell Mathematica to do this?

There are many ways which one could exploit. For example let's use calculate o1''[t] eliminating o1'[t], o2'[t], o3'[t]. We can use Eliminate:

#1/.( Eliminate[ Join[ D[Halp, t], Halp], {#2, #3}] //ToRules)&[o1''[t], o2''[t], o3''[t]]

2 o1[t] o2[t]^2 + 2 o1[t] o2[t] o3[t] - 2 o2[t]^2 o3[t] + 2 o1[t] o3[t]^2 - 2 o2[t] o3[t]^2


Changing the order in the square bracket we can calculate o2''[t] and o3''[t].

For the third order derivative of e.g. o1[t] we can eliminate other dependent variables in many different ways, let's point out one of them:

#1 /.(
Eliminate[ Join[ D[Halp, {t, 2}], D[Halp, t], Halp], {##2}]
//ToRules) &[o1'''[t], o2'''[t], o3'''[t], o2''[t], o3''[t], o1''[t]]

  2 o1[t]^2 o2[t]^2 + 4 o1[t] o2[t]^3 - 4 o1[t]^2 o2[t] o3[t] + 8 o1[t] o2[t]^2 o3[t]
- 4 o2[t]^3 o3[t] + 2 o1[t]^2 o3[t]^2 + 8 o1[t] o2[t] o3[t]^2 - 10 o2[t]^2 o3[t]^2
+ 4 o1[t] o3[t]^3 - 4 o2[t] o3[t]^3


In the above we eliminated all dependent variables starting from the second position in the square bracket (note very useful sign ##n i.e. SlotSequence ) and calculated o1'''[t]. Of course you can try to get rid of different dependent variables. In general this is a difficult issue to decide which dependent variables could be eliminated and how they could be represented by another variables. I suggest to take a closer look at GroebnerBasis with the MonomialOrder -> EliminationOrder option.

GroebnerBasis[ polys, vars, elims, MonomialOrder -> EliminationOrder]


If given equations are relatively simple you can't observe a common problem in differential elimination - so called intermediate expresssion swell. For interesting mathematical issues related to the problem I'd suggest to study e.g. a recent monograph Involution by Werner Seiler (Springer 2009). However I couldn't mention interesting packages in Mathematica related to differential elimination (there are such packages e.g. in Maple).

I would simply define rules for derivatives of o1, o2 and o3 and avoid the use of Eliminate or GroebnerBasis:

o1' = Function[t, o1[t]*o2[t]+o1[t]*o3[t]-o2[t]*o3[t]];
o2' = Function[t, o1[t]*o2[t]+o2[t]*o3[t]-o1[t]*o3[t]];
o3' = Function[t, o1[t]*o3[t]+o2[t]*o3[t]-o1[t]*o2[t]];

Derivative[n_?Positive][o1] := simplifyFunction @ Derivative[n-1][o1']
Derivative[n_?Positive][o2] := simplifyFunction @ Derivative[n-1][o2']
Derivative[n_?Positive][o3] := simplifyFunction @ Derivative[n-1][o3']

simplifyFunction[HoldPattern@Function[t_, body_]] := Function[t, Evaluate @ Simplify @ body]


Example:

o1'''[t]
o2'''[t]
o3'''[t]


2 o1[t]^2 (o2[t] - o3[t])^2 - 2 o2[t] o3[t] (2 o2[t]^2 + 5 o2[t] o3[t] + 2 o3[t]^2) + 4 o1[t] (o2[t]^3 + 2 o2[t]^2 o3[t] + 2 o2[t] o3[t]^2 + o3[t]^3)

2 (2 o1[t]^3 (o2[t] - o3[t]) - 2 o1[t] (o2[t] - o3[t])^2 o3[t] + o2[t] o3[t]^2 (o2[t] + 2 o3[t]) + o1[t]^2 (o2[t]^2 + 4 o2[t] o3[t] - 5 o3[t]^2))

2 (-2 o1[t]^3 (o2[t] - o3[t]) - 2 o1[t] o2[t] (o2[t] - o3[t])^2 + o2[t]^2 o3[t] (2 o2[t] + o3[t]) + o1[t]^2 (-5 o2[t]^2 + 4 o2[t] o3[t] + o3[t]^2))

• Would you recommend avoiding (if possible) the use of GroebnerBasis and Eliminate in general (e.g. because of their efficiency), or is it just a neat shortcut which sometimes appear to be useful? Nonetheless we cannot completely get rid of GroebnerBasis in differential elimination. +1. – Artes Apr 20 '18 at 9:14
• @Artes I'm just saying that for a system of equations of the form y'[t] == f[y[t]], where both y and f are vectors, there is no need to use Eliminate. – Carl Woll Apr 20 '18 at 14:51
• Ok, fair enough. – Artes Apr 20 '18 at 14:56