Work done on a harmonically vibrating deformable body, with damping

Question

Explanation

The following equation may be obtained by applying Finite Element Analysis to a deformable body considering linear elasticity and a harmonic input force f(t) = Fe^jωt. We assume that the steady-state response of the displacementes will be u(t) = Ue^jωt. Here italic letters are vectors, bold face letters are matrices, and normal letters are scalars.

(K + jωC - ω²M)U = F       (1)

where K, C and M are the stiffness, damping and mass matrices, respectively. All matrices are symmetric.

Now, when no damping is considered, the amplitude of the work of the external forces done on the system can be calculated simply as:

W_{no damping} = (1/2)U^T(K - ω²M)U   ∈ ℝ       (2)

The question

The question then is: what is the amplitude of the work done when there is damping?

Although I am using FE notation here, answers using continuum mechanics notation are welcome! (or any other, for that matter)

My attempt

The main issue here, I believe, is that while in a conservative system (no damping), the work of the external forces always oscilates around 0, in a non-conservative system, there is a net linear increase in work done, which is equal to the energy dissipated by the damping terms. But I do not know how to calculate this.

I tried explicitly calculating W = ∫ f(t)^Tu(t) dt for t ∈ [t₀ , t₁], but I don't know if I should actually take the conjugate transpose of f(t) instead of just the transpose, nor if that equation is even correct, since it yields a complex value because U ∈ ℂ^{n x 1}. The integral yielded:

W_{with dampnig} = (1/2)U^T (K + jωC - ω²M)U (e^2jωt₁ - e^2jωt₀)       (3)

I am having difficulties finding references for this in FEA using complex values.

---

A few tests, based on suggested answers

Suppose $f(t)=Fe^{j\omega t}$ and $u(t)=Ue^{j\omega t}$, with $F, U \in \mathbb{C}^{N \times 1}$. Then $\frac{du}{dt}(t)=j\omega Ue^{j\omega t}$.

We have that $\left(\mathbf{K} + j\omega\mathbf{C} - \omega^2\mathbf{M}\right)U = F$

Case A:

$$\text{W}=\int_0^{\frac{2\pi}{\omega}}f(t)^T\frac{du}{dt}(t)\,dt = \int_0^{\frac{2\pi}{\omega}} j\omega F^TUe^{2j\omega t} \,dt = F^TU\frac{j\omega}{2j\omega}\left.e^{2j\omega t}\right|_0^{\frac{2\pi}{\omega}}=$$ $$=\frac{F^TU}{2}\left(e^{4\pi j}-e^0\right)=0$$

Case B:

$$\text{W}=\int_0^{\frac{2\pi}{\omega}}f(t)^H\frac{du}{dt}(t)\,dt = \int_0^{\frac{2\pi}{\omega}} j\omega F^HUe^{0} \,dt = 2\pi jF^HU=2\pi jU^H\left(\mathbf{K} - j\omega\mathbf{C} - \omega^2\mathbf{M}\right)U$$ $$\text{W} = -2\pi\omega U^H\mathbf{C}U+2\pi jU^H\left(\mathbf{K}-\omega^2\mathbf{M}\right)U$$ Since all matrices are real and symmetric, the vector, matrix, vector products all result in real values.

Case C:

$$\text{W}=\int_0^{\frac{2\pi}{\omega}}\Re\left(f(t)\right)^T\Re\left(\frac{du}{dt}(t)\right)\,dt$$ Now we separate, $U=U_R+jU_I$, with $U_R, U_I\in\mathbb{R}^{N\times1}$. And we also have to explicitly write $e^{j\omega t}$ using Euler's formula, so that: $$u(t)=Ue^{j\omega t}=U_R\cos{\omega t}-U_I\sin{\omega t} +j\left( U_R\sin{\omega t} + U_I\cos{\omega t} \right)$$ $$\dfrac{du}{dt}(t)=-\omega\left( U_R\sin{\omega t} + U_I\cos{\omega t} \right)+j\omega\left(U_R\cos{\omega t}-U_I\sin{\omega t}\right)$$ After a lot of algebraic manuplation, one can find that the terms multiplying $\mathbf{K}$ and $\mathbf{M}$ vanish after the integration. We then obtain: $$\text{W}=\int_0^{\frac{2\pi}{\omega}}\omega^2\left(U_R^T\mathbf{C}U_R\sin^2{\omega t}+U_I^T\mathbf{C}U_I\cos^2{\omega t}+U_R^T\mathbf{C}U_I\sin{2\omega t}\right)\,dt$$ Knowing that $\int_0^{\frac{2\pi}{\omega}}\sin^2{\omega t}\,dt=\int_0^{\frac{2\pi}{\omega}}\cos^2{\omega t}\,dt=\frac{\pi}{\omega}$ and $\int_0^{\frac{2\pi}{\omega}}\sin{2\omega t}\,dt=0$: $$\text{W}=\pi\omega\left(U_R^T\mathbf{C}U_R+U_I^T\mathbf{C}U_I\right)=\pi\omega U^H\mathbf{C}U$$

Answer 1

Writing the equations for steady state forced response using complex variables is very convenient mathematically, but interpreting the results is not always "trivial", as you have discovered.

The (physically real and time-dependent) displacements and forces are actually the real parts of the complex expressions for $f(t)$ and $U(t)$.

The complex expressions are not completely determined, since they do not specify the absolute phase angle of any of the quantities at time $t=0$. They do specify the relative phase angles of the components of force, displacement, velocity, acceleration, etc, but everything can be multiplied by an arbitrary constant $e^{i\phi}$ representing a phase shift.

Since the choice of when to set $t = 0$ is arbitrary, only the relative phases have any physical meaning, so the fact that the absolute phases are not determined by the complex equation of motion doesn't matter.

If you use the real parts and calculate the work $$\text{W} = \int \Re(f(t))^T \Re(u(t))\,dt$$ over one cycle of the harmonic motion, things should become clearer. Once you have the "correct answer" from the above, you can write it more simply without using the real parts of $f$ and $U$.

Of course if there is no damping term in the equation of motion, then integrated over a complete cycle the result will always be $\text{W} = 0$.

Another way to think about this is to eliminate the arbitrary phase shift by arbitrarily fixing one component of either $F$ or $u$ to be a real value. This is equivalent to choosing a fixed value of $\phi$ for the arbitrary constant $e^{i\phi}$, or (equivalently) choosing the instant when $t = 0$.

Answer 2

1. Going back to definitions in 1D : $W=\int F(u).du$ where $du$ is the infinitesimal displacement along application of force $F$. You may also write $W= \int F.v dt$ with $v=du/dt$, that is $\underline{v}=j \omega\underline{u}$ in harmonic regime using complex notations. You should stick to these definitions. For what you want to do I'd go with the second definition.
2. Since ${F}=\underline{Z}{v}$, with $Z$ complex $v$ is out of phase by $\phi= -arg (\underline{Z})$ with the excitation force. The instantaneous power has, similarly to electrical circuits, a $\cos^{2}(\omega t)$ term with a net contribution over a period (the active power loss due to damping) and a $\cos(\omega t)\sin(\omega t)$ with zero net contribution over a period (elastic energy restored).The real and imaginary part of $\underline {F}\underline{v^{\star}}$ will respectively give you the correct active and reactive power amplitude.