Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

System of linear second order differential equations

  1. Oct 10, 2009 #1
    I think I know how to solve
    [tex]\frac{d\vec{x}}{dt}= A \vec{x}[/tex]

    where A is a constant nXn matrix. We just compute the eigenvalues and the corresponding eigenvectors.

    But how do we solve
    [tex]\frac{d^2\vec{x}}{dt^2}= A \vec{x}[/tex]

    Can we say straight away that the solution is (following that of one dependent variable)
    [tex]\vec{x}(t) = exp(-Mt) \vec{x}_1+ exp(Mt) \vec{x}_2 [/tex]

    where M2=A and x1 and x2 are constant vectors.
  2. jcsd
  3. Oct 10, 2009 #2


    User Avatar
    Science Advisor
    Homework Helper

    Hi matematikawan! :smile:

    Try putting y = dx/dt + Mx :wink:
  4. Oct 10, 2009 #3
    This is great. How did you came out with the substitution y = dx/dt + Mx tiny-tim ?. I have try searching the internet for system of DE. Most of them are about first order DE.
    With your substitution I obtained

    Using http://en.wikipedia.org/wiki/Matrix_exponential" [Broken]
    [tex]\vec{y} = exp(Mt) \vec{c_1}[/tex]

    [tex]\frac{d\vec{x}}{dt} = -M\vec{x} + exp(Mt) \vec{c_1}[/tex]

    [tex] \vec{x}=exp(-Mt)\vec{c_2} + exp(-Mt)\int exp(2Mt) \vec{c_1}[/tex]

    [tex] \vec{x}=exp(-Mt)\vec{c_2} + exp(-Mt) (2M)^{-1} exp(2Mt) \vec{c_1}[/tex]

    So not exactly

    [tex]\vec{x}(t) = exp(-Mt) \vec{x}_1+ exp(Mt) \vec{x}_2 [/tex]

    as I predicted.
    Last edited by a moderator: May 4, 2017
  5. Oct 10, 2009 #4
    Or it could possible that

    [tex] \vec{x}=exp(-Mt)\vec{c_2} + exp(-Mt) exp(2Mt) (2M)^{-1}\vec{c_1}[/tex]

    Then it is of the form

    [tex]\vec{x}(t) = exp(-Mt) \vec{x}_1+ exp(Mt) \vec{x}_2 [/tex]

  6. Oct 10, 2009 #5


    User Avatar
    Science Advisor
    Homework Helper

    Hi matematikawan! :smile:

    (btw, we usually write the constant before the exponential, and we usually use A and B if we're expecting more than one :wink:)
    That's right …

    now write it dx/dt + Mx = BeMt,

    and use the particular solution and general solution method (do you know what that is?) …

    the g s gives you Ce-Mt, and the p s of course is (B/M - 1)eMt, or just BeMt. :wink:

    And I got it from the "D" method …

    we write D instead of d/dt, so the original equation is (D2 - A)x = 0,

    then (D + A)(D - A)x = 0,

    so either (D + A)x = 0 or (D - A)x = 0 …

    this works for any polynomial in D (it's called the characteristic polynomial of the original equation) …

    google "characteristic polynomial", or see the recurrence relation page in the PF Library for some details. :smile:
  7. Oct 11, 2009 #6
    That is a clever trick. I have seen in some old DE textbooks where they use similar 'D operator' technique to solve ODE. But that works only for DE with one dependent variable.

    In this case of system of DE

    [tex]\frac{d^2\vec{x}}{dt^2}= A \vec{x}[/tex]

    I'm quite lucky to be able to find M such that M^2=A. To find M in matlab we just use the command sqrtm(A).

    How about these system of DE
    [tex]\frac{d^3\vec{x}}{dt^3}= A \vec{x}[/tex]

    [tex]\frac{d^2\vec{x}}{dt^2} + \frac{d\vec{x}}{dt}= A \vec{x}[/tex]

    How are we going to 'factorise' the so called characteristic polynomial ?
  8. Oct 12, 2009 #7


    User Avatar
    Science Advisor
    Homework Helper

    (erm :redface: … it's not "so-called" … that really is its name!)

    For the first one, use the cube root of A (I don't know whether there's a matlab command for that: if not, is there one for lnA or expA?).

    The second is easy: the roots are -(1/2)I ±√((1/4)I + A2). :smile:
  9. Oct 12, 2009 #8
    Oh I'm sorry. :shy: I didn't mean to offend anyone. May be I'm not using proper words there. English is not my first language. So it is characteristic polynomial. Really sorry.

    But by the way when I Google 'characteristic polynomial'. This term gives the meaning as the polynomial whose roots are the eigenvalues of a given matrix A, i.e. det(A-aI)=0.

    In matlab we use the command expm(A) for matrix exponential. I'm not aware of commands for cube root or log of A. I'm trying to solve second order pde numerically. I discretise the variables x and t. Solution for the third order system of DE is not that important for me at the moment compare to that of second order.
    It looks like that we can factorise the characteristic quadratic by completing the square.

    [tex]D^2+D-A =\left(D + \frac{1}{2}I\right)^2 - \frac{1}{4} I - A [/tex]
    [tex]= \left(D + \frac{1}{2}I - \sqrt{\frac{1}{4}I+A} \right) \left(D + \frac{1}{2}I + \sqrt{\frac{1}{4}I+A} \right) [/tex]

    Thank you for your help. :smile:
  10. Oct 12, 2009 #9

    D H

    User Avatar
    Staff Emeritus
    Science Advisor

    Alternatively, form a state vector comprising [itex]\vec x[/itex] and [itex]d\vec x/dt[/itex].

    \vec u = \left [\vec x \\ d\vec x/dt \right]


    \dot{\vec u} =
    \bmatrix \boldsymbol 0_{N\times N} & \boldsymbol 1_{N\times N} \\
    \boldsymbol A & \boldsymbol 0_{N\times N} \endbmatrix \vec u

    where [itex]\boldsymbol 1[/itex] is the NxN identity matrix and [itex]\boldsymbol 0[/itex] is the NxN zero matrix. Denoting

    [tex] \boldsymbol B \equiv \bmatrix
    \boldsymbol 0_{N\times N} & \boldsymbol 1_{N\times N} \\
    \boldsymbol A & \boldsymbol 0_{N\times N} \endbmatrix[/tex]

    This is of course a constant matrix. Then

    [tex]\dot{\vec u} = \boldsymbol B \vec u[/tex]

    and tada! A linear first order equation with a constant matrix.
  11. Oct 13, 2009 #10
    Of course, we could then solve the equation with Runge-Kutta method, thus avoiding the need for matrix exponential. Although we have double the matrix size but their entries are sparse. Anyway I'm looking for numerical answer. Analytical answer may be suitable for theoretical purposes.

    Thanks D H, especially for formulating the equation with nice block matrices. Is there any chance of having B as
    [tex] \boldsymbol B \equiv \bmatrix
    \boldsymbol A & \boldsymbol 0_{N\times N} \\
    \boldsymbol 0_{N\times N} & \boldsymbol 1_{N\times N} \endbmatrix[/tex]

    I have try a few permutation but without success.
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook