In Harmonic Oscillator via Discrete Lagrangian (https://community.wolfram.com/groups/-/m/t/3161109), we numerically integrated the undamped harmonic oscillator via the discrete Lagrangian, comparing the numerical solution favorably to the analytical solution. The purpose of that exercise was to build confidence in discrete-Lagrangian methods. We’re interested in these methods at large because they promise far better conservation of momentum and energy and scale to very large systems. The damped harmonic oscillator (DHO) is much more interesting than the undamped, but not much more difficult. DHO is both elementary and fundamental. DHO appears everywhere in physics and engineering. One can spend a lifetime studying its parameter spaces and relating it to stability theory, Bode plots, Laplace transforms, and more. In this notebook, we extend Kobilarov and Marsden, Reference 3, to add a non-conservative damping force following the advice of References 1 and 2.

Analytical Solution

In[]:=
ClearAll[L,F,V,T,Ld];

Lagrangian

Potential energy, entailing a spring constant κ as a free variable:
In[]:=
V[q_]:=
1
2
κ
2
q
Lagrangian, entailing the kinetic energy
1
2
m
2
v
with mass
m
as a free variable:
In[]:=
T[v_]:=
1
2
m
2
v
;​​L[q_,v_,t_]:=T[v]-V[q]

Modified Euler-Lagrange Equation

The Lagrange-D'Alembert Principle adds a force term to the ordinary Euler-Lagrange Equation, here with a minus sign distributed across the equation versus the wiki page:

t
∂L
∂v
-
∂L
∂q
+
∂F
∂v
=0
(
1
)

First Term of the Modified Euler-Lagrange Equation

The first term is a total derivative of the partial derivative of the Lagrangian
L
with respect to its second parameter, generalized velocity
v
. The total derivative entails a derivative of mass with respect to time via the chain rule. We insert an explicit constant-mass rule. Non-constant mass would be useful for, say, the rocket equation, but is not pertinent here.
In[]:=
Dt[D[L[q,v,t],v],t]/.{Dt[m,t]->0}
Out[]=
mDt[v,t]

Second Term of the Euler-Lagrange Equations

No explicit assumption of constant spring-force coefficient
κ
is necessary because this term is just a partial derivative of
L
with respect to
q
and does not entail the chain rule.
In[]:=
D[L[q,v,t],q]
Out[]=
-qκ

Force Term

Choosing Rayleigh's Dissipation Function,
F=
b
2
2
v
, with
b
, the dissipation constant, a third free variable, akin to κ, the spring constant we already have, produces the following:
In[]:=
F[q_,v_,t_]:=
b
2
2
v

Continuous Differential Equation

needing no prose:
In[]:=
deqn=deqn=((Dt[D[L[q,v,t],v],t]-D[L[q,v,t],q]+D[F[q,v,t],v])/.{Dt[m,t]->0}/.{v->Dt[q,t]}/.{q->x[t]})==0
Out[]=
κx[t]+b
′
x
[t]+m
′′
x
[t]0

Damped Solution

Mathematica easily solves the equation and we pick the only pertinent part of the expression returned by DSolve:
In[]:=
dsoln=(DSolve[deqn,x,t])[[1,1,2]]
Out[]=
Function{t},
1
2
t-
b
m
-
2
b
-4mκ
m


1
+
1
2
t-
b
m
+
2
b
-4mκ
m


2

For plotting, we must be a bit careful, because if
2
b
is smaller than
4mκ
, we get complex arguments to the exponentials, yielding oscillation, but that’s the interesting part! Only the real parts of the solution are physical.
Convince yourself that as
2
b
varies from 0 to
4mκ
, the system transitions from undamped, to underdamped (showing some oscillation), to critically damped, to overdamped (showing no oscillation). See this page on damping.
In[]:=
Manipulate[Plot[Re[dsoln[t]/.{κ->1,m->1,b->bManip,C[2]->1,C[1]->0}],{t,0,25},GridLines->Automatic],​​{bManip,0,4,Appearance->"Open"},SaveDefinitions->True]
Out[]=
​
bManip
0

Discrete Solution

Reference 2 gives us the forced discrete Euler-Lagrange Equations:
References 2 and 3 are very difficult and long. The derivations dive very deeply into differential geometry, topology, and more. No shallow reading can give us a recipe for the force terms
-
F
d
and
+
F
d
. However, from our success with the undamped simple harmonic oscillator, we can hope for an easy extension.
Start with quadrature of the unforced discrete Lagrangian according to the third, unnumbered equation on page 643 of Reference 3, and introducing another free variable, the time step
h
:
Guess that
and write
Pick the first one to start:
Total success!

References

1
.
Marsden and West, 2001, Discrete Mechanics and Variational Integrators, http://www.cds.caltech.edu/~marsden/bib/2001/09-MaWe2001/MaWe2001.pdf
2
.
Matthew West, Variational Integrators, Caltech Thesis, May 2004, https://thesis.library.caltech.edu/2492/1/west_thesis.pdf
3
.
Marin Kobilarov and Jerrold E. Marsden, Discrete Geometric Optimal Control on Lie Groups, IEEE Transactions on Robotics, Vol. 27, No. 4, August 2011.

CITE THIS NOTEBOOK

Discrete damped harmonic oscillator​
by Brian Beckman​
Wolfram Community, STAFF PICKS, May 13, 2024
​https://community.wolfram.com/groups/-/m/t/3174169