# Optimal damping of vibrating systems: dependence on initial conditions

K. Lelas<EMAIL_ADDRESS>I. Nakić<EMAIL_ADDRESS>
###### Abstract

Common criteria used for measuring performance of vibrating systems have one
thing in common: they do not depend on initial conditions of the system. In
some cases it is assumed that the system has zero initial conditions, or some
kind of averaging is used to get rid of initial conditions. The aim of this
paper is to initiate rigorous study of the dependence of vibrating systems on
initial conditions in the setting of optimal damping problems. We show that,
based on the type of initial conditions, especially on the ratio of potential
and kinetic energy of the initial conditions, the vibrating system will have
quite different behavior and correspondingly the optimal damping coefficients
will be quite different. More precisely, for single degree of freedom systems
and the initial conditions with mostly potential energy, the optimal damping
coefficient will be in the under-damped regime, while in the case of the
predominant kinetic energy the optimal damping coefficient will be in the
over-damped regime. In fact, in the case of pure kinetic initial energy, the
optimal damping coefficient is $+\infty$! Qualitatively, we found the same
behavior in multi degree of freedom systems with mass proportional damping. We
also introduce a new method for determining the optimal damping of vibrating
systems, which takes into account the peculiarities of initial conditions and
the fact that, although in theory these systems asymptotically approach
equilibrium and never reach it exactly, in nature and in experiments they
effectively reach equilibrium in some finite time.

###### keywords:

viscous damping , optimal damping , multi-degree of freedom , initial
conditions

††journal: Journal of Sound and Vibration

[1] organization=Faculty of Textile Technology, University of Zagreb,
country=Croatia [2] organization=Department of Mathematics, Faculty of
Science, University of Zagreb, country=Croatia

##  1 Introduction

If we have an multi-degree of freedom (MDOF) linear vibrating system, i.e. a
system of coupled damped oscillators, how to determine damping coefficients
that ensure optimal evanescence of free vibrations? In the literature one
finds several different criteria, typically based on frequency domain analysis
of the system, although there are also approaches based on time domain
analysis [1]. The tools used for designing the criteria include modal analysis
[2], transfer functions [3], $H_{2}$ and $H_{\infty}$ norms coming from
systems theory [4, 5] and spectral techniques [6]. A general overview of the
optimization tools for structures analysis can be found in e.g. [7]. Another
optimization criterion used is to take as optimal the damping coefficients
that minimize the (zero to infinity) time integral of the energy of the
system, averaged over all possible initial conditions corresponding to the
same initial energy [8]. This criterion was investigated widely, mostly by
mathematicians in the last two decades, more details can be found, e.g., in
references [8, 9, 10, 11, 12].

However, what is common to all these criteria is that they implicitly or
explicitly ignore the dependence of the dynamics of the system on the initial
conditions. Sometimes this is suitable, e.g. for systems with continuous
excitation, but for some, like e.g. earthquake models, it make sense to study
the free vibrations of the system with non-zero initial conditions.

Specifically, the dependence of the energy integral on the initial conditions
has not been investigated, as far as we are aware, and therefore it is not
clear how much information about the behavior of vibrating systems is lost by
taking the average over all initial conditions or by assuming zero initial
conditions and it is not clear how well the optimal damping obtained in this
way works for a specific choice of initial conditions, e.g. for an experiment
with initial conditions such that the initial energy consists only of
potential energy, etc. We have chosen to study the particular criterion of
minimizing time integral of the energy as in this case it is straightforward
to modify it to take into account the initial conditions: instead of averaging
over all possible initial conditions, we study the dependence of the time
integral of the energy of the system on initial conditions. Specifically, for
criteria based in frequency domain analysis, which are designed for forced
vibrations, it is not clear how to take into account the non-zero initial
conditions.

We will explore this dependence by considering free vibrations of single
degree of freedom (SDOF), two-degree of freedom (2-DOF) and MDOF vibrating
systems with mass proportional damping (MPD). In particular, for a SDOF,
averaging over all initial conditions gives the critical damping as optimal
[8, 10], and we show, by considering the minimization of the energy integral
without averaging over initial conditions, that damping coefficients
approximately $30\%$ less than critical to infinite are obtained as optimal,
depending on the initial conditions. We systematize all our results with
respect to the relationship between initial potential and initial kinetic
energy, e.g., for initial conditions with initial potential energy grater than
initial kinetic energy the optimal damping coefficient is in the under-damped
regime, while for initial conditions with initial kinetic energy grater than
initial potential energy we find the optimal damping deep in the over-damped
regime. We also consider the minimization of the energy integral averaged over
a subset of initial conditions and obtain a significant dependence of the
optimal damping coefficient on the selected subset. Qualitatively, we find the
same behavior in 2-DOF and MDOF systems as well.

Furthermore, we show that the minimization of the energy integral for certain
types of initial conditions does not give a satisfactory optimal damping
coefficient. Specifically, for SDOF systems, the obtained optimal damping
coefficient does not distinguish between two initial states with the same
magnitude of initial displacement and initial velocity, but which differ in
the relative sign of initial displacement and initial velocity. These initial
conditions differ significantly in the rate of energy dissipation as a
function of the damping coefficient, i.e. it is not realistic for one damping
coefficient to be optimal for both of these initial conditions. The same is
true for each individual mode of MDOF systems with respect to the signs of
initial displacements and velocities, expressed via modal coordinates. Another
disadvantage of this criterion is that, for initial conditions with purely
kinetic initial energy, it gives an infinite optimal damping coefficient,
which is not practical for experiments. Also, the energy integral is
calculated over the entire time, due to the fact that these systems
asymptotically approach equilibrium and never reach it exactly, but in nature
and experiments they effectively reach equilibrium in some finite time.

We introduce a new method for determining the optimal damping of MDOF systems,
which practically solves the aforementioned problems and gives optimal damping
coefficients that take into account the peculiarities of each initial
condition and the fact that these systems effectively reach equilibrium in
some finite time. We take that the system has effectively reached equilibrium
when its energy drops to some small fraction of the initial energy, e.g., to
the energy resolution of the measuring device with which we observe the
system. Our method is based on the determination of the damping coefficient
for which the energy of the system drops to that desired energy level the
fastest. In this paper, we focus on MPD, but, as we comment at the end of the
paper, our conclusions and the method we present applies to other types of
vibrating systems as well.

In this paper we focus on mass proportional damping so that we could
analytically perform a modal analysis and present ideas in the simplest
possible way. Everything we have done can be done in a similar fashion
analytically for the case of Rayleigh damping [13] as well as for tuned mass
damper [14, 15] . Also, it is possible to carry out this kind of analysis
numerically for systems with damping that does not allow analytical treatment.
This will be the subject of our further research.

The rest of the paper is organized as follows: Section 2 is devoted to SDOF
systems, in particular minimization of the energy integral and optimal damping
is studied for the chosen set of initial conditions. In Section 3 we analyze
2-DOF systems with MPD. MDOF systems with MPD are the subject of Section 4. In
Section 5 we propose a new optimization criterion and analyze its properties.
Section 6 summarizes important findings of the paper.

##  2 SDOF systems

Free vibrations of a SDOF linear vibrating system can be described by the
equation

$\ddot{x}(t)+2\gamma\dot{x}(t)+\omega_{0}^{2}x(t)=0,\,x(0)=x_{0},\,\dot{x}(0)=v_{0},$
(1)

where $x(t)$ denotes the displacement from the equilibrium position (set to
$x=0$) as a function of time, the dots denote time derivatives, $\gamma>0$ is
the damping coefficient (sometimes called damping ratio), $\omega_{0}$ stands
for the undamped oscillator angular frequency (sometimes called the natural
frequency of the oscillator) and $(x_{0},v_{0})$ encode the initial conditions
[16, 17]. The physical units of the displacement $x(t)$ depend on the system
being considered. For example, for a mass on a spring (or a pendulum) in
viscous fluid, when it is usually called _elongation_ , it is measured in
$[m]$, while for an RLC circuit it could either be voltage, or current, or
charge. In contrast, the units of $\gamma$ and $\omega_{0}$ are $[s^{-1}]$ for
all systems described with the SDOF model. The form of the solution to the
equation (1) depends on the relationship between $\gamma$ and $\omega_{0}$,
producing three possible regimes [16, 17]: under-damped ($\gamma<\omega_{0}$),
critically damped ($\gamma=\omega_{0}$) and over-damped ($\gamma>\omega_{0}$)
regime.

Here we would like to point out that, although it is natural to classify the
solution of SDOF into three regimes depending on the value of $\gamma$, we can
actually take one form of the solution as a unique solution valid for all
values of $\gamma>0$, $\gamma\neq\omega_{0}$,

$x(t)=e^{-\gamma t}\left(x_{0}\cos(\omega t)+\frac{v_{0}+\gamma
x_{0}}{\omega}\sin(\omega t)\right)\,,$ (2)

where $\omega=\sqrt{\omega_{0}^{2}-\gamma^{2}}$ is the (complex valued) damped
angular frequency. In order to describe the critically damped regime, one can
take the limit $\gamma\rightarrow\omega_{0}$ of the solution (2) to obtain the
general solution of the critically damped regime

$x_{c}(t)=e^{-\omega_{0}t}\left(x_{0}+(v_{0}+\omega_{0}x_{0})t\right)\,.$ (3)

Therefore, in order to calculate the energy and the time integral of the
energy, we do not need to perform separate calculations for all three regimes,
but a single calculation using the displacement given by (2) and the velocity
given by

$\dot{x}(t)=e^{-\gamma t}\left(v_{0}\cos(\omega t)-\frac{\gamma
v_{0}+\omega_{0}^{2}x_{0}}{\omega}\sin(\omega t)\right)\,.$ (4)

For simplicity, in this section we will refer to the quantity

$E(t)=\dot{x}(t)^{2}+\omega_{0}^{2}x(t)^{2}$ (5)

as the _energy_ of the system, and to the quantities $E_{K}(t)=\dot{x}(t)^{2}$
and $E_{P}(t)=\omega_{0}^{2}x(t)^{2}$ as the _kinetic energy_ and _potential
energy_ of the system respectively. The connection of the quantity (5) to the
usual expressions for the energy is straightforward, e.g., for a mass $m$ on a
spring in viscous fluid

$\mathcal{E}(t)=\frac{m}{2}E(t)\,,$ (6)

and similarly for other systems described with the SDOF model. Using (2) and
(4) in (5), we obtain

$E(t)=e^{-2\gamma t}\left(E_{0}\cos^{2}(\omega
t)+\gamma\left(\omega_{0}^{2}x_{0}^{2}-v_{0}^{2}\right)\frac{\sin(2\omega
t)}{\omega}+\left(E_{0}(\omega_{0}^{2}+\gamma^{2})+4\omega_{0}^{2}\gamma
x_{0}v_{0}\right)\frac{\sin^{2}(\omega t)}{\omega^{2}}\right)$ (7)

for the energy of the system, where $E_{0}=v_{0}^{2}+\omega_{0}^{2}x_{0}^{2}$
is the initial energy given to the system at $t=0$. Accordingly,
$E_{0K}=v_{0}^{2}$ is the initial kinetic energy and
$E_{0P}=\omega_{0}^{2}x_{0}^{2}$ is the initial potential energy. Expression
(7) is valid for both under-damped and over-damped regimes, and to obtain the
energy of the critically damped regime we take the
$\gamma\rightarrow\omega_{0}$ limit of the energy (7), and obtain

$E_{c}(t)=e^{-2\omega_{0}t}\left(E_{0}+2\omega_{0}\left(\omega_{0}^{2}x_{0}^{2}-v_{0}^{2}\right)t+2\omega_{0}^{2}\left(E_{0}+2\omega_{0}x_{0}v_{0}\right)t^{2}\right)\,.$
(8)

###  2.1 Minimization of the energy integral and optimal damping in dependence
of initial conditions

We consider the SDOF system with initially energy $E_{0}$. All possible
initial conditions that give this energy can be expressed in polar coordinates
with constant radius $r=\sqrt{E_{0}}$ and angle
$\theta=\arctan\left(\frac{v_{0}}{\omega_{0}x_{0}}\right)$, i.e. we have

$\begin{split}\omega_{0}x_{0}=r\cos\theta\\\ v_{0}=r\sin\theta\,.\end{split}$
(9)

In Fig. 1 we sketch the circle given by (9), i.e. given by all possible
initial conditions with the same energy $E_{0}$. For clarity of the
exposition, here we comment on a few characteristic points of the circle
presented in Fig. 1:

  * 1.

Initial conditions $\omega_{0}x_{0}=\pm\sqrt{E_{0}}$ and $v_{0}=0$, i.e. with
purely potential initial energy (and zero initial kinetic energy), correspond
to two points on the circle with $\theta=\\{{0,\pi\\}}$.

  * 2.

Initial conditions $\omega_{0}x_{0}=\pm\sqrt{E_{0}/2}$ and
$v_{0}=\pm\sqrt{E_{0}/2}$, i.e. with initial potential energy equal to initial
kinetic energy, correspond to four points on the circle with
$\theta=\\{\pi/4,3\pi/4,5\pi/4,7\pi/4\\}$.

  * 3.

Initial conditions $\omega_{0}x_{0}=0$ and $v_{0}=\pm\sqrt{E_{0}}$, i.e. with
purely kinetic initial energy (and zero initial potential energy), correspond
to two points on the circle with $\theta=\\{{\pi/2,3\pi/2\\}}$.

