# Solution of a system of ODE with Dsolve

I am trying to solve following ODE with Dsolve.

a1[t_] = 2 - 3*Exp[-2*t];
a2[t_] = 1 + Exp[-3*t];
a3[t_] = 4 + Exp[-5*t]


$$B=\begin{bmatrix}a1[t] & a2[t] & 0\\a2[t] & a1[t] - a2[t] & 0\\ a3[t] & 0 & 1 \end{bmatrix}$$

  X[t_] = {x[t], y[t], z[t]}
system = X'[t] == B.X[t]


which gives

{x'[t],y'[t],z'[t]} == {(2 - 3 E^(-2 t)) x[t] + (1 + 3 E^(-3 t)) y[t],
(1 +E^(-3 t)) x[t] + (1 -E^(-3 t) - 3 E^(-2 t)) y[t], (4+E^(-5 t))x[t]+z[t]}

sol = DSolve[system, {x, y, z}, t]


Mathematica gives output as an input. Can you give any recommendation?

DSolve[x'[t],y'[t],z'[t]} == {(2 - 3 E^(-2 t)) x[t] + (1 + E^(-3 t)) y[t],
(1 + E^(-3 t)) x[t] + (1 - E^(-3 t) - 3 E^(-2 t)) y[t],
(4 + E^(-5 t)) x[t] + z[t]}, {x, y, z}, t]

• Yours code give me: "DSolve::deqx: Supplied equations are not differential or integral equations of the given functions."? – Mariusz Iwaniuk Oct 12 at 17:57
• hi. thank you for reply.i edited the post to become more clear. – cabri61 Oct 12 at 18:38
• I do not think there is analytical solution to this. Try numerical. – Nasser Oct 12 at 18:40
• Maple cannot do it, either. – Mariusz Iwaniuk Oct 12 at 18:41
• @bbgodfrey.Yes Maple is better for solving differential equations. 12000.org/my_notes/kamek/mma_12_maple_2019/KEse1.htm#x3-20001 – Mariusz Iwaniuk Oct 13 at 13:01

With help, DSolve can obtain symbolic solutions for {x, y}, after which z can be represented in terms of an integral. It is convenient to start by labeling the three equations,

{eqx, eqy, eqz} = Thread[{x'[t], y'[t], z'[t]} ==
{(2 - 3 E^(-2 t)) x[t] + (1 + E^(-3 t)) y[t],
(1 + E^(-3 t)) x[t] + (1 - E^(-3 t) - 3 E^(-2 t)) y[t],
(4 + E^(-5 t)) x[t] + z[t]}];


Then, just as in the question,

DSolveValue[{eqx, eqy, eqz}, {x, y, z}, t]


returns unevaluated. Now, since z appears only in eqz, we might try

DSolveValue[{eqx, eqy}, {x, y}, t]


but it too returns unevaluated. Nonetheless, eliminating y from {eqx, eqy} leads to an equation that can be solved for x.

Solve[{eqx, eqy}, {y[t], y'[t]}] // Flatten // Simplify;
Simplify[D[eqx, t] /. %];
sx = DSolveValue[%, x, t]

(* Function[{t}, E^(1/6 E^(-3 t) (1 + Sqrt + 9 E^t)) (E^t)^(3/2 - Sqrt/2) C +
(E^(E^(-3 t)/6 - 1/6 Sqrt E^(-3 t) + (3 E^(-2 t))/2) (E^t)^(3/2 + Sqrt/2) C)
/Sqrt] *)


Then, y is obtained algebraically by

Solve[Simplify[eqx /. x -> sx], y[t]][[1, 1]] // Simplify // Values;
sy = Function[{t}, Evaluate[%]]
(* Function[{t}, -(1/10) E^(1/6 E^(-3 t) (1 - Sqrt + 9 E^t + 6 E^(3 t) t))
(E^t)^(1/2 - Sqrt/2)
(5 (1 + Sqrt) E^(1/3 Sqrt E^(-3 t)) C + (-5 + Sqrt) (E^t)^Sqrt C)] *)


To verify the correctness of these solutions, insert them into the original equations.

Simplify[{eqx, eqy} /. x -> sx /. y -> sy]
(* {True, True} *)


Finally, solve eqz for z.

DSolveValue[eqz /. x -> sx, z[t], t] // Flatten // Simplify;
sz = Function[{t}, Evaluate[%]]
(* Function[{t}, E^t (C + Inactive[Integrate][
1/5 E^(-(1/6) E^(-3 K) (-1 + Sqrt - 9 E^K + 30 E^(3 K) K))
(E^K)^(1/2 - Sqrt/2) (1 + 4 E^(5 K)) (5 E^(1/3 Sqrt E^(-3 K)) C +
Sqrt (E^K)^Sqrt C), {K, 1, t}])] *)