Up One Level Sage Tutorial Contents Index
Previous: Up: 2.10 Calculus Next:


2.10.3 Systems of DEs using Laplace transforms

In this section, we provide a few details which are useful to teaching a lower level ordinary differential equations course using Sage.

The displacement from equilibrium (respectively) for a coupled spring attached to a wall on the left 

|------\/\/\/\/\---|mass1|----\/\/\/\/\/----|mass2|
         spring1               spring2
is modeled by the system of 2nd order ODEs

$\displaystyle m_1x_1'' + (k_1+k_2)x_1 - k_2x_2 = 0, \ m_2x_2''+ k_2(x_2-x_1) = 0, $

where $ x_1$ denotes the displacement from equilibrium of mass 1, denoted $ m_1$ , $ x_2$ denotes the displacement from equilibrium of mass 2, denoted $ m_2$ , and $ k_1$ , $ k_2$ are the respective spring constants.

Example: Use Sage to solve the above problem with $ m_1 = 2$ , $ m_2 = 1$ , $ k_1 = 4$ , $ k_2 = 2$ , $ x_1(0) = 3$ , $ x_1'(0) = 0$ , $ x_2(0) = 3$ , $ x_2'(0) = 0$ .

Solution: Take the Laplace transform of the first DE (for simplicity of notation, let $ x = x_1$ , $ y = x_2$ ): 

sage: de1 = maxima("2*diff(x(t),t, 2) + 6*x(t) - 2*y(t)")
sage: lde1 = de1.laplace("t","s"); lde1
2*(-?%at('diff(x(t),t,1),t=0)+s^2*?%laplace(x(t),t,s)-x(0)*s)-2*?%laplace(y(t),t,s)+6*?%laplace(x(t),t,s)
This says $ -2x'(0) + 2s^2*X(s) - 2sx(0) - 2Y(s) + 2X(s) = 0$ (where the Laplace transform of a lower case function is the upper case function). Take the Laplace transform of the second DE: 

sage: de2 = maxima("diff(y(t),t, 2) + 2*y(t) - 2*x(t)")
sage: lde2 = de2.laplace("t","s"); lde2
-?%at('diff(y(t),t,1),t=0)+s^2*?%laplace(y(t),t,s)+2*?%laplace(y(t),t,s)-2*?%laplace(x(t),t,s)-y(0)*s
This says $ s^2Y(s) + 2Y(s) - 2X(s) - 3s = 0$ . Solve these two equations: 

sage: var('s X Y')
(s, X, Y)
sage: eqns = [(2*s^2+6)*X-2*Y == 6*s, -2*X +(s^2+2)*Y == 3*s] 
sage: solve(eqns, X,Y)
[[X == (3*s^3 + 9*s)/(s^4 + 5*s^2 + 4), 
  Y == (3*s^3 + 15*s)/(s^4 + 5*s^2 + 4)]]
This says $ X(s) = (3s^3 + 9s)/(s^4 + 5s^2 + 4)$ , $ Y(s) = (3s^3 + 15s)/(s^4 + 5s^2 + 4)$ . Take inverse Laplace transforms to get the answer: 

sage: var('s t')
(s, t)
sage: inverse_laplace((3*s^3 + 9*s)/(s^4 + 5*s^2 + 4),s,t)
cos(2*t) + 2*cos(t)
sage: inverse_laplace((3*s^3 + 15*s)/(s^4 + 5*s^2 + 4),s,t)
4*cos(t) - cos(2*t)
Therefore, $ x_1(t) = \cos(2t) + 2\cos(t)$ , $ x_2(t) = 4\cos(t) - \cos(2t)$ . Using matplotlib, this can be plotted parametrically using (do not type the ...) 

sage: t = var('t')
sage: P = parametric_plot( (cos(2*t) + 2*cos(t), 4*cos(t) - cos(2*t) ),\
...   0, 2*pi, rgbcolor=hue(0.6) )
(type show(P) to view this). The individual components can be plotted using 

sage: t = var('t')
sage: Px = plot( cos(2*t) + 2*cos(t) , 0, 2*pi, rgbcolor=hue(0.6) )
sage: Py = plot( 4*cos(t) - cos(2*t) , 0, 2*pi, rgbcolor=hue(0.6) )

For an interactive ``zoomable'' plot, try (do not type the ...) 

sage: maxima.plot2d('cos(2*x) + 2*cos(x)','[x,0,1]',\
...   '[plot_format,openmath]') # not tested
if you have tcl/tk installed.

REFERENCES: Nagle, Saff, Snider, Fundamentals of DEs, 6th ed, Addison-Wesley, 2004. (see §5.5).

Up One Level Sage Tutorial Contents Index
Previous: Up: 2.10 Calculus Next:
Release 2008.05.04, documentation updated on May 4, 2008.
See About this document... for information on suggesting changes.