$\theta$$r=\sqrt{E_{0}}$$(\omega_{0}x_{0},v_{0})$$\omega_{0}x_{0}$$v_{0}$
Figure 1: Sketch of all possible initial conditions with the same initial
energy $E_{0}$ in the $(\omega_{0}x_{0},v_{0})$ coordinate system. Square of
the coordinates corresponds to initial potential energy
$E_{0P}=\omega_{0}^{2}x_{0}^{2}$ and initial kinetic energy $E_{0K}=v_{0}^{2}$
respectively. This representation gives us a useful visualization, e.g.: all
initial conditions with $E_{0P}>E_{0K}$ are represented by two arcs, i.e.
points with $\theta\in(-\pi/4,\pi/4)\cup(3\pi/4,5\pi/4)$ (blue dotted arcs);
initial conditions with $E_{0K}=E_{0}$ and $E_{0P}=0$ are represented by two
points on a circle with $\theta=\\{\pi/2,3\pi/2\\}$ (two red filled circles);
etc.

Using (9) in (7) and (8), we obtain the energy of the under-damped and over-
damped regime

$E(t,\theta)=E_{0}e^{-2\gamma t}\left(\cos^{2}(\omega t)+\gamma\cos
2\theta\frac{\sin(2\omega
t)}{\omega}+\left(\omega_{0}^{2}+\gamma^{2}+2\omega_{0}\gamma\sin
2\theta\right)\frac{\sin^{2}(\omega t)}{\omega^{2}}\right)\,,$ (10)

and the energy of the critically damped regime

$E_{c}(t,\theta)=E_{0}e^{-2\omega_{0}t}\left(1+2\omega_{0}(\cos
2\theta)t+2\omega_{0}^{2}\left(1+\sin 2\theta\right)t^{2}\right)\,,$ (11)

as functions of $\theta$, instead of $x_{0}$ and $v_{0}$. Now we integrate
energy (10) over all time, i.e.

$I(\gamma,\theta)=\int_{0}^{\infty}E(t)dt\,,$ (12)

and obtain

$I(\gamma,\theta)=\frac{E_{0}}{2\omega_{0}}\left(\frac{\omega_{0}^{2}+\gamma^{2}}{\gamma\omega_{0}}+\frac{\gamma}{\omega_{0}}\cos
2\theta+\sin 2\theta\right)\,.$ (13)

Integral (13) is valid for all three regimes, i.e. for any $\gamma>0$.

We note here that the energy (see (7) and (8)) is invariant to a simultaneous
change of the signs of the initial conditions, i.e. to the change
$(x_{0},v_{0})\rightarrow(-x_{0},-v_{0})$ (but not to $x_{0}\rightarrow-x_{0}$
or $v_{0}\rightarrow-v_{0}$ separately). This change of signs corresponds to
the change in angle $\theta\rightarrow\theta+\pi$, therefore, functions
$\eqref{Et2}$, (11) and (13) are all periodic in $\theta$ with period $\pi$.

In Fig. 2 we show the integral (13) for
$\gamma\in\left[0.1\omega_{0},3\omega_{0}\right]$ for three different initial
conditions, i.e. for $\theta=\\{{0,\pi/4,\pi/2\\}}$. We can see that
$I(\gamma,\theta=0)$ (red solid curve), with purely potential initial energy
and zero initial kinetic energy, attains minimum for $\gamma=0.707\omega_{0}$
(rounded to three decimal places), i.e. well in the under-damped regime. For
the initial condition with equal potential and kinetic energy,
$I(\gamma,\theta=\pi/4)$ (black dotted curve) attains minimum for
$\gamma=\omega_{0}$, i.e. at the critical damping condition. Interestingly,
for the initial condition with purely kinetic energy and zero potential
energy, $I(\gamma,\theta=\pi/2)$ (blue dashed curve) has no minimum in the
displayed range of damping coefficients, therefore here we explicitly show
this function

$I(\gamma,\theta=\pi/2)=\frac{E_{0}}{2\gamma}\,,$ (14)

and it is clear that (14) has no minimum. This is easy to understand from a
physical point of view, i.e. if all the initial energy is kinetic, the higher
the damping coefficient, the faster the energy dissipation will be.

$0.1$$0.71$$1$$2$$3$$0$$1$$2$$3$$4$$5$$\gamma/\omega_{0}$$I\left(\gamma,\theta\right)\omega_{0}/E_{0}$$\theta=0$$\theta=\pi/4$$\theta=\pi/2$
Figure 2: Integral (13) for three initial conditions
$\theta=\\{0,\pi/4,\pi/2\\}$.

If we consider the optimal damping as the one for which the integral (13) is
minimal, we can easily determine the optimal damping coefficient
$\gamma_{\mathrm{opt}}(\theta)$ from the condition

$\frac{\partial
I(\gamma,\theta)}{\partial\gamma}\bigg{|}_{\gamma_{\mathrm{opt}}}=0\,,$ (15)

and we obtain

$\gamma_{\mathrm{opt}}(\theta)=\sqrt{\frac{1}{2\cos^{2}\theta}}\omega_{0}\,.$
(16)

In Fig. 3 we show the optimal damping coefficient (16) for
$\theta\in\left[0,2\pi\right]$ (function (16) has a period $\pi$, but here we
choose this interval for completeness), and here we comment on the shown
results with respect to the relationship between initial potential energy
($E_{0P}=\omega_{0}^{2}x_{0}^{2}$) and initial kinetic energy
($E_{0K}=v_{0}^{2}$) for any given initial condition, i.e. for any $\theta$:

  * 1.

Initial conditions with $E_{0P}>E_{0K}$ correspond to the set
$\theta\in\left(-\pi/4,\pi/4\right)\cup\left(3\pi/4,5\pi/4\right)$. For these
initial conditions, optimal damping coefficients (16) are in the under-damped
regime, i.e.
$\gamma_{\mathrm{opt}}\in\left[\sqrt{2}\omega_{0}/2,\omega_{0}\right)$, with
the minimum value $\gamma_{\mathrm{opt}}=\sqrt{2}\omega_{0}/2=0.707\omega_{0}$
(rounded to three decimal places) obtained for $\theta=\\{0,\pi\\}$, i.e. for
two initial conditions with $E_{0}=E_{0P}$ and $E_{0K}=0$.

  * 2.

Initial conditions with $E_{0P}=E_{0K}$ correspond to four points
$\theta=\\{\pi/4,3\pi/4,5\pi/4,7\pi/4\\}$ with optimal damping coefficient
(16) equal to critical damping, i.e. $\gamma_{\mathrm{opt}}=\omega_{0}$.

  * 3.

Initial conditions with $E_{0K}>E_{0P}$ correspond to the set
$\theta\in\left(\pi/4,3\pi/4\right)\cup\left(5\pi/4,7\pi/4\right)$. For these
initial conditions, optimal damping coefficients (16) are in the over-damped
regime, i.e. $\gamma_{\mathrm{opt}}\in\left(\omega_{0},\infty\right)$, where
$\gamma_{\mathrm{opt}}$ diverges for $\theta=\\{\pi/2,3\pi/2\\}$, i.e. for two
initial conditions with $E_{0K}=E_{0}$ and $E_{0P}=0$.

$0$$0.25$$0.5$$0.75$$1$$1.25$$1.5$$1.75$$2$$0$$1$$2$$3$$4$$\gamma_{\mathrm{opt}}=\omega_{0}$$\theta/\pi$$\gamma_{\mathrm{opt}}\left(\theta\right)/\omega_{0}$
Figure 3: Optimal damping coefficient (16) (solid red curve) as a function of
all possible initial conditions, i.e. for $\theta\in\left[0,2\pi\right]$.
Below the dashed horizontal line, optimal damping coefficients are in the
under-damped regime, above the line in the over-damped regime, and in the
critically damped regime at the crossing points of the line and the solid red
curve.

Before closing this subsection, we would like to point out two more ways in
which we can write relation (16) that will prove useful when dealing with MDOF
systems. The ratio of the initial potential energy to the initial total energy
is

$\beta=\frac{E_{0P}}{E_{0}}=\cos^{2}\theta\,,$ (17)

where we used first of the relations (9) and $E_{0P}=\omega_{0}^{2}x_{0}^{2}$.
Using (17), optimal damping coefficient (16) can be written as a function of
the fraction of potential energy in the initial total energy, i.e.

$\gamma_{\mathrm{opt}}(\beta)=\sqrt{\frac{1}{2\beta}}\omega_{0}\,.$ (18)

Thus, from (18) one can simply see that $\gamma_{\mathrm{opt}}$ is in the
under-damped regime for $\beta\in(1/2,1]$, in the critically damped regime for
$\beta=1/2$ and in the over-damped regime for $\beta\in[0,1/2)$. Using
$\beta=\omega_{0}^{2}x_{0}^{2}/E_{0}$ in (18) we can express the optimal
damping coefficient in yet another way, as a function of the initial
displacement $x_{0}$, i.e.

$\gamma_{\mathrm{opt}}(x_{0})=\sqrt{\frac{E_{0}}{2x_{0}^{2}}}=\sqrt{\frac{v_{0}^{2}+\omega_{0}^{2}x_{0}^{2}}{2x_{0}^{2}}}\,,$
(19)

where $x_{0}\in[-\sqrt{E_{0}}/\omega_{0},\sqrt{E_{0}}/\omega_{0}]$ and for
$v_{0}$ the condition $v_{0}^{2}=E_{0}-\omega_{0}^{2}x_{0}^{2}$ holds. One of
the benefits of relation (19) is that it can be seen most directly that the
optimal damping coefficient does not distinguish initial conditions $(\pm
x_{0},\pm v_{0})$ and $(\pm x_{0},\mp v_{0})$, which is a shortcoming of this
optimization criterion, because the energy as a function of time is not the
same for those two types of initial conditions (see (7) and (8)) and the
energy decay may differ significantly depending on which of those initial
conditions is in question. We will deal with these and other issues of energy
integral minimization as an optimal damping criterion in the subsection 4.2.

###  2.2 Minimization of the energy integral averaged over a set of initial
conditions and optimal damping in dependence of the chosen set

Now we calculate the average of the integral (12) over a set of initial
conditions with $\theta\in\left[\phi_{1},\phi_{2}\right]$, i.e.

$\overline{I}(\gamma,\phi_{1},\phi_{2})=\frac{1}{\phi_{2}-\phi_{1}}\int_{\phi_{1}}^{\phi_{2}}I(\gamma,\theta)d\theta\,,$
(20)

and we obtain

$\overline{I}(\gamma,\phi_{1},\phi_{2})=\frac{E_{0}}{2\omega_{0}}\left(\frac{\omega_{0}^{2}+\gamma^{2}}{\gamma\omega_{0}}+\frac{\gamma}{2\omega_{0}(\phi_{2}-\phi_{1})}(\sin
2\phi_{2}-\sin 2\phi_{1})+\frac{1}{2(\phi_{2}-\phi_{1})}(\cos 2\phi_{1}-\cos
2\phi_{2})\right)\,.$ (21)

In Fig. 4 we show averaged integral (21) for three different sets of initial
conditions. For the set of initial conditions with $\phi_{1}=-\pi/4$ and
$\phi_{2}=\pi/4$, i.e. with $E_{0P}\geq E_{0K}$ (where the equality holds only
at the end points of the set), minimum of the averaged integral (solid red
curve) is at $\gamma=0.781\omega_{0}$ (rounded to three decimal places). For
the set of initial conditions with $\phi_{1}=\pi/4$ and $\phi_{2}=3\pi/4$,
i.e. with $E_{0K}\geq E_{0P}$ (where the equality holds only at the end points
of the set), minimum of the averaged integral (dashed blue curve) is at
$\gamma=1.658\omega_{0}$ (rounded to three decimal places). For the set of
mixed initial conditions with $\phi_{1}=-\pi/4$ and $\phi_{2}=3\pi/4$, i.e.
with $E_{0P}>E_{0K}$ and $E_{0K}>E_{0P}$ points equally present in the set,
minimum of the averaged integral (dotted black curve) is at the critical
damping condition $\gamma=\omega_{0}$.

$0.1$$0.78$$1$$1.66$$2$$3$$0$$1$$2$$3$$4$$5$$\gamma/\omega_{0}$$\overline{I}\left(\gamma,\phi_{1},\phi_{2}\right)\omega_{0}/E_{0}$$\phi_{1}=-\pi/4$
, $\phi_{2}=\pi/4$$\phi_{1}=\pi/4$ , $\phi_{2}=3\pi/4$$\phi_{1}=-\pi/4$ ,
$\phi_{2}=3\pi/4$ Figure 4: Averaged integral (21) for three sets of initial
conditions.

If we consider the optimal damping as the one for which the averaged integral
(21) is minimal, we can easily determine the optimal damping coefficient
$\overline{\gamma}_{\mathrm{opt}}(\phi_{1},\phi_{2})$ form the condition

$\frac{\partial\overline{I}(\gamma,\phi_{1},\phi_{2})}{\partial\gamma}\bigg{|}_{\overline{\gamma}_{\mathrm{opt}}}=0\,,$
(22)

and we obtain

$\overline{\gamma}_{\mathrm{opt}}(\phi_{1},\phi_{2})=\sqrt{\frac{2(\phi_{2}-\phi_{1})}{2(\phi_{2}-\phi_{1})+\sin
2\phi_{2}-\sin 2\phi_{1}}}\omega_{0}\,.$ (23)

We note here that averaged integral (21) and optimal damping coefficient (23)
are not periodic functions in variables $\phi_{1}$ and $\phi_{2}$, if we keep
one variable fixed and change the other. But they are periodic, with period
$\pi$, if we change both variables simultaneously.

In Fig. 5 we show the optimal damping coefficient (23) as a function of
$\phi_{2}$ with fixed $\phi_{1}=0$, and the results shown can be summarized as
follows:

  * 1.

For $\phi_{1}=0$ and
$\phi_{2}\in\left[0,\pi/2\right)\cup\left(\pi,3\pi/2\right)$, the optimal
damping coefficient (23) is in the under-damped regime. In this case, integral
(20) is averaged over sets that have more points corresponding to initial
conditions with $E_{0P}>E_{0K}$, in comparison to the points corresponding to
initial conditions with $E_{0K}>E_{0P}$.

  * 2.

For $\phi_{1}=0$ and $\phi_{2}=\\{\pi/2,\pi,3\pi/2,2\pi\\}$, the optimal
damping coefficient (23) is equal to critical damping. In this case, integral
(20) is averaged over sets that have equal amount of points corresponding to
initial conditions with $E_{0P}>E_{0K}$ and initial conditions with
$E_{0K}>E_{0P}$.

  * 3.

For $\phi_{1}=0$ and
$\phi_{2}\in\left(\pi/2,\pi\right)\cup\left(3\pi/2,2\pi\right)$, the optimal
damping coefficient (23) is in the over-damped regime. In this case, integral
(20) is averaged over sets that have more points corresponding to initial
conditions with $E_{0K}>E_{0P}$, in comparison to the points corresponding to
initial conditions with $E_{0P}>E_{0K}$.

$0$$0.25$$0.5$$0.75$$1$$1.25$$1.5$$1.75$$2$$0.7$$0.8$$0.9$$1$$1.1$$\overline{\gamma}_{\mathrm{opt}}=\omega_{0}$$\phi_{2}/\pi$$\overline{\gamma}_{\mathrm{opt}}\left(0,\phi_{2}\right)/\omega_{0}$
Figure 5: Optimal damping coefficient (23) (solid red curve) as a function of
$\phi_{2}\in\left[0,2\pi\right]$ for fixed $\phi_{1}=0$. Below the dashed
horizontal line, optimal damping coefficients are in the under-damped regime,
above the line in the over-damped regime, and in the critically damped regime
at the crossing points of the line and the solid red curve.

##  3 2-DOF systems with MPD

Figure 6: Schematic figure of a 2-DOF system.

Here we consider 2-DOF system shown schematically in Fig. 6. The corresponding
equations of motion are

$\begin{split}m_{1}\ddot{x}_{1}(t)=-c_{1}\dot{x}_{1}(t)-k_{1}x_{1}(t)-k_{2}\left(x_{1}(t)-x_{2}(t)\right)\,,\\\
m_{2}\ddot{x}_{2}(t)=-c_{2}\dot{x}_{2}(t)-k_{3}x_{2}(t)+k_{2}\left(x_{1}(t)-x_{2}(t)\right)\,.\end{split}$
(24)

We will consider MPD [18], i.e. masses $\\{m_{1},m_{2}\\}$, spring constants
$\\{k_{1},k_{2},k_{3}\\}$, and dampers $\\{c_{1},c_{2}\\}$ can in general be
mutually different but the condition $c_{1}/m_{1}=c_{2}/m_{2}$ holds. In this
case we can use modal analysis [2, 16] and the system of equations (24) can be
written via modal coordinates [16] as

$\begin{split}\ddot{q}_{1}(t)+2\gamma\dot{q}_{1}(t)+\omega_{01}^{2}q_{1}(t)=0\,\\\
\ddot{q}_{2}(t)+2\gamma\dot{q}_{2}(t)+\omega_{02}^{2}q_{2}(t)=0\,,\end{split}$
(25)

where $q_{i}(t)$ and $\omega_{0i}$, with $i=\\{1,2\\}$, denote the modal
coordinates and undamped modal frequencies of the two modes, while
$\gamma=c_{i}/2m_{i}$ is the damping coefficient. In the analysis that we will
carry out in this subsection, we will not need the explicit connection of
modal coordinates $q_{i}(t)$ and mass coordinates, i.e. displacements
$x_{i}(t)$, and we will deal with this in the next subsection when considering
a specific example with given masses, springs and dampers. Similarly as in
section 2 (see (2)), the general solution for the $i$-th mode can be written
as

$q_{i}(t)=e^{-\gamma t}\left(q_{0i}\cos(\omega_{i}t)+\frac{\dot{q}_{0i}+\gamma
q_{0i}}{\omega_{i}}\sin(\omega_{i}t)\right)\,,$ (26)

where $\omega_{i}=\sqrt{\omega_{0i}^{2}-\gamma^{2}}$ is the damped modal
frequency, and $q_{i}(0)\equiv q_{0i}$ and $\dot{q}_{i}(0)\equiv\dot{q}_{0i}$
are the initial conditions of the $i$-th mode. Thus, the reasoning and the
results presented in section 2, with some adjustments, can by applied for the
analysis of the 2-DOF system we are considering here.

The energy of the system is

$E(t)=\sum_{i=1}^{2}\frac{m_{i}\dot{x}_{i}(t)^{2}}{2}+\frac{k_{1}x_{1}(t)^{2}}{2}+\frac{k_{3}x_{2}(t)^{2}}{2}+\frac{k_{2}(x_{1}(t)-x_{2}(t))^{2}}{2}\,,$
(27)

and we take that the modal coordinates are normalised so that (27) can be
written as

$E(t)=\sum_{i=1}^{2}E_{i}(t)=\sum_{i=1}^{2}\left(\dot{q}_{i}(t)^{2}+\omega_{0i}^{2}q_{i}(t)^{2}\right)\,$
(28)

where $E_{i}(t)$ in (28) denotes the energy of the $i$-th mode. Total energy
at $t=0$, i.e. the initial energy, is given by

$E_{0}=\sum_{i=1}^{2}E_{0i}=\sum_{i=1}^{2}\left(E_{0Ki}+E_{0Pi}\right)=\sum_{i=1}^{2}\left(\dot{q}_{0i}^{2}+\omega_{0i}^{2}q_{0i}^{2}\right)\,,$
(29)

where $E_{0i}$ denotes the initial energy of the $i$-th mode,
$E_{0Ki}=\dot{q}_{0i}^{2}$ and $E_{0Pi}=\omega_{0i}^{2}q_{0i}^{2}$ denote
initial kinetic and initial potential energy of the $i$-th mode.

All possible initial conditions with the same initial energy (29) can be
expressed similarly as in the SDOF case (see (9) and Fig. 1) but with two
pairs of polar coordinates, one pair for each mode. For the $i$-th mode we
have radius $r_{i}=\sqrt{E_{0i}}$ and angle
$\theta_{i}=\arctan\left(\frac{\dot{q}_{0i}}{\omega_{0i}q_{0i}}\right)$, i.e.
we can write

$\begin{split}\omega_{0i}q_{0i}=r_{i}\cos\theta_{i}\\\
\dot{q}_{0i}=r_{i}\sin\theta_{i}\,.\end{split}$ (30)

Thus, each initial condition with energy $E_{0}=E_{01}+E_{02}$ can be
represented by points on two circles with radii $r_{1}=\sqrt{E_{01}}$ and
$r_{2}=\sqrt{E_{02}}$, for which condition $r_{1}^{2}+r_{2}^{2}=E_{0}$ holds,
and with angles $\theta_{1}$ and $\theta_{2}$ that tell us how initial
potential and initial kinetic energy are distributed within the modes. Using
relation (10) for SDOF systems, we can write the energy of the $i$-th mode in
polar coordinates (30) as

$E_{i}(t)=E_{0i}e^{-2\gamma t}\left(\cos^{2}(\omega_{i}t)+\gamma\cos
2\theta_{i}\frac{\sin(2\omega_{i}t)}{\omega_{i}}+\left(\omega_{0i}^{2}+\gamma^{2}+2\omega_{0i}\gamma\sin
2\theta_{i}\right)\frac{\sin^{2}(\omega_{i}t)}{\omega_{i}^{2}}\right)\,$ (31)

for the under-damped ($\gamma<\omega_{0i}$) and over-damped
($\gamma>\omega_{0i}$) regime, and the energy of the $i$-th mode in the
critically damped regime is obtained analogously using the relation (11).

Consequently, the integral of the energy (28) over the entire time, for some
arbitrary initial condition, is simply calculated using relation (13). For
each individual mode, we obtain

$I(\gamma,\\{E_{0i}\\},\\{\theta_{i}\\})=\sum_{i=1}^{2}\int_{0}^{\infty}E_{i}(t)dt=\sum_{i=1}^{2}\frac{E_{0i}}{2\omega_{0i}}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}}+\frac{\gamma}{\omega_{0i}}\cos
2\theta_{i}+\sin 2\theta_{i}\right)\,.$ (32)

Furthermore, initial energy of the $i$-th mode can be written as
$E_{0i}=a_{i}^{2}E_{0}$, where coefficient $a_{i}^{2}\in[0,1]$ denotes the
fraction of the initial energy of the $i$-th mode in the total initial energy.
Coefficients of the two modes satisfy $a_{1}^{2}+a_{2}^{2}=1$ and therefore
can be parameterized as

$\begin{split}a_{1}=\cos\psi\\\ a_{2}=\sin\psi,\end{split}$ (33)

where $\psi\in[0,\pi/2]$. Taking (33) into account, we can write (32) as

$I(\gamma,\psi,\theta_{1},\theta_{2})=E_{0}\sum_{i=1}^{2}\frac{a_{i}^{2}}{2\omega_{0i}}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}}+\frac{\gamma}{\omega_{0i}}\cos
2\theta_{i}+\sin 2\theta_{i}\right)\,.$ (34)

If we consider the optimal damping coefficient as the one for which the
integral (34) is minimal, we can easily determine the optimal damping
coefficient form the condition

$\frac{\partial
I(\gamma,\psi,\theta_{1},\theta_{2})}{\partial\gamma}\bigg{|}_{\gamma_{\mathrm{opt}}}=0\,,$
(35)

and we obtain

$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})=\sqrt{\frac{\omega_{01}^{2}\omega_{02}^{2}}{2\omega_{02}^{2}\cos^{2}\psi\cos^{2}\theta_{1}+2\omega_{01}^{2}\sin^{2}\psi\cos^{2}\theta_{2}}}\,.$
(36)

It is easy to see that, for any fixed $\psi$, the function (36) has smallest
magnitude for $\cos^{2}\theta_{1}=\cos^{2}\theta_{2}=1$, which corresponds to
the initial conditions with initial energy comprised only of potential energy
distributed within the two modes, i.e $E_{0}=E_{0P1}+E_{0P2}$. In that case we
can write the denominator of (36) as

$f(\psi)=\sqrt{2\omega_{02}^{2}\cos^{2}\psi+2\omega_{01}^{2}\sin^{2}\psi}=\sqrt{2(\omega_{02}^{2}-\omega_{01}^{2})\cos^{2}\psi+2\omega_{01}^{2}}\,,$
(37)

where we used $\sin^{2}\psi=1-\cos^{2}\psi$. Since $\omega_{01}<\omega_{02}$,
the function (37) has maximum for $\psi=0$. Thus, the minimum value of the
optimal damping coefficient (36) is $\sqrt{2}\omega_{01}/2$, and it is
obtained for $\psi=0$ and $\theta_{1}=\\{0,\pi\\}$, which corresponds to the
initial conditions with initial energy comprised only of potential energy in
the first mode, i.e. $E_{0}=E_{0P1}$. On the other hand, for any fixed $\psi$,
the function (36) has singularities for
$\cos^{2}\theta_{1}=\cos^{2}\theta_{2}=0$, which corresponds to the initial
conditions with initial energy comprised only of kinetic energy. Thus, the
range of the optimal damping coefficient (36) is

$\gamma_{\mathrm{opt}}\in\left[\sqrt{2}\omega_{01}/2,+\infty\right)\,.$ (38)

Now we calculate the average of the integral (34) over a set of all initial
conditions, we obtain

$\overline{I}(\gamma)=\frac{1}{2\pi^{3}}\int_{0}^{\pi/2}d\psi\int_{0}^{2\pi}d\theta_{1}\int_{0}^{2\pi}d\theta_{2}\,I(\gamma,\psi,\theta_{1},\theta_{2})=\frac{E_{0}}{4}\sum_{i=1}^{2}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}^{2}}\right)\,,$
(39)

and from the condition

$\frac{\partial\overline{I}(\gamma)}{\partial\gamma}\bigg{|}_{\overline{\gamma}_{\mathrm{opt}}}=0\,,$
(40)

we find that the optimal damping coefficient with respect to the averaged
integral (39) is given by

$\overline{\gamma}_{\mathrm{opt}}=\sqrt{\frac{2\omega_{01}^{2}\omega_{02}^{2}}{\omega_{01}^{2}+\omega_{02}^{2}}}\,.$
(41)

In order to more easily analyze the behavior of the damping coefficient (36)
with regard to the distribution of the initial potential energy within the
modes and its relationship with the damping coefficient (41), similarly as in
subsection 2.1 (see (17) and (18)), we define the ratio of the initial
potential energy of the $i$-th mode and the total initial energy, i.e.

$\beta_{i}=\frac{E_{0Pi}}{E_{0}}\,.$ (42)

Since the initial potential energy satisfies $E_{0P}=E_{0P1}+E_{0P2}\leq
E_{0}$, we have $\beta_{i}\in[0,1]$ and the condition
$0\leq(\beta_{1}+\beta_{2})\leq 1$ holds. Taking
$E_{0Pi}=E_{0i}\cos^{2}\theta_{i}$ (see (30)) and $E_{0i}=a_{i}^{2}E_{0}$ with
(33) into account, we have

$\begin{split}\beta_{1}=\cos^{2}\psi\cos^{2}\theta_{1}\\\
\beta_{2}=\sin^{2}\psi\cos^{2}\theta_{2}.\end{split}$ (43)

Using (43), relation (36) can be written as

$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})=\sqrt{\frac{\omega_{01}^{2}\omega_{02}^{2}}{2\omega_{02}^{2}\beta_{1}+2\omega_{01}^{2}\beta_{2}}}\,.$
(44)

For clarity, we will repeat briefly, the minimum value of (44) is
$\sqrt{2}\omega_{01}/2$, obtained for $\beta_{1}=1$ and $\beta_{2}=0$ (or in
terms of the angles in (36), for $\psi=0$ and $\theta_{1}=\\{0,\pi\\}$), while
$\gamma_{\mathrm{opt}}\rightarrow+\infty$ for $\beta_{1}=\beta_{2}=0$ (or in
terms of the angles in (36), for any $\psi$ with
$\theta_{1}=\\{\pi/2,3\pi/2\\}$ and $\theta_{2}=\\{\pi/2,3\pi/2\\}$). The
benefit of relation (44) is that we expressed (36) through two variables
instead of three, i.e. this way we lost information about the signs of the
initial conditions and about distribution of initial kinetic energy within the
modes, but the optimal damping coefficient (36) does not depend on those signs
anyway, due to the squares of trigonometric functions in variables
$\theta_{1}$ and $\theta_{2}$, and, for a fixed distribution of initial
potential energy within the modes, the optimal damping coefficient (36) is
constant for different distributions of initial kinetic energy within the
modes. By looking at relations (44) and (41), it is immediately clear that
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})=\overline{\gamma}_{\mathrm{opt}}$
for

$\omega_{02}^{2}\beta_{1}+\omega_{01}^{2}\beta_{2}=\frac{\omega_{01}^{2}+\omega_{02}^{2}}{4},$
(45)

while
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})<\overline{\gamma}_{\mathrm{opt}}$
if the left hand side of relation (45) is greater than the right hand side,
and
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})>\overline{\gamma}_{\mathrm{opt}}$
if the left hand side of relation (45) is smaller than the right hand side.

Again, similarly as in subsection 2.1 (see (19)), using
$\beta_{i}=\omega_{0i}^{2}q_{0i}^{2}/E_{0}$ we can express the optimal damping
coefficient (44) as a function of the initial modal coordinates as well, i.e.

$\gamma_{\mathrm{opt}}(q_{01},q_{02})=\sqrt{\frac{E_{0}}{2q_{01}^{2}+2q_{02}^{2}}}\,,$
(46)

where $q_{0i}\in[-\sqrt{E_{0}}/\omega_{0i},\sqrt{E_{0}}/\omega_{0i}]$ and the
condition $0\leq(\omega_{01}^{2}q_{01}^{2}+\omega_{02}^{2}q_{02}^{2})\leq
E_{0}$ holds. We can express condition (45) in terms of initial modal
coordinates, i.e.
$\gamma_{\mathrm{opt}}(\\{q_{0i}\\})=\overline{\gamma}_{\mathrm{opt}}$ for

$\frac{q_{01}^{2}+q_{02}^{2}}{E_{0}}=\frac{\omega_{01}^{2}+\omega_{02}^{2}}{4\omega_{01}^{2}\omega_{02}^{2}},$
(47)

while $\gamma_{\mathrm{opt}}(\\{q_{0i}\\})<\overline{\gamma}_{\mathrm{opt}}$
if the left hand side of relation (47) is greater than the right hand side,
and $\gamma_{\mathrm{opt}}(\\{q_{0i}\\})>\overline{\gamma}_{\mathrm{opt}}$ if
the left hand side of relation (47) is smaller than the right hand side.

We note here that we did not use explicit values of the undamped modal
frequencies $\omega_{01}$ and $\omega_{02}$ in the analysis so far, and
relations presented so far are valid for any 2-DOF system with MPD. In the
next subsection, we provide a more detailed quantitative analysis using an
example with specific values of modal frequencies.

###  3.1 Quantitative example

Here we consider the 2-DOF system as the one shown schematically in Fig. 6,
but with $m_{1}=m_{2}=m$, $k_{1}=k_{2}=k_{3}=k$ and $c_{1}=c_{2}=c$. The
corresponding equations of motion are

$\begin{split}m\ddot{x}_{1}(t)=-c\dot{x}_{1}(t)-kx_{1}(t)-k\left(x_{1}(t)-x_{2}(t)\right)\,,\\\
m\ddot{x}_{2}(t)=-c\dot{x}_{2}(t)-kx_{2}(t)+k\left(x_{1}(t)-x_{2}(t)\right)\,.\end{split}$
(48)

For completeness, we will investigate here the behavior of the optimal damping
coefficient given by the minimization of the energy integral for different
initial conditions, and its relationship with the optimal damping coefficient
given by the minimization of the averaged energy integral, in all three
coordinate systems that we introduced in the previous subsection and
additionally in the coordinate system defined by the initial displacements of
the masses. System of equations (48) can be easily recast to the form (25)
with the modal coordinates

$\begin{split}q_{1}(t)=\sqrt{\frac{m}{4}}\left(x_{1}(t)+x_{2}(t)\right)\\\
q_{2}(t)=\sqrt{\frac{m}{4}}\left(x_{1}(t)-x_{2}(t)\right)\,,\end{split}$ (49)

and with the natural (undamped) frequencies of the modes
$\omega_{01}=\omega_{0}$ and $\omega_{02}=\sqrt{3}\omega_{0}$, where
$\omega_{0}=\sqrt{k/m}$. Normalisation factors $\sqrt{m/4}$ in (49) ensure
that our expression (28) for the energy of the system corresponds to energy
expressed over the displacements and velocities of the masses, i.e.

$E(t)=\sum_{i=1}^{2}\left(\dot{q}_{i}(t)^{2}+\omega_{0i}^{2}q_{i}(t)^{2}\right)=\sum_{i=1}^{2}\left(\frac{m\dot{x}_{i}(t)^{2}}{2}+\frac{kx_{i}(t)^{2}}{2}\right)+\frac{k(x_{1}(t)-x_{2}(t))^{2}}{2}\,.$
(50)

Using the specific values of undamped modal frequencies of this system,
relations (36), (41) and (44) become

$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})=\sqrt{\frac{3}{6\cos^{2}\psi\cos^{2}\theta_{1}+2\sin^{2}\psi\cos^{2}\theta_{2}}}\omega_{0}\,,$
(51) $\overline{\gamma}_{\mathrm{opt}}=\frac{\sqrt{6}}{2}\omega_{0}\,,$ (52)
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})=\sqrt{\frac{3}{6\beta_{1}+2\beta_{2}}}\omega_{0}\,.$
(53)

Since $\omega_{01}=\omega_{0}$, the range of (53) is
$\gamma_{\mathrm{opt}}\in[\sqrt{2}\omega_{0}/2,+\infty)$ (see (38)).

As examples of the behavior of the damping coefficient (51) as a function of
the angles $\\{\psi,\theta_{1},\theta_{2}\\}$ and its relationship with the
damping coefficient (52), in Fig. 7 we show
$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})/\overline{\gamma}_{\mathrm{opt}}$
for $\psi=\\{\pi/3,\pi/6\\}$ and $\theta_{i}\in[0,\pi]$.

Figure 7: Ratio
$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})/\overline{\gamma}_{\mathrm{opt}}$
of the optimal damping coefficients (51) and (52) for $\psi=\\{\pi/3,\pi/6\\}$
and $\theta_{i}\in[0,\pi]$. (a) For $\psi=\pi/3$, the total initial energy
$E_{0}$ is distributed within the modes as $E_{01}=E_{0}/4$ and
$E_{02}=3E_{0}/4$. (b) For $\psi=\pi/6$, the total initial energy is
distributed within the modes as $E_{01}=3E_{0}/4$ and $E_{02}=E_{0}/4$.
Singularities for $\theta_{1}=\theta_{2}=\pi/2$ are indicated by infinity
symbols, and the points around singularities for which
$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})/\overline{\gamma}_{\mathrm{opt}}>3$
are removed on both figures (central white areas). Black lines, on both
figures, indicate the points for which
$\gamma_{\mathrm{opt}}(\psi,\theta_{1},\theta_{2})/\overline{\gamma}_{\mathrm{opt}}=1$.
On both figures, ratio attains minimum for the corner points, i.e. for
$(\theta_{1},\theta_{2})=\\{(0,0),(0,\pi),(\pi,0),(\pi,\pi)\\}$.

In Fig. 8 we show ratio of the damping coefficient (53) and the damping
coefficient (52), i.e.
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})/\overline{\gamma}_{\mathrm{opt}}$.

Figure 8: Ratio
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})/\overline{\gamma}_{\mathrm{opt}}$
of the optimal damping coefficients (53) and (52) for $\beta_{1}\in[0,1]$,
$\beta_{2}\in[0,1]$ and the constraint $0\leq(\beta_{1}+\beta_{2})\leq 1$.
Singularity for $\beta_{1}=\beta_{2}=0$ is indicated by the filled red circle,
and the points near singularity, for which
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})/\overline{\gamma}_{\mathrm{opt}}>3$,
are removed, thus, a small white triangle is formed with the right angle at
the origin. Black line indicates the points for which
$\gamma_{\mathrm{opt}}(\beta_{1},\beta_{2})/\overline{\gamma}_{\mathrm{opt}}=1$.
The minimum value of the ratio is at the point $(\beta_{1},\beta_{2})=(1,0)$.

If the initial energy is comprised only of potential energy, in terms of
initial modal coordinates we have
$E_{0}=\omega_{0}^{2}q_{01}^{2}+3\omega_{0}^{2}q_{02}^{2}$, thus, the initial
modal coordinates satisfy

$\begin{split}q_{01}\sqrt{\frac{\omega_{0}^{2}}{E_{0}}}\in\left[-1,1\right]\,,\\\
q_{02}\sqrt{\frac{\omega_{0}^{2}}{E_{0}}}\in\left[-\sqrt{3}/3,\sqrt{3}/3\right]\,,\\\
0\leq\frac{\omega_{0}^{2}}{E_{0}}\left(q_{01}^{2}+3q_{02}^{2}\right)\leq
1\,.\end{split}$ (54)

Furthermore, we can write the optimal damping coefficient (46) as

$\gamma_{\mathrm{opt}}(q_{01},q_{02})=\sqrt{\frac{E_{0}}{2\omega_{0}^{2}(q_{01}^{2}+q_{02}^{2})}}\omega_{0}\,,$
(55)

and the condition (47) as

$\frac{\omega_{0}^{2}}{E_{0}}\left(q_{01}^{2}+q_{02}^{2}\right)=\frac{1}{3}.$
(56)

In Fig. 9(a) we show the ratio of (55) and (52), i.e.
$\gamma_{\mathrm{opt}}(q_{01},q_{02})/\overline{\gamma}_{\mathrm{opt}}$. The
domain of this function consists of points inside and on the ellipse, i.e. it
is given by (54). Similarly as before, singularity at $(q_{01},q_{02})=(0,0)$
is indicated by the infinity symbol, and the points for which
$\gamma_{\mathrm{opt}}(q_{01},q_{02})/\overline{\gamma}_{\mathrm{opt}}>3$ are
removed. For points inside the circle we have
$\gamma_{\mathrm{opt}}(q_{01},q_{02})/\overline{\gamma}_{\mathrm{opt}}>1$, and
for points outside the circle we have
$\gamma_{\mathrm{opt}}(q_{01},q_{02})/\overline{\gamma}_{\mathrm{opt}}<1$.
Minimum values of this ratio are $\sqrt{3}/3\approx 0.58$, obtained for the
points
$(q_{01},q_{02})=\\{(-\sqrt{E_{0}}/\omega_{0},0),(\sqrt{E_{0}}/\omega_{0},0)\\}$.

Figure 9: (a) Ratio
$\gamma_{\mathrm{opt}}(q_{01},q_{02})/\overline{\gamma}_{\mathrm{opt}}$ of the
optimal damping coefficients (55) and (52). (b) Ratio
$\gamma_{\mathrm{opt}}(x_{01},x_{02})/\overline{\gamma}_{\mathrm{opt}}$ of the
optimal damping coefficients (57) and (52). Singularities, at points $(0,0)$
on the both figures, are denoted by infinity symbols, and the points near
singularities, for which
$\gamma_{\mathrm{opt}}/\overline{\gamma}_{\mathrm{opt}}>3$, are removed. Black
circles on both figures indicate the points for which
$\gamma_{\mathrm{opt}}/\overline{\gamma}_{\mathrm{opt}}=1$.

Using (49) we can write the optimal damping coefficient (55) in terms of
initial displacements $x_{i}(0)\equiv x_{0i}$ as

