# Solving System of Differential Equations and Plotting

I'm new to Mathematica and I'm struggling with this exercise of Initial Value Problem:
Two equations with conditions

• $$x_1$$'=-2$$x_1$$-3$$x_2$$+cos(5t)
• $$x_2$$'=-$$x_1$$-5$$x_2$$+sin(10t)
• $$x_1$$(0)=1
• $$x_2$$(0)=-1

Then i need to build the graph of R(t)= √($$x_1$$^2+$$x_2$$^2)
For 0≤t≤10
So far I have this code, but the graph is appearing empty.

solution = NDSolve[
{x1'[t] == -2 x1[t] - 3 x2[t] + Cos[5 t],
x2'[t] == -x1[t] - 5 x2[t] + Sin[10 t],
x1[0] == 1, x2[0] == -1},
{x1, x2},
{t, 0, 10}
]
X = x1[t] /. solution[[1, 1]];
Y = x2[t] /. solution[[1, 2]];
z[t_] := Sqrt[X[t]^2 + Y[t]^2];
Plot[z[t], {t, 0, 10}]


Thanks for any help provided.

• Use Plot[Evaluate[Sqrt[x1[t]^2 + x2[t]^2] /. solution], {t, 0, 10}, PlotRange -> All] – Alex Trounev Apr 15 at 20:22
• X = x1 /. solution[[1, 1]];Y = x2 /. solution[[1, 2]]; will fix your problem. – m_goldberg Apr 15 at 21:32
• Welcome to Mathematica.SE! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Take the tour! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! – Michael E2 Apr 18 at 10:25

Starting at around $$t=10.25$$, the solution calculated with NDSolve[] starts to blow up:

 Plot[Evaluate[Sqrt[x1[t]^2 + x2[t]^2] /. solution], {t, 0, 10.5}, PlotRange -> All]


(Thank you @Alex Trounev for clarification on how to plot this with Evaluate)

 Plot[Evaluate[{x1[t], x2[t]} /. solution], {t, 0, 10.5}, PlotRange -> All]


 ParametricPlot[Evaluate[{x1[t], x2[t]} /. solution], {t, 0, 10.5},
PlotRange -> All, Axes -> False, Frame -> True,
FrameLabel -> {Subscript[x, 1], Subscript[x, 2]}]


Better to use DSolve[]:

solution =
DSolve[{x1'[t] == -2 x1[t] - 3 x2[t] + Cos[5 t],
x2'[t] == -x1[t] - 5 x2[t] + Sin[10 t], x1[0] == 1,
x2[0] == -1}, {x1, x2}, {t, 0, 10}]

Plot[Evaluate[Sqrt[x1[t]^2 + x2[t]^2] /. solution], {t, 0, 24.5}, PlotRange -> All]


 Plot[Evaluate[{x1[t], x2[t]} /. solution], {t, 0, 25}, PlotRange -> All]


 ParametricPlot[Evaluate[{x1[t], x2[t]} /. solution], {t, 0, 25},
PlotRange -> All, Axes -> False, Frame -> True,
FrameLabel -> {Subscript[x, 1], Subscript[x, 2]}]


• If you put the same upper limit of t in both the DSolve and the Plot lines of code (for example 12), you will not see a steep increase after 10.5...If these limits are not the same I doubt your solution is mathematically.correct. – Sigis K Apr 18 at 10:55
• This is true for NDSolve, but I don't think for DSolve. Since DSolve finds the exact solution, it does not depend on the range of $t$ (unless, there are multiple solutions in different regions). Nevertheless, you are right, if you ask for NDSolve in a certain range, that's all you get. Going further in $t$ it will extrapolate, causing the large errors. – mjw Apr 18 at 15:48
Clear["Global*"]

eqns = {x1'[t] == -2 x1[t] - 3 x2[t] + Cos[5 t],
x2'[t] == -x1[t] - 5 x2[t] + Sin[10 t], x1[0] == 1, x2[0] == -1};

solution = NDSolve[eqns, {x1, x2}, {t, 0, 10}];
X = x1[t] /. solution[[1, 1]];
Y = x2[t] /. solution[[1, 2]];


As defined, X and Y are implicit functions of t they do not take t as an explicit argument

z[t_] := Sqrt[X^2 + Y^2];

Plot[z[t], {t, 0, 10}, PlotRange -> All]


DSolve will provide an exact solution

solExact = DSolve[eqns, {x1, x2}, t][[1]];


Verifying the solution,

eqns /. solExact // Simplify

(* {True, True, True, True} *)

\$Assumptions = t >= 0;

z[t_] = Norm[{x1[t], x2[t]} /. solExact // Simplify] // Simplify

(* (1/41974802)E^(-Sqrt[
21] t) \[Sqrt]((19573853 E^(1/2 (-7 + Sqrt[21]) t) +
938127 Sqrt[21] E^(1/2 (-7 + Sqrt[21]) t) +
19573853 E^(1/2 (-7 + 3 Sqrt[21]) t) -
938127 Sqrt[21] E^(1/2 (-7 + 3 Sqrt[21]) t) -
487764 E^(Sqrt[21] t) Cos[5 t] + 3314860 E^(Sqrt[21] t) Cos[10 t] +
948430 E^(Sqrt[21] t) Sin[5 t] - 1592372 E^(Sqrt[21] t) Sin[10 t])^2 +
4 (9755223 E^(1/2 (-7 + Sqrt[21]) t) -
4189868 Sqrt[21] E^(1/2 (-7 + Sqrt[21]) t) +
9755223 E^(1/2 (-7 + 3 Sqrt[21]) t) +
4189868 Sqrt[21] E^(1/2 (-7 + 3 Sqrt[21]) t) +
1151665 E^(Sqrt[21] t) Cos[5 t] + 325290 E^(Sqrt[21] t) Cos[10 t] +
3590485 E^(Sqrt[21] t) Sin[5 t] + 432171 E^(Sqrt[21] t) Sin[10 t])^2) *)

Plot[z[t], {t, 0, 10}, PlotRange -> All]


Here's an alternative approach using NDSolveValue:

{X, Y, z[t_]} =
NDSolveValue[{x1'[t] == -2 x1[t] - 3 x2[t] + Cos[5 t],
x2'[t] == -x1[t] - 5 x2[t] + Sin[10 t], x1[0] == 1, x2[0] == -1},
{x1, x2, Sqrt[x1[t]^2 + x2[t]^2]}, {t, 0, 10}];
Plot[z[t], {t, 0, 10}, PlotRange -> All]
`