Preface

This is a tutorial made solely for the purpose of education and it was designed for students taking Applied Math 0330. It is primarily for students who have very little experience or have never used Mathematica before and would like to learn more of the basics for this computer algebra system. As a friendly reminder, don't forget to clear variables in use and/or the kernel.

Finally, the commands in this tutorial are all written in bold black font, while Mathematica output is in normal font. This means that you can copy and paste all commands into Mathematica, change the parameters and run them. You, as the user, are free to use the scripts for your needs to learn the Mathematica program, and have the right to distribute this tutorial and refer to this tutorial as long as this tutorial is accredited appropriately.

Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to Mathematica tutorial for the second course APMA0330
Return to Mathematica tutorial for the first course APMA0340
Return to the main page for the course APMA0340
Return to the main page for the course APMA0330
Return to Part IV of the course APMA0330

Complex Roots

An exponential solution \( y = C\, e^{\lambda\,t} , \) where C is an arbitrary real number and λ is a complex or real number, to the homogeneous constant coefficient linear differential equation

\[ a_n y^{(n)} + a_{n-1} y^{(n-1)} + \cdots + a_1 y' + a_0 \,y = 0 , \qquad a_n \ne 0, \]
is called a modal solution and \( C\, e^{\lambda\,t} \) is called a mode of the system. If we assign to the left-hand side a linear differential operator
\[ L\left[ \texttt{D} \right] = a_n \texttt{D}^n + a_{n-1} \texttt{D}^{n-1} + \cdots + a_1 \,\texttt{D} + a_0 , \qquad \texttt{D} = \frac{{\text d}}{{\text d}t} , \]
then the exponential function/mode is a solution to the homogeneous linear equation \( L\left[ \texttt{D} \right] y = 0 \) exactly when λ is a root of the characteristic equation
\[ L\left( \lambda \right) = a_n \lambda^n + a_{n-1} \lambda^{n-1} + \cdots + a_1 \,\lambda + a_0 =0 . \]
If the above polynomial equation has a simple complex root \( \lambda = \alpha + {\bf j} \beta , \) then its complex conjugate \( \overline{\lambda} = \alpha - {\bf j} \beta \) also must be a root because we consider only real-valued differential equations and all coefficients \( a_0 , a_1 , \cdots , a_n \) are assumed to be real numbers. In this case, the characteristic polynomial can be factored:
\[ L\left( \lambda \right) = a_n \left[ \left( \lambda - \alpha \right)^2 + \beta^2 \right] M \left( \lambda \right) , \qquad M \left( \alpha \pm {\bf j} \beta \right) \ne 0 . \]
In the above equation, \( L\left( \lambda \right) \) is the real-valued polynomial of degree n and \( M\left( \lambda \right) \) is the real-valued polynomial of degree n-2. If the characteristic polynomial has complex pair roots \( \lambda = \alpha \pm {\bf j} \beta \) of multiplicity m, then it can be factored:
\[ L\left( \lambda \right) = a_n \left[ \left( \lambda - \alpha \right)^2 + \beta^2 \right]^m M \left( \lambda \right) , \qquad M \left( \alpha \pm {\bf j} \beta \right) \ne 0 , \]
where \( M\left( \lambda \right) \) is the real-valued polynomial of degree n-2m.

Now we go from the general case to analyze solutions of the second order linear differential equation with constant coefficients:

\[ a \, \ddot{y} + b\, \dot{y} + c \,y =0 , \qquad a\ne 0, \quad \ddot{y} = \frac{{\text d}^2 y}{{\text d} t^2} , \quad \dot{y} = \frac{{\text d} y}{{\text d} t} . \]
Assume that the characteristic polynomial \( L\left( \lambda \right) = a\,\lambda^2 + b\,\lambda + c \) has a pair of complex roots:
\[ \lambda = \frac{-b \pm {\bf j} \sqrt{4ac-b^2}}{2a} , \qquad \mbox{where} \quad 4ac-b^2 >0 . \]
Recall that j is a unit vector in positive vertical direction on the complex plane so that \( {\bf j}^2 =-1 . \) Upon introducing two real numbers
\[ \alpha = -\frac{b}{2a} , \qquad \beta = \frac{\sqrt{4ac-b^2}}{2a} , \qquad \mbox{where} \quad 4ac-b^2 >0 , \]
we see that the characteristic equation \( a\,\lambda^2 + b\,\lambda + c =0 \) has two complex conjugate roots \( \lambda = \alpha \pm {\bf j} \beta . \) So the differential equation \( L\left[ \texttt{D} \right] y = 0 \) has two exponential solutions
\[ z_1 (t) = e^{\left( \alpha + {\bf j} \beta \right) t} = e^{\alpha\,t} \, e^{{\bf j} \beta\,t} , \qquad z_2 (t) = \overline{z_1 (t)} = e^{\left( \alpha - {\bf j} \beta \right) t} = e^{\alpha\,t} \, e^{-{\bf j} \beta\,t} . \]
We use the letter z here to indicate the functions are complex valued. The general solution is a linear combination of these two basic solutions:
\[ y(t) = C_1 z_1 (t) +C_2 z_2 (t) = C_1 e^{\alpha\,t} \, e^{{\bf j} \beta\,t} + C_2 e^{\alpha\,t} \, e^{-{\bf j} \beta\,t} = e^{\alpha\,t} \left[ C_1 e^{{\bf j} \beta\,t} + C_2 e^{-{\bf j} \beta\,t} \right] . \]
But, because the differential equation \( L\left[ \texttt{D} \right] y = 0 \) has real coefficients, we were expecting real-valued solutions. This could be achieved only when one arbitrary constant \( C_2 = \overline{C_1} \) is complex conjugate of another one. So upon introducing two real arbitrary constants A and B such that \( 2\,C_1 = A + {\bf j} B , \quad 2\,C_2 = A - {\bf j} B , \) we get
\[ y(t) = e^{\alpha\,t} \left[ \frac{1}{2} \left( A + {\bf j} B \right) e^{{\bf j} \beta\,t} + \frac{1}{2} \left( A - {\bf j} B \right) e^{-{\bf j} \beta\,t} \right] = e^{\alpha\,t} \, \frac{A}{2} \left[ e^{{\bf j} \beta\,t} + e^{-{\bf j} \beta\,t} \right] + e^{\alpha\,t} \, \frac{B}{2} \left[ {\bf j}\,e^{{\bf j} \beta\,t} -{\bf j}\, e^{-{\bf j} \beta\,t} \right] . \]
Now we recall Euler's formula
\[ e^{{\bf j}\theta} = \cos \theta + {\bf j}\,\sin \theta, \qquad \Longrightarrow \qquad \cos\theta = \frac{1}{2} \left[ e^{{\bf j}\theta} + e^{-{\bf j}\theta} \right] , \quad \sin\theta = \frac{1}{2{\bf j}} \left[ e^{{\bf j}\theta} + e^{-{\bf j}\theta} \right] . \]
This allows us to rewrite the general solution as
\[ y(t) = e^{\alpha\,t} \, A\, \cos \left( \beta\,t \right) - e^{\alpha\,t} \, B\, \sin \left( \beta\,t \right) \]
that involve only real-valued functions. Here A and B are arbitrary real constants and the real number β is referred to as pseudo-frequency.

Theorem: If \( z(t) \) is a complex-valued solution to

\[ L\left[ \texttt{D} \right] z = a_n \texttt{D}^n z + a_{n-1} \texttt{D}^{n-1} z + \cdots + a_1 \,\texttt{D} \,z+ a_0 \, z =0 , \qquad \texttt{D} = \frac{{\text d}}{{\text d}t} , \]
where all coefficients \( a_0 , a_1 , \ldots , a_n \) are real numbers, then the real and imaginary parts of z are also solutions. ■

Example: Consider second order differential equation

\[ y'' + a^2 y =0 , \qquad a>0. \]
The corresponding characteristic equation, \( \lambda^2 + a^2 =0 \) has two complex conjugate roots \( \lambda = \pm a{\bf j} , \) where j is the unit vector on the complex plane in positive vertical direction. We use Mathematica to show that the fundamental set of solutions consists of sine and cosine functions.
Clear[y]
y[t_]=c1 Cos[a*t]+c2*Sin[a*t]
Simplify[y''[t]+a^2*y[t]==0]

In this syntax, you define cosine and sine functions by Cos[ ] and Sin[ ]. When you enter the general solution, you can verify that this is indeed the solution.

Clear[y]
y[t_]=c1 Cos[a*t]+c2*Sin[a*t]

Out[7]= c1 Cos[a t] + c2 Sin[a t]
Simplify[y''[t]+a^2*y[t]==0]
Out[8]= True
Now we plot family of solutions for \( a =1 . \)
DSolve[y''[x] + y[x] == 0, y[x], x]
Out[14]= {{y[x] -> C[1] Cos[x] + C[2] Sin[x]}}
toplot1 =
Table[C[2] Cos[x] - C[1] Sin[x] /. {C[2] -> i, C[1] -> 0}, {i, -2.5, 2.5, 0.5}];
toplot2 =
Table[C[2] Cos[x] - C[1] Sin[x] /. {C[2] -> 0, C[1] -> i}, {i, -2.5, 2.5, 0.5}];
plot1 = Plot[Evaluate[toplot1], {x, 0, 4*Pi},
DisplayFunction -> Identity]

plot2 = Plot[Evaluate[toplot2], {x, 0, 4*Pi}, DisplayFunction -> Identity]

Example: We verify that the function \( y(x) = C_1 \,\cos 2x + C_2 \,\sin 2x + x^3 -3x/2 \) is a solution of the differential equation \( y'' + 4\,y = 4\,x^3 .\)

y[x_] = c1 Cos[2 x] + c2 Sin[2 x] + x^3 - 3*x/2;
Simplify[y''[x] + 4 y[x] == 4 x^3]
Out[2]= True
ics = {y[0] == 1, y'[0] == -1};
c1c2 = Solve[ics, {c1, c2}];
SolRule[x_] = y[x] /. c1c2[[1]]
SolRule''[x] + 4 SolRule[x] == 4 x^3 && SolRule[0] == 1 && SolRule'[0] == -1
Simplify[%]
Out[7]= True

To check linearly dependence:

Clear[y];
y={Cos[2 x], Sin[2 x]}
w={y,D[y,x]}
Det[w]
Out[11]= 2 Cos[2 x]^2 + 2 Sin[2 x]^2
y={x,3 x}
w={y,D[y,x]}
Det[w]
Out[15]= 0
ExpToTrig[Exp[I 2 x]]
Out[9]= Cos[2 x] + I Sin[2 x]
y={E^(a x) Cos[2x], E^(a x) Sin[2x]};
w = {y, D[y, x]};
Det[w] (* or *)
Simplify[Det[w], Trig -> True]
Out[22]= 2 E^(2 a x)
which is not zero. So these two functions are linearly independent.

Example: Consider the differential equation

\[ y'' -4\,y' + 5\, y =0 . \]
Its general solution is
\[ y= e^{2x} \left( C_1 \cos x + C_2 \sin x \right) . \]
DSolve[y''[x]-4 y'[x]+5 y[x]==0,y[x],x]
Out[1]= {{y[x] -> E^(2 x) C[2] Cos[x] + E^(2 x) C[1] Sin[x]}}
soln[x_] = Expand[y[x]/.%[[1]]]
Out[2]= E^(2 x) C[2] Cos[x] + E^(2 x) C[1] Sin[x]
y1[x_] = Coefficient[soln[x],C[1]]
Out[3]= E^(2 x) Sin[x]
y2[x_] = Coefficient[soln[x], C[2]]
Out[4]= E^(2 x) Cos[x]
Wronskian[{y1[x], y2[x]}, x] (* built-in Wronskian *)
Out[5]= -E^(4 x)
DSolve[{y''[x]-4 y'[x]+5 y[x]==0,y[0]==1},y[x],x]
Out[6]= {y[x] -> E^(2 x) (Cos[x] + C[1] Sin[x])}}

IVP:

DSolve[{y''[x] - 4 y'[x] + 5 y[x] == 0, y[0] == 1}, y[x], x];
soln[x_] = Expand[y[x] /. %[[1]]] /. C[1] -> c1
curves = Table[soln[x], {c1, -2, 2}];
Plot[Evaluate[curves], {x, -2, 2}, PlotRange -> {-10, 10}]
Out[11]= E^(2 x) Cos[x] + c1 E^(2 x) Sin[x]
Out[13]=

DSolve[{y''[x] - 4 y'[x] + 5 y[x] == 0, y'[0] == 1}, y[x], x];
soln[x_] = Expand[y[x] /. %[[1]]] /. C[1] -> c2
curves = Table[soln[x], {c2, -2, 2}];
Plot[Evaluate[curves], {x, -0.5, 2}, PlotRange -> {-10, 10}]
Out[15]= 1/2 E^(2 x) Cos[x] - 1/2 c2 E^(2 x) Cos[x] + c2 E^(2 x) Sin[x]
Out[17]=

 

Return to Mathematica page

Return to the main page (APMA0330)
Return to the Part 1 (Plotting)
Return to the Part 2 (First Order ODEs)
Return to the Part 3 (Numerical Methods)
Return to the Part 4 (Second and Higher Order ODEs)
Return to the Part 5 (Series and Recurrences)
Return to the Part 6 (Laplace Transform)
Return to the Part 7 (Boundary Value Problems)