$\gamma_{\mathrm{opt}}(x_{01},x_{02})=\sqrt{\frac{E_{0}}{m(x_{01}^{2}+x_{02}^{2})}}=\sqrt{\frac{E_{0}}{m\omega_{0}^{2}(x_{01}^{2}+x_{02}^{2})}}\omega_{0}\,.$
(57)

If the initial energy is comprised only of potential energy, in terms of
initial displacements we have
$E_{0}=m\omega_{0}^{2}(x_{01}^{2}+x_{02}^{2}-x_{01}x_{02})$, thus, the initial
displacements of the masses satisfy

$\begin{split}x_{0i}\sqrt{\frac{m\omega_{0}^{2}}{E_{0}}}\in\left[-1,1\right]\,,\\\
0\leq\frac{m\omega_{0}^{2}}{E_{0}}\left(x_{01}^{2}+x_{02}^{2}-x_{01}x_{02}\right)\leq
1\,,\end{split}$ (58)

and the condition (56) is now

$\frac{m\omega_{0}^{2}}{E_{0}}\left(x_{01}^{2}+x_{02}^{2}\right)=\frac{2}{3}.$
(59)

In Fig. 9(b) we show the ratio of (57) and (52), i.e.
$\gamma_{\mathrm{opt}}(x_{01},x_{02})/\overline{\gamma}_{\mathrm{opt}}$. The
domain of this function consists of points given by (58). Similarly as before,
singularity at $(x_{01},x_{02})=(0,0)$ is indicated by the infinity symbol,
and the points for which
$\gamma_{\mathrm{opt}}(x_{01},x_{02})/\overline{\gamma}_{\mathrm{opt}}>3$ are
removed. For points inside the circle
$\gamma_{\mathrm{opt}}(x_{01},x_{02})/\overline{\gamma}_{\mathrm{opt}}>1$, and
for points outside the circle
$\gamma_{\mathrm{opt}}(x_{01},x_{02})/\overline{\gamma}_{\mathrm{opt}}<1$.
Minimum values of this ratio are $\sqrt{3}/3\approx 0.58$, obtained for the
points
$(x_{01},x_{02})=\left(\pm\sqrt{\frac{E_{0}}{m\omega_{0}^{2}}},\pm\sqrt{\frac{E_{0}}{m\omega_{0}^{2}}}\right)$.

##  4 MDOF systems with MPD

Figure 10: Schematic figure of a MDOF system with $N$ degrees of freedom.

Here we consider the MDOF system with $N$ degrees of freedom shown
schematically in Fig. 10. As in the section 3, we will consider MPD, i.e.
masses $\\{m_{1},m_{2},...,m_{N}\\}$, spring constants
$\\{k_{1},k_{2},...,k_{N+1}\\}$, and dampers $\\{c_{1},c_{2},...,c_{N}\\}$ can
in general be mutually different but the condition $c_{i}/m_{i}=2\gamma$ holds
for any $i=\\{1,...,N\\}$, where $\gamma$ is the damping coefficient.
Therefore, the reasoning we presented in section 3 can be applied here, with
the main difference that now the system has $N$ modes instead of two. Again,
we can write each initial condition over polar coordinates, as in the 2-DOF
case (see (30)), only now we have $N$ pairs of polar coordinates instead of
two.

The energy of each mode is given by (31), and consequently, the integral of
the total energy over the entire time, for some arbitrary initial condition,
is simply calculated similarly as in (32), i.e.

$I(\gamma,\\{E_{0i}\\},\\{\theta_{i}\\})=\sum_{i=1}^{N}\int_{0}^{\infty}E_{i}(t)dt=\sum_{i=1}^{N}\frac{E_{0i}}{2\omega_{0i}}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}}+\frac{\gamma}{\omega_{0i}}\cos
2\theta_{i}+\sin 2\theta_{i}\right)\,,$ (60)

where, again, $E_{i}(t)$ is the energy of the $i$-th mode, $E_{0i}$ is the
initial energy of the $i$-th mode. Thus, each initial condition with energy
$E_{0}=\sum_{i=1}^{N}E_{0i}$ is represented by points on $N$ circles with
radii $r_{i}=\sqrt{E_{0i}}$, for which condition
$\sum_{i=1}^{N}r_{i}^{2}=E_{0}$ holds, and with angles $\theta_{i}$ that tell
us how initial potential and initial kinetic energy are distributed within the
modes.

Similarly as before, initial energy of the $i$-th mode can be written as
$E_{0i}=a_{i}^{2}E_{0}$, where coefficient $a_{i}^{2}\in[0,1]$ denotes the
fraction of the initial energy of the $i$-th mode in the total initial energy
$E_{0}$, and the condition

$\sum_{i=1}^{N}a_{i}^{2}=1$ (61)

holds. Relation (61) defines a sphere embedded in $N$-dimensional space and we
can express the coefficients $a_{i}$ over $N$-dimensional spherical
coordinates ($N-1$ independent coordinates, i.e. angles, since the radius is
equal to one), but for the sake of simplicity we will not do that here and we
will stick to writing the expressions as a functions of the coefficients
$a_{i}$. Thus, we can write (60) as

$I(\gamma,\\{a_{i}\\},\\{\theta_{i}\\})=\sum_{i=1}^{N}\int_{0}^{\infty}E_{i}(t)dt=E_{0}\sum_{i=1}^{N}\frac{a_{i}^{2}}{2\omega_{0i}}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}}+\frac{\gamma}{\omega_{0i}}\cos
2\theta_{i}+\sin 2\theta_{i}\right)\,.$ (62)

We differentiate relation (62) by $\gamma$ and equate it to zero and get

$\gamma_{\mathrm{opt}}(\\{a_{i}\\},\\{\theta_{i}\\})=\left(\sum_{i=1}^{N}\frac{2a_{i}^{2}\cos^{2}\theta_{i}}{\omega_{0i}^{2}}\right)^{-1/2}\,$
(63)

as the optimal damping coefficient for which integral (62) is minimal. For any
fixed set of coefficients $\\{a_{i}\\}$, the smallest magnitude of the
function (63) is obtained for $\cos^{2}\theta_{i}=1$ $\forall i$, which
corresponds to the initial conditions with initial energy comprised only of
potential energy distributed within the modes, i.e
$E_{0}=\sum_{i=1}^{N}E_{0Pi}$. In that case the denominator of (63) is

$f(\\{a_{i}\\})=\left(\sum_{i=1}^{N}\frac{2a_{i}^{2}}{\omega_{0i}^{2}}\right)^{1/2}\,$
(64)

and using $a_{1}^{2}=1-\sum_{i=2}^{N}a_{i}^{2}$ (see (61)) we can write (64)
as

$f(\\{a_{i}\\})=\left(\frac{2}{\omega_{01}^{2}}+\sum_{i=2}^{N}2a_{i}^{2}\left(\frac{1}{\omega_{0i}^{2}}-\frac{1}{\omega_{01}^{2}}\right)\right)^{1/2}\,.$
(65)

Since $\omega_{01}<\omega_{0i}$ for any $i\geq 2$, each term in the sum of
relation (65) is negative, and we can conclude that the function (65) has
maximum for the set $\\{a_{i}\\}=\\{1,0,...,0\\}$. Thus, the minimum value of
the optimal damping coefficient (63) is $\sqrt{2}\omega_{01}/2$, and it is
obtained for $a_{1}=1$ and $\theta_{1}=\\{0,\pi\\}$, which corresponds to the
initial conditions with initial energy comprised only of potential energy in
the first mode, i.e. $E_{0}=E_{0P1}$. On the other hand, for any fixed set
$\\{a_{i}\\}$, the function (63) has singularities for $\cos^{2}\theta_{i}=0$
$\forall i$. Thus, the range of the optimal damping coefficient (63) is

$\gamma_{\mathrm{opt}}\in\left[\sqrt{2}\omega_{01}/2,+\infty\right)\,.$ (66)

In A we have calculated the average of the integral (62) over a set of all
initial conditions and obtained

$\overline{I}(\gamma)=\frac{E_{0}}{2N}\sum_{i=1}^{N}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}^{2}}\right)\,.$
(67)

We differentiate relation (67) by $\gamma$ and equate it to zero and obtain

$\overline{\gamma}_{\mathrm{opt}}=N^{1/2}\left(\sum_{i=1}^{N}\frac{1}{\omega_{0i}^{2}}\right)^{-1/2}\,$
(68)

as the optimal damping coefficient with respect to the averaged integral (67).

Since the ratio of the initial potential energy of the $i$-th mode and the
total initial energy is

$\beta_{i}=\frac{E_{0Pi}}{E_{0}}=a_{i}^{2}\cos^{2}\theta_{i}\,,$ (69)

where $\beta_{i}\in[0,1]$ and the condition $0\leq\sum_{i=1}^{N}\beta_{i}\leq
1$ holds, we can write (63) as a function of the distribution of the initial
potential energy over the modes, i.e.

$\gamma_{\mathrm{opt}}(\\{\beta_{i}\\})=\left(\sum_{i=1}^{N}\frac{2\beta_{i}}{\omega_{0i}^{2}}\right)^{-1/2}\,.$
(70)

The minimum value of (70) is $\sqrt{2}\omega_{01}/2$, obtained for
$\beta_{1}=1$ and $\beta_{i}=0$ for $i\geq 2$, while
$\gamma_{\mathrm{opt}}\rightarrow+\infty$ for $\beta_{i}=0$ $\forall i$. Using
$\beta_{i}=\omega_{0i}^{2}q_{0i}^{2}/E_{0}$, we can write (63) as a function
of initial modal coordinates as well, i.e.

$\gamma_{\mathrm{opt}}(\\{q_{0i}\\})=\sqrt{\frac{E_{0}}{2\sum_{i=1}^{N}q_{0i}^{2}}}\,,$
(71)

where $q_{0i}\in[-\sqrt{E_{0}}/\omega_{0i},\sqrt{E_{0}}/\omega_{0i}]$ and the
condition $0\leq\sum_{i=1}^{N}\omega_{0i}^{2}q_{0i}^{2}\leq E_{0}$ holds.

###  4.1 Quantitative example

Here we consider the MDOF system as the one shown schematically in Fig. 10 but
with $m_{i}=m$, $c_{i}=c$ for $i=\\{1,...,N\\}$, and with $k_{i}=k$ for
$i=\\{1,...,N+1\\}$. Such a system without damping, i.e with $c_{i}=0$
$\forall i$, is a standard part of the undergraduate physics/mechanics courses
[17]. Therefore, for the MDOF system with $N$ degrees of freedom we are
considering here, the undamped modal frequencies are [16, 17]

$\omega_{0i}=2\omega_{0}\sin\left(\frac{i\pi}{2(N+1)}\right)\,,\textrm{with}\,\,i=\\{1,...,N\\}\,,$
(72)

and where $\omega_{0}=\sqrt{k/m}$. In Fig. 11(a) we show undamped modal
frequencies $\omega_{01}$, $\omega_{0N}$ and damping coefficient
$\overline{\gamma}_{\mathrm{opt}}$, i.e. (68), calculated with (72), as
functions of $N$. We clearly see that the coefficient
$\overline{\gamma}_{\mathrm{opt}}$ is in the over-damped regime from the
perspective of the first mode, and in the under-damped regime from the
perspective of highest mode, for any $N>1$, and in the case $N=1$ all three
values match. In Fig. 11(b) we show ratios
$\overline{\gamma}_{\mathrm{opt}}/\omega_{01}$ and
$\omega_{0N}/\overline{\gamma}_{\mathrm{opt}}$ and we see that both ratios
increase with increasing $N$.

Figure 11: (a) Undamped modal frequencies $\omega_{01}$ (blue circles),
$\omega_{0N}$ (red x’s) and the damping coefficient
$\overline{\gamma}_{\mathrm{opt}}$ (black squares) as functions of the number
of the masses $N$. (b) Ratios $\overline{\gamma}_{\mathrm{opt}}/\omega_{01}$
(blue line) and $\omega_{0N}/\overline{\gamma}_{\mathrm{opt}}$ (red line),
shown as solid lines due to the high density of the shown points.

We show in B that the following limits hold

$\lim_{N\to+\infty}\overline{\gamma}_{\mathrm{opt}}(N)=0,$ (73)
$\lim_{N\to+\infty}\frac{\overline{\gamma}_{\mathrm{opt}}(N)}{\omega_{01}(N)}=+\infty,$
(74)
$\lim_{N\to+\infty}\frac{\omega_{0N}(N)}{\overline{\gamma}_{\mathrm{opt}}(N)}=+\infty.$
(75)

We note here that these limit values do not correspond to the transition from
a discrete to a continuous system, but simply tell us the behavior of these
quantities with respect to the increase in the number of masses, i.e. with
respect to the increase in the size of the discrete system.

From everything that has been said so far, it is clear that the damping
coefficient $\overline{\gamma}_{\mathrm{opt}}$, obtained by minimizing the
energy integral averaged over all initial conditions that correspond to the
same initial energy, cannot be considered generally as optimal and that, by
itself, it says nothing about optimal damping of the system whose dynamics
started with some specific initial condition. Damping coefficient (63), which
is given by the minimization of the energy integral for a specific initial
condition, is of course a better choice for optimal damping of an MDOF system,
than the damping coefficient $\overline{\gamma}_{\mathrm{opt}}$, if we want to
consider how the system dissipates energy the fastest for a particular initial
condition, but, as we argue in the section 4.2, this damping coefficient also
has some obvious deficiencies.

###  4.2 Issues with the minimum of the energy integral as a criterion for
optimal damping

We can ask, for example, whether in an experiment, with known initial
conditions, in which an MDOF system is excited to oscillate, a damping
coefficient (63) would be the best choice if we want that the system settles
down in equilibrium as soon as possible? Here, in three points, we explain why
we think the answer to that question is negative:

  * 1.

From relation (62), we see that, due to the term $\sin 2\theta_{i}$, the
energy integral is sensitive to changes $\theta_{i}\rightarrow-\theta_{i}$ and
$\theta_{i}\rightarrow\pi-\theta_{i}$, which correspond to changes of initial
conditions $(q_{0i},\dot{q}_{0i})\rightarrow(q_{0i},-\dot{q}_{0i})$ and
$(q_{0i},\dot{q}_{0i})\rightarrow(-q_{0i},\dot{q}_{0i})$. When we
differentiate (62) to determine $\gamma$ for which the energy integral has a
minimum, the term $\sin 2\theta_{i}$ cancels and as a result the coefficient
(63) is not sensitive to this change in initial conditions. Such changes in
the initial conditions lead to significantly different situations. For
example, if $q_{0i}>0$ and $\dot{q}_{0i}>0$, the $i$-th mode in the critical
and over-damped regime (i.e. for $\gamma\geq\omega_{0i}$) will never reach the
equilibrium position, while for $q_{0i}>0$ and $\dot{q}_{0i}<0$, and $i$-th
mode initial kinetic energy grater than initial potential energy, it can go
through the equilibrium position once, depending on the magnitude of the
damping coefficient, and there will be the smallest damping coefficient in the
over-damped regime for which no crossing occurs and for which the solution
converges to equilibrium faster than for any other damping coefficient [19].
Therefore, the damping coefficient considered optimal would have to be
sensitive to this change in initial conditions.

  * 2.

Damping coefficient (63) has singularities for $\cos\theta_{i}=0$ $\forall i$,
i.e. for initial conditions for which all initial energy is kinetic. For such
initial conditions, the higher the damping coefficient, the higher and faster
the dissipation. In other words, the higher the damping coefficient, the
faster the energy integral decreases. Therefore, coefficient (63) diverges for
that type of initial conditions. This would actually mean that, for this
initial conditions, it is optimal to take the damping coefficient as high as
possible, but in principle this corresponds to a situation in which all modes
are highly over-damped, i.e. all masses reach their maximum displacements in a
very short time and afterwards they begin to return to the equilibrium
position almost infinitely slowly. Figuratively speaking, it is as if we
immersed the system in concrete. This issue has recently been addressed in the
context of free vibrations of SDOF [19] and was already noticed in [20].
Therefore, simply taking the highest possible damping coefficient, as
suggested by relation (63) for this type of initial conditions, is not a good
option.

  * 3.

The damping coefficient (63) is determined on the basis of the energy integral
over the entire time and therefore it does not take into account that in
nature and experiments these systems effectively return to the equilibrium
state for some finite time.

Because of the above points, in the next section we provide a new approach to
determine the optimal damping of MDOF systems.

##  5 Optimal damping of an MDOF system: a new perspective

From a theoretical perspective, systems with viscous damping asymptotically
approach the equilibrium state and never reach it exactly. In nature and in
experiments, these systems reach the equilibrium state which is not an exact
zero energy state, but rather a state in which the energy of the system has
decreased to the level of the energy imparted to the system by the surrounding
noise, or to the energy resolution of the measuring apparatus. Following this
line of thought, we will define a system to be in equilibrium for times
$t>\tau$ such that

$\frac{E(\tau)}{E_{0}}=10^{-\delta}\,,$ (76)

where $E(\tau)$ is the energy of the system at $t=\tau$, $E_{0}$ is the
initial energy, and $\delta>0$ is a dimensionless parameter that defines what
fraction of the initial energy is left in the system. This line of thought has
recently been used to determine the optimal damping of SDOF systems [19], and
here we extend it to MDOF systems. Therefore, in what follows, we will
consider as optimal the damping coefficient for which the systems energy drops
to some energy level of interest, e.g. to the energy resolution of the
experiment, the fastest and we will denote it with $\tilde{\gamma}$.

###  5.1 Optimal damping of the $i$-th mode of a MDOF system with MPD

Here we will consider the behavior of the energy of the $i$-th mode of the
MDOF system with MPD and determine the optimal damping coefficient
$\tilde{\gamma}_{i}$ of the $i$-th mode with respect to criterion (76). For
any MDOF system with $N\geq 1$ degrees of freedom with MPD, each mode behaves
as a SDOF system studied in section 2, with the damping coefficient $\gamma$
and the undamped (natural) frequency $\omega_{0i}$. Thus (see relation (31)),
the ratio of the energy of the $i$-th mode, $E_{i}(\gamma,t)$, and initial
energy of the $i$-th mode, $E_{0i}$, is given by

$\frac{E_{i}(\gamma,t)}{E_{0i}}=e^{-2\gamma
t}\left(\cos^{2}(\omega_{i}t)+\gamma\cos
2\theta_{i}\frac{\sin(2\omega_{i}t)}{\omega_{i}}+\left(\omega_{0i}^{2}+\gamma^{2}+2\omega_{0i}\gamma\sin
2\theta_{i}\right)\frac{\sin^{2}(\omega_{i}t)}{\omega_{i}^{2}}\right)\,$ (77)

for the under-damped ($\gamma<\omega_{oi}$) and over-damped
($\gamma>\omega_{0i}$) regime. We will repeat here briefly for clarity,
$\omega_{i}=\sqrt{\omega_{0i}^{2}-\gamma^{2}}$ is the damped angular frequency
and $\theta_{i}$ is the polar angle which determines the initial conditions
$q_{0i}$ and $\dot{q}_{0i}$ of the $i$-th mode and the distribution of the
initial energy within the mode, i.e. initial potential and initial kinetic
energy of the $i$-th mode are $E_{0Pi}=E_{0i}\cos^{2}\theta_{i}$ and
$E_{0Ki}=E_{0i}\sin^{2}\theta_{i}$ respectively. Energy to initial energy
ratio for the $i$-th mode in the critically damped regime
($\gamma=\omega_{0i}$) is simply obtained by taking
$\gamma\rightarrow\omega_{0i}$ limit of the relation (77), and we obtain

$\frac{E_{i}(\gamma=\omega_{0i},t)}{E_{0i}}=e^{-2\omega_{0i}t}\left(1+2\omega_{0i}(\cos
2\theta_{i})t+2\omega_{0i}^{2}\left(1+\sin 2\theta_{i}\right)t^{2}\right)\,.$
(78)

In relations (77) and (78), we explicitly show that the energy depends on the
damping coefficient and time, because in what follows we will plot these
quantities as functions of these two variables for fixed initial conditions,
i.e. fixed $\theta_{i}$. We will investigate the behavior for several types of
initial conditions, which of course will not cover all possible types of
initial conditions, but will give us a sufficiently clear picture of the
determination and behavior of the optimal damping with respect to the initial
conditions and the equilibrium state defined with condition (76).

####  5.1.1 Initial energy of the $i$-th mode comprised only of potential
energy

In Fig. 12 we show the base $10$ logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition
$\theta_{i}=0$, which corresponds to the initial energy of the $i$-th mode
comprised only of potential energy. Four black contour lines denote points
with $E_{i}(\gamma,t)/E_{0i}=\\{10^{-3},10^{-4},10^{-5},10^{-6}\\}$
respectively, as indicated by the numbers placed to the left of each contour
line. Each contour line has a unique point closest to the $\gamma$ axis, i.e.
corresponding to the damping coefficient $\tilde{\gamma}_{i}$ for which that
energy level is reached the fastest. As an example, we draw arrow in Fig. 12
that points to the coordinates
$(\gamma,t)=(0.840\omega_{0i},5.15\omega_{0i}^{-1})$, i.e. to the tip of the
contour line with points corresponding to $E_{i}(\gamma,t)=10^{-4}E_{0i}$.
Thus, for the initial condition $\theta_{i}=0$,
$\tilde{\gamma}_{i}=0.840\omega_{0i}$ is the optimal damping coefficient for
the $i$-th mode to reach this energy level the fastest, and it does so at the
instant $\tau_{i}=5.15\omega_{0i}^{-1}$. In Table 1 we show results for other
energy levels corresponding to contour lines shown in Fig. 12. Here, and in
the rest of the paper, we have rounded the results for the damping coefficient
to three decimal places, and for the time to two decimal places.

Figure 12: The base $10$ logarithm of the ratio (77), i.e. $\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition $\theta_{i}=0$. For this initial condition, initial energy of the $i$-th mode is comprised only of potential energy. Four black contour lines denote points with $E_{i}(\gamma,t)/E_{0i}=\\{10^{-3},10^{-4},10^{-5},10^{-6}\\}$ respectively, as indicated by the numbers placed to the left of each contour line. As an example of determining the optimal damping for which the system reaches the desired energy level the fastest, i.e. with respect to the condition (76), we draw the arrow that points to the coordinates $(\gamma,t)=(0.840\omega_{0i},5.15\omega_{0i}^{-1})$ for which the $i$-th mode reaches the level $E_{i}(\gamma,t)/E_{0i}=10^{-4}$ the fastest. Thus, $\tilde{\gamma}_{i}=0.840\omega_{0i}$ is the optimal damping coefficient to reach this energy level the fastest. Optimal values for other energy levels, denoted with contour lines, are given in Table 1.  $E_{i}(\gamma,t)/E_{0i}$ | $\tilde{\gamma}_{i}\left[\omega_{0i}\right]$ | $\tau_{i}[\omega_{0i}^{-1}]$
---|---|---
$10^{-3}$ |  $0.769$ | $4.18$
$10^{-4}$ |  $0.840$ | $5.15$
$10^{-5}$ |  $0.885$ | $6.16$
$10^{-6}$ |  $0.915$ | $7.20$
Table 1: Optimal damping coefficient $\tilde{\gamma}_{i}$ for which the energy
of the $i$-th mode drops to the level $10^{-\delta}E_{0i}$ the fastest, with
the initial condition $\theta_{i}=0$.

Consider now, for example, a thought experiment in which we excite a MDOF
system so that it vibrates only in the first mode and that all initial energy
was potential, i.e. $E_{01}=E_{0}$ and $\theta_{1}=0$. Furthermore, suppose
that the system has effectively returned to equilibrium when its energy drops
below $10^{-6}E_{0}$, due to the resolution of the measuring apparatus. It is
clear form the Table 1 that $\tilde{\gamma}_{1}=0.915\omega_{01}$ would be
optimal in such a scenario. In the same scenario, optimal damping coefficient
given by the minimization of the energy integral, i.e. (63), would be
$\gamma_{\mathrm{opt}}=\sqrt{2}\omega_{01}/2=0.707\omega_{01}$, thus, a very
bad choice in the sense that this damping coefficient would not be optimal
even in an experiment with a significantly poorer energy resolution (see Table
1). This simple example illustrates that, from a practical point of view, one
has to take into account both the initial conditions and the resolution of the
measuring apparatus in order to determine the optimal damping coefficient.

####  5.1.2 Initial energy of the $i$-th mode comprised only of kinetic energy

In Fig. 13(a) and (b) we show the base $10$ logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition
$\theta_{i}=\pi/2$, which corresponds to the initial energy of the $i$-th mode
comprised only of kinetic energy. In Fig. 13(b) we show results for larger
data span than in Fig. 13(a), and only contour line for points corresponding
to $E_{i}(\gamma,t)=10^{-3}E_{0i}$. The left arrow in Fig. 13(b) indicates the
same coordinates as the arrow in Fig. 13(a), and the right arrow in Fig. 13(b)
points to the coordinates
$(\gamma,t)=(13.316\omega_{0i},4.66\omega_{0i}^{-1})$ with
$E_{i}(\gamma,t)=10^{-3}E_{0i}$. Thus, for $\gamma>13.316\omega_{0i}$ the
system comes sooner to the energy level $10^{-3}E_{0i}$ than for
$\gamma=0.722\omega_{0i}$, but these highly over-damped damping coefficients
would correspond to restricting the system to infinitesimal displacements from
equilibrium, after which the system returns to the equilibrium state
practically infinitely slowly [19]. Thus, for this initial condition we take
the damping coefficient in the under-damped regime, i.e.
$\tilde{\gamma}_{i}=0.722\omega_{0i}$, as optimal for reaching the level
$E_{i}(\gamma,t)=10^{-3}E_{0i}$ the fastest. For all energy levels the
behaviour is qualitatively the same, and the results are given in Table 2.

Figure 13: The base $10$ logarithm of the ratio (77), i.e. $\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition $\theta_{i}=\pi/2$. For this initial condition, initial energy of the $i$-th mode is comprised only of kinetic energy. (a) Four black contour lines denote points with $E_{i}(\gamma,t)/E_{0i}=\\{10^{-3},10^{-4},10^{-5},10^{-6}\\}$ respectively, and the arrow points to the coordinates $(\gamma,t)=(0.722\omega_{0i},4.66\omega_{0i}^{-1})$, with $E_{i}(\gamma,t)/E_{0i}=10^{-3}$, for which this level of energy is reached in shortest time for the shown data span. (b) Contour line for points with $E_{i}(\gamma,t)=10^{-3}E_{0i}$ is shown for larger data span, left arrow points to the coordinates $(\gamma,t)=(0.722\omega_{0i},4.66\omega_{0i}^{-1})$, and the right arrow to the coordinates $(\gamma,t)=(13.316\omega_{0i},4.66\omega_{0i}^{-1})$, both with $E_{i}(\gamma,t)/E_{0i}=10^{-3}$. Thus, for $\gamma>13.316\omega_{0i}$ energy level $10^{-3}E_{0i}$ is reached faster than for $\gamma=0.7223\omega_{0i}$. See text for details. $E_{i}(\gamma,t)/E_{0i}$ | $\tilde{\gamma}_{i}\left[\omega_{0i}\right]$ | $\tau_{i}[\omega_{0i}^{-1}]$
---|---|---
$10^{-3}$ |  $0.722$ | $4.66$
$10^{-4}$ |  $0.794$ | $5.50$
$10^{-5}$ |  $0.852$ | $6.42$
$10^{-6}$ |  $0.892$ | $7.40$
Table 2: Optimal damping coefficient $\tilde{\gamma}_{i}$ for which the energy
of the $i$-th mode drops to the level $10^{-\delta}E_{0i}$ the fastest, with
the initial condition $\theta_{i}=\pi/2$.

Consider now, for example, a thought experiment in which we excite a MDOF
system so that it vibrates only in the first mode and that all initial energy
was kinetic, i.e. $E_{01}=E_{0}$ and $\theta_{1}=\pi/2$. Furthermore, suppose
that the system has effectively returned to equilibrium when its energy drops
below $10^{-6}E_{0}$, due to the resolution of the measuring apparatus. It is
clear form the Table 2 that $\tilde{\gamma}_{1}=0.892\omega_{01}$ would be
optimal in such a scenario. In the same scenario, optimal damping coefficient
given by the minimization of the energy integral, i.e. (63), would be
$\gamma_{\mathrm{opt}}=+\infty$.

Here we note that if in such an experiment we can set the damping coefficient
to be in the over-damped regime in the first part of the motion, i.e. when the
system is moving from the equilibrium position to the maximum displacement,
and in the under-damped regime in the second part of the motion, i.e. when the
system moves from the position of maximum displacement back towards the
equilibrium position, then the fastest way to achieve equilibrium would be to
take the largest experimentally available over-damped coefficient in the first
part of the motion, and the under-damped coefficient optimised like in 5.1.1
in the second part of the motion, with the fact that we have to carry out the
optimization with respect to the energy left in the system at the moment when
the system reached the maximum displacement and with respect to the energy
resolution of the experiment.

####  5.1.3 Initial energy of the $i$-th mode comprised of potential and
kinetic energy

In Fig. 14(a) we show the base $10$ logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition
$\theta_{i}=\pi/3$, which corresponds to the initial energy of the $i$-th mode
comprised of kinetic energy $E_{0Ki}=3E_{0i}/4$ and potential energy
$E_{0Pi}=E_{0i}/4$, with both initial normal coordinate and velocity positive,
i.e. with $q_{0i}>0$ and $\dot{q}_{0i}>0$. The results for optimal damping are
obtained by the same method as in 5.1.1 and are given in Table 4 for data
shown in Fig. 14(a), and in Table 4 for data shown in Fig. 14(b). We see that
the energy dissipation strongly depends on the relative sign between $q_{0i}$
and $\dot{q}_{0i}$. It was recently shown, for free vibrations of SDOF, that
for an initial condition with initial kinetic energy greater than initial
potential energy and opposite signs between $x_{0}$ and $v_{0}$, an optimal
damping coefficient can be found in the over-damped regime [19], thus, the
same is true when we consider any mode of a MDOF system with MPD.

Figure 14: The base $10$ logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, (a) for initial condition
$\theta_{i}=\pi/3$, and (b) for initial condition $\theta_{i}=-\pi/3$. For
both initial conditions, initial energy of the $i$-th mode is comprised of
kinetic energy $E_{0Ki}=3E_{0i}/4$ and potential energy $E_{0Pi}=E_{0i}/4$.
For $\theta_{i}=\pi/3$ initial normal coordinate and velocity are of the same
signs, i.e. $q_{0i}>0$ and $\dot{q}_{0i}>0$. For $\theta_{i}=-\pi/3$ initial
normal coordinate and velocity are of the opposite signs, i.e. $q_{0i}>0$ and
$\dot{q}_{0i}<0$.

$E_{i}(\gamma,t)/E_{0i}$ | $\tilde{\gamma}_{i}\left[\omega_{0i}\right]$ | $\tau_{i}[\omega_{0i}^{-1}]$
---|---|---
$10^{-3}$ |  $0.751$ | $4.66$
$10^{-4}$ |  $0.825$ | $5.58$
$10^{-5}$ |  $0.875$ | $6.55$
$10^{-6}$ |  $0.908$ | $7.58$
Table 3: Optimal damping coefficient $\tilde{\gamma}_{i}$ for which the energy
of the $i$-th mode drops to the level $10^{-\delta}E_{0i}$ the fastest, with
the initial condition $\theta_{i}=\pi/3$.

$E_{i}(\gamma,t)/E_{0i}$ | $\tilde{\gamma}_{i}\left[\omega_{0i}\right]$ | $\tau_{i}[\omega_{0i}^{-1}]$
---|---|---
$10^{-3}$ |  $1.075$ | $1.87$
$10^{-4}$ |  $1.112$ | $2.42$
$10^{-5}$ |  $1.135$ | $3.02$
$10^{-6}$ |  $1.145$ | $3.64$
Table 4: Optimal damping coefficient $\tilde{\gamma}_{i}$ for which the energy
of the $i$-th mode drops to the level $10^{-\delta}E_{0i}$ the fastest, with
the initial condition $\theta_{i}=-\pi/3$.

Consider now, for example, a thought experiment in which we excite a MDOF
system so that it vibrates only in the first mode and that $75\%$ of initial
energy was kinetic and $25\%$ of initial energy was potential, and with
$q_{01}>0$ and $\dot{q}_{01}>0$, i.e. $E_{01}=E_{0}$ and $\theta_{1}=\pi/3$.
Furthermore, suppose that the system has effectively returned to equilibrium
when its energy drops below $10^{-6}E_{0}$, due to the resolution of the
measuring apparatus. It is clear form the Table 4 that
$\tilde{\gamma}_{1}=0.908\omega_{01}$ would be optimal in such a scenario. In
the same scenario, but with $q_{01}>0$ and $\dot{q}_{01}<0$, i.e. for
$\theta_{1}=-\pi/3$, we see from Table 4 that
$\tilde{\gamma}_{1}=1.145\omega_{01}$ would be optimal. Optimal damping
coefficient given by the minimization of the energy integral, i.e. (63), is
insensitive to the change of the sign of $\dot{q}_{01}$, and it would be
$\gamma_{\mathrm{opt}}=\sqrt{2}\omega_{01}=1.414\omega_{01}$ in both cases.

We note here, that for the initial conditions of the $i$-th mode with initial
kinetic energy much grater than initial potential energy, i.e.
$E_{0Ki}>>E_{0Pi}$, and with opposite signs of initial displacement and
velocity, i.e.
$\mathop{\mathrm{sgn}}(q_{0i})\neq\mathop{\mathrm{sgn}}(\dot{q}_{0i})$, the
optimal damping coefficient is going to be deep in the over-damped regime and
dissipation of initial energy will happen in a very short time. If, for any
reason, this is not desirable in some particular application, one can always
find damping coefficient in the under-damped regime, with that same initial
condition, which can serve as an alternative. As an example of such a
situation, in Fig. 15 we show the base 10 logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition
$\theta_{i}=-9\pi/20$, which corresponds to the initial energy of the $i$-th
mode comprised of kinetic energy $E_{0Ki}\approx 0.98E_{0i}$ and potential
energy $E_{0Pi}\approx 0.02E_{0i}$, with $q_{0i}>0$ and $\dot{q}_{0i}<0$. In
Fig. 15 we see that the $i$-th mode will reach the energy level
$10^{-6}E_{0i}$ the fastest for $\gamma=3.222\omega_{0i}$, and in case, e.g.,
that such damping coefficient is difficult to achieve experimentally, another
choice for the optimal damping coefficient can be $\gamma=0.883\omega_{0i}$.

Figure 15: The base $10$ logarithm of the ratio (77), i.e.
$\log\left(E_{i}(\gamma,t)/E_{0i}\right)$, for initial condition
$\theta_{i}=-9\pi/20$. Black contour line denotes the points with
$E_{i}(\gamma,t)=10^{-6}E_{0i}$. Left arrow points to the coordinates
$(\gamma,t)=(0.883\omega_{0i},7.30\omega_{0i}^{-1})$ for which level
$10^{-6}E_{0i}$ is reached the fastest in the under-damped regime, and the
right arrow points to the coordinates
$(\gamma,t)=(3.222\omega_{0i},0.87\omega_{0i}^{-1})$ for which the same level
is reached the fastest in the over-damped regime.

###  5.2 Optimal damping of a MDOF system with MPD

If all modes of a MDOF system with $N$ defrees of freedom are excited, the
ratio of the energy of the system, $E(\gamma,t)$, and initial energy of the
system, $E_{0}$, is given by

$\frac{E(\gamma,t)}{E_{0}}=\sum_{i=1}^{N}\frac{E_{0i}}{E_{0}}e^{-2\gamma
t}\left(\cos^{2}(\omega_{i}t)+\gamma\cos
2\theta_{i}\frac{\sin(2\omega_{i}t)}{\omega_{i}}+\left(\omega_{0i}^{2}+\gamma^{2}+2\omega_{0i}\gamma\sin
2\theta_{i}\right)\frac{\sin^{2}(\omega_{i}t)}{\omega_{i}^{2}}\right)\,,$ (79)

where the set of all initial energies of the modes, i.e. $\\{E_{0i}\\}$, and
the set of all polar angles, i.e. $\\{\theta_{i}\\}$, determines the initial
condition of the whole system. Since for MPD the damping of the system as a
whole is determined by only one damping coefficient $\gamma$, we can calculate
the base $10$ logarithm of the ratio (79), but using a unique units for
$\gamma$, $t$ and $\omega_{0i}$ for all modes, and from these data determine
the optimal damping coefficient $\tilde{\gamma}$, for which the system will
come to equilibrium in the sense of the condition (76) the fastest, in the
same way as in subsubsections 5.1.1-5.1.3 where we showed how to determine the
optimal damping of individual modes. One practical choice for the units might
be $\omega_{01}$ for $\gamma$ and for $\omega_{0i}$ $\forall i$, and
$\omega_{01}^{-1}$ for $t$. This way, we have the easiest insight into the
relationship between the first mode and the optimal damping coefficient that
we want to determine, in the sense that we can easily see whether the first
mode is under-damped, over-damped or critically damped in relation to it,
which is important since the first mode is often the dominant mode. If we
apply this to the 2-DOF system studied in 3.1, and we obtain

$\frac{E(\gamma,t)}{E_{0}}=\sum_{i=1}^{2}\frac{E_{0i}}{E_{0}}e^{-2\gamma
t}\left(\cos^{2}(\omega_{i}t)+\gamma\cos
2\theta_{i}\frac{\sin(2\omega_{i}t)}{\omega_{i}}+\left(\omega_{0i}^{2}+\gamma^{2}+2\omega_{0i}\gamma\sin
2\theta_{i}\right)\frac{\sin^{2}(\omega_{i}t)}{\omega_{i}^{2}}\right)\,,$ (80)

where $\omega_{01}=\omega_{0}$, $\omega_{02}=\sqrt{3}\omega_{0}$,
$\omega_{1}=\sqrt{\omega_{0}^{2}-\gamma^{2}}$,
$\omega_{2}=\sqrt{3\omega_{0}^{2}-\gamma^{2}}$ and we take that the damping
coefficient is in $\omega_{0}$ units, while the time is in $\omega_{0}^{-1}$
units. We are now in a position to determine the optimal damping of this 2-DOF
system for different initial conditions. Again, we will not investigate all
possible types of the initial conditions, but two qualitatively different
ones, one with initial energy comprised only of potential energy, and the
other with initial energy comprised only of kinetic energy. These two examples
will give us a picture of the procedure for determining the optimal damping
coefficient $\tilde{\gamma}$ for this 2-DOF system. The same procedure for
determining the optimal damping can be in principle carried out for any MDOF
system with MPD, with any initial condition.

####  5.2.1 Optimal damping of the 2-DOF system with initial energy comprised
only of potential energy

Here we choose initial condition with $E_{01}=E_{02}=E_{0}/2$ and
$\theta_{1}=\theta_{2}=0$, i.e. with initial potential energy distributed
equally between the two modes and zero initial kinetic energy. In Fig. 16 we
show the base $10$ logarithm of the ratio (80), i.e.
$\log\left(E(\gamma,t)/E_{0}\right)$, for the chosen initial condition. In
Table 5 we show results for other energy levels corresponding to contour lines
shown in Fig. 16. For this initial condition, optimal damping coefficient
given by the minimization of the energy integral, i.e. (63), is
$\gamma_{\mathrm{opt}}=\sqrt{3/4}\omega_{0}=0.866\omega_{0}$.

Figure 16: The base $10$ logarithm of the ratio (80), i.e. $\log\left(E(\gamma,t)/E_{0}\right)$, for initial condition $E_{01}=E_{02}=E_{0}/2$ and $\theta_{1}=\theta_{2}=0$. For this initial condition, initial energy of the 2-DOF system is comprised only of potential energy distributed equally between the modes. Four black contour lines denote points with $E(\gamma,t)/E_{0}=\\{10^{-3},10^{-4},10^{-5},10^{-6}\\}$ respectively, as indicated by the numbers placed to the left of each contour line. As an example of determining the optimal damping for which the system reaches the desired energy level the fastest, i.e. with respect to the condition (76), we draw the arrow that points to the coordinates $(\gamma,t)=(0.859\omega_{0},5.37\omega_{0}^{-1})$ for which the energy of the system reaches the level $E(\gamma,t)/E_{0}=10^{-4}$ the fastest. Thus, $\tilde{\gamma}=0.859\omega_{0}$ is the optimal damping coefficient to reach this energy level the fastest.  $E(\gamma,t)/E_{0}$ | $\tilde{\gamma}\left[\omega_{0}\right]$ | $\tau[\omega_{0}^{-1}]$
---|---|---
$10^{-3}$ |  $0.817$ | $4.36$
$10^{-4}$ |  $0.859$ | $5.37$
$10^{-5}$ |  $0.893$ | $6.27$
$10^{-6}$ |  $0.924$ | $7.55$
Table 5: Optimal damping coefficient $\tilde{\gamma}$ for which the energy of
the system drops to the level $10^{-\delta}E_{0}$ the fastest, with the
initial condition $E_{01}=E_{02}=E_{0}/2$ and $\theta_{1}=\theta_{2}=0$.

####  5.2.2 Optimal damping of the 2-DOF system with initial energy comprised
only of kinetic energy

Here we choose initial condition with $E_{01}=E_{02}=E_{0}/2$ and
$\theta_{1}=\theta_{2}=\pi/2$, i.e. with initial kinetic energy distributed
equally between the two modes and zero initial potential energy. In Fig. 17(a)
and (b) we show the base $10$ logarithm of the ratio (80), i.e.
$\log\left(E(\gamma,t)/E_{0}\right)$, for the chosen initial condition. In
Table 6 we show results for other energy levels corresponding to contour lines
shown in Fig. 17(a). For this initial condition, optimal damping coefficient
given by the minimization of the energy integral, i.e. (63), is
$\gamma_{\mathrm{opt}}=+\infty$.

Figure 17: The base $10$ logarithm of the ratio (80), i.e. $\log\left(E(\gamma,t)/E_{0}\right)$, for initial condition $E_{01}=E_{02}=E_{0}/2$ and $\theta_{1}=\theta_{2}=\pi/2$. For this initial condition, initial energy of the 2-DOF system is comprised only of kinetic energy distributed equally between the modes. (a) Four black contour lines denote points with $E(\gamma,t)/E_{0}=\\{10^{-3},10^{-4},10^{-5},10^{-6}\\}$ respectively, and the arrow points to the coordinates $(\gamma,t)=(0.783\omega_{0},4.60\omega_{0}^{-1})$, with $E(\gamma,t)/E_{0}=10^{-3}$, for which this level of energy is reached in shortest time for the shown data span. (b) Contour line for points with $E(\gamma,t)=10^{-3}E_{0}$ is shown for larger data span, left arrow points to the coordinates $(\gamma,t)=(0.783\omega_{0},4.60\omega_{0}^{-1})$, and the right arrow to the coordinates $(\gamma,t)=(15.927\omega_{0},4.60\omega_{0}^{-1})$, both with $E_{i}(\gamma,t)/E_{0i}=10^{-3}$. Thus, for $\gamma>15.927\omega_{0}$ energy level $10^{-3}E_{0}$ is reached faster than for $\gamma=0.783\omega_{0}$. $E(\gamma,t)/E_{0}$ | $\tilde{\gamma}\left[\omega_{0}\right]$ | $\tau[\omega_{0}^{-1}]$
---|---|---
$10^{-3}$ |  $0.783$ | $4.60$
$10^{-4}$ |  $0.838$ | $5.72$
$10^{-5}$ |  $0.861$ | $6.47$
$10^{-6}$ |  $0.909$ | $7.78$
Table 6: Optimal damping coefficient $\tilde{\gamma}$ for which the energy of
the system drops to the level $10^{-\delta}E_{0}$ the fastest, with the
initial condition $E_{01}=E_{02}=E_{0}/2$ and $\theta_{1}=\theta_{2}=\pi/2$.

##  6 Conclusion

The main message of this paper is that the dissipation of the initial energy
in vibrating systems strongly depends on the initial conditions with which the
dynamics of the system started, and ideally it would be optimal to always
adjust the damping to the initial conditions. We show that having some
information about the type of initial conditions can also help in devising the
optimal damping. Although the paper is devoted to the case of mass
proportional damping, the same phenomena also appears for other types of
damping. The new method we propose for optimizing the damping is general and
we expect that it will be applicable also to the types of vibrating systems we
did not study in this paper.

##  7 Acknowledgments

We are grateful to Bojan Lončar for making schematic figures of 2-DOF and MDOF
systems, i.e. Fig. 6 and 10, according to our sketches. This work was
supported by the QuantiXLie Center of Excellence, a project co-financed by the
Croatian Government and European Union through the European Regional
Development Fund, the Competitiveness and Cohesion Operational Programme
(Grant No. KK.01.1.1.01.0004).

The authors have no conflicts to disclose.

## References

  * [1] R. A. Rojas, E. Wehrle, R. Vidoni, Optimal design for the passive control of vibration based on limit cycles, Shock and Vibration 2019 (2019).
  * [2] Z.-F. Fu, J. He, Modal analysis, Elsevier, 2001.
  * [3] I. Takewaki, Optimal damper placement for minimum transfer functions, Earthquake Engineering & Structural Dynamics 26 (11) (1997) 1113–1124.
  * [4] Y. Cheung, W. Wong, ${H}_{\infty}$ and ${H}_{2}$ optimizations of a dynamic vibration absorber for suppressing vibrations in plates, Journal of Sound and Vibration 320 (1) (2009) 29–42.
  * [5] N. Alujević, G. Zhao, B. Depraetere, P. Sas, B. Pluymers, W. Desmet, ${H}_{2}$ optimal vibration control using inertial actuators and a comparison with tuned mass dampers, Journal of Sound and Vibration 333 (18) (2014) 4073–4083.
  * [6] P. Lancaster, U. Prells, Inverse problems for damped vibrating systems, Journal of Sound and Vibration 283 (3) (2005) 891–914.
  * [7] R. T. Haftka, Z. Gürdal, Elements of structural optimization, Vol. 11, Springer Science & Business Media, 2012.
  * [8] K. Veselić, On linear vibrational systems with one dimensional damping II, Integral Equations and Operator Theory 13 (1990) 883–897.
  * [9] S. J. Cox, I. Nakić, A. Rittmann, K. Veselić, Lyapunov optimization of a damped system, Systems & control letters 53 (3-4) (2004) 187–194.
  * [10] I. Nakic, Optimal damping of vibrational systems, PhD thesis, Fernuniversität Hagen (2003).
  * [11] N. Truhar, K. Veselić, An efficient method for estimating the optimal dampers’ viscosity for linear vibrating systems using Lyapunov equation, SIAM journal on matrix analysis and applications 31 (1) (2009) 18–39.
  * [12] K. Veselić, Damped oscillations of linear systems: a mathematical introduction, Vol. 2023, Springer Science & Business Media, 2011.
  * [13] J. T. Katsikadelis, Dynamic analysis of structures, Academic press, 2020.
  * [14] M. Gutierrez Soto, H. Adeli, Tuned mass dampers, Archives of Computational Methods in Engineering 20 (2013) 419–431.
  * [15] L. Zuo, S. Nayfeh, Minimax optimization of multi-degree-of-freedom tuned-mass dampers, Journal of Sound and Vibration 272 (3) (2004) 893–908.
  * [16] M. D. Grigoriu, Linear Dynamical Systems, Springer, 2021.
  * [17] F. S. Crawford, Waves : Berkeley physics course / volume 3, McGraw-Hill, New York, USA, 1968.
  * [18] A. Muravyov, S. G. Hutton, A geometrical interpretation of damping for discrete classically damped systems, International journal of mechanical sciences 40 (1) (1998) 83–86.
  * [19] K. Lelas, N. Poljak, D. Jukić, Damped harmonic oscillator revisited: the fastest route to equilibrium (2023).  arXiv:2301.08222.
  * [20] I. Nakić, Z. Tomljanović, N. Truhar, Mixed control of vibrational systems, ZAMM-Journal of Applied Mathematics and Mechanics/Zeitschrift für Angewandte Mathematik und Mechanik 99 (9) (2019) e201800328.

##  Appendix A Average of the integral (62) over a set of all initial
conditions

For reader’s convenience, we will repeat the integral (62) here

$I(\gamma,\\{a_{i}\\},\\{\theta_{i}\\})=E_{0}\sum_{i=1}^{N}\frac{a_{i}^{2}}{2\omega_{0i}}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}}+\frac{\gamma}{\omega_{0i}}\cos
2\theta_{i}+\sin 2\theta_{i}\right)\,.$ (81)

In order to calculate the average of (81) over a set of all initial
conditions, one has to integrate (81) over all coefficients $a_{i}$, which
satisfy $\sum_{i=1}^{N}a_{i}^{2}=1$ and $a_{i}\in[-1,1]$, and over all angles
$\theta_{i}\in[0,2\pi]$. Due to $\int_{0}^{2\pi}\cos
2\theta_{i}d\theta_{i}=\int_{0}^{2\pi}\sin 2\theta_{i}d\theta_{i}=0$, terms
with sine and cosine functions don’t contribute to the average of (81).
Integration over all possible coefficients $a_{i}$ amounts to calculating the
average of $a_{i}^{2}$ over a sphere of radius one embedded in $N$ dimensional
space. If we were to calculate the average of the equation of a sphere
$\sum_{i=1}^{N}a_{i}^{2}=1$ over a sphere defined by that equation, we would
get

$\sum_{i=1}^{N}\overline{a_{i}^{2}}=1\,,$ (82)

where $\overline{a_{i}^{2}}$ denotes the average of $a_{i}^{2}$ over a sphere.
Due to the symmetry of the sphere and the fact that we are integrating over
the whole sphere, contribution of each $\overline{a_{i}^{2}}$ in the sum (82)
has to be the same, so we can easily conclude that

$\overline{a_{i}^{2}}=\frac{1}{N}\,,$ (83)

for any $i$. Thus, the average of (81) over all possible initial conditions is

$\overline{I}(\gamma)=\frac{E_{0}}{2N}\sum_{i=1}^{N}\left(\frac{\omega_{0i}^{2}+\gamma^{2}}{\gamma\omega_{0i}^{2}}\right)\,.$
(84)

##  Appendix B Limit values (73), (74) and (75)

For reader’s convenience, we repeat here (68) and (72)

$\overline{\gamma}_{\mathrm{opt}}=N^{1/2}\left(\sum_{i=1}^{N}\frac{1}{\omega_{0i}^{2}}\right)^{-1/2}\,$
(85)
$\omega_{0i}=2\omega_{0}\sin\left(\frac{i\pi}{2(N+1)}\right)\,,\textrm{with}\,\,i=\\{1,...,N\\}\,.$
(86)

Using (86), we can write (85) as

$\overline{\gamma}_{\mathrm{opt}}=2\omega_{0}N^{1/2}\left(\sum_{i=1}^{N}\frac{1}{\sin^{2}\zeta_{i}}\right)^{-1/2}\,,$
(87)

where $\zeta_{i}=\frac{i\pi}{2(N+1)}$. Using the fact that $\sin x<x$ for
$0<x<\pi/2$, we obtain

$\overline{\gamma}_{\mathrm{opt}}<2\omega_{0}N^{1/2}\left(\frac{4(N+1)^{2}}{\pi^{2}}\sum_{i=1}^{N}\frac{1}{i^{2}}\right)^{-1/2}=2\omega_{0}N^{1/2}\frac{\pi}{2(N+1)}\left(\sum_{i=1}^{N}\frac{1}{i^{2}}\right)^{-1/2}.$

Now taking $N\to\infty$ and using the well-known formula
$\sum_{i=1}^{\infty}\frac{1}{i^{2}}=\frac{\pi^{2}}{6}$, we obtain (73).

Now we focus on the limit (74). We will use the following well-known
inequality $\sin x>x/2$ for $0<x<\pi/2$ (this can be easily seen by, e.g.
using the fact that $\sin$ is a concave function on $[0,\pi/2]$). From (87) it
follows

$\frac{\overline{\gamma}_{\mathrm{opt}}}{\omega_{01}}=N^{1/2}(\sin\zeta_{1})^{-1}\left(\sum_{i=1}^{N}\frac{1}{\sin^{2}\zeta_{i}}\right)^{-1/2}>N^{1/2}\zeta_{1}^{-1}\cdot\frac{1}{2}\left(\sum_{i=1}^{N}\frac{1}{\zeta_{i}^{2}}\right)^{-1/2}=\frac{1}{2}N^{1/2}\left(\sum_{i=1}^{N}\frac{1}{i^{2}}\right)^{-1/2},$

hence we obtain (74).

The limit (75) is also easy to prove. Since

$\lim_{N\to+\infty}\omega_{0N}=\lim_{N\to+\infty}2\omega_{0}\sin\left(\frac{N\pi}{2(N+1)}\right)=2\omega_{0}\,$
(88)

and we already showed (73), it is easy to conclude that

$\lim_{N\to+\infty}\frac{\omega_{0N}}{\overline{\gamma}_{\mathrm{opt}}}=+\infty\,,$
(89)

i.e. the limit (75) holds.