## 2 Materials and Methods
### Numerical Simulation by FEM  
#### 2.3.1 Elasto-Dynamic Equilibrium Equations and Approximation with Finite Elements  
For this purpose, we will consider a viscoelastic, homogeneous, and isotropic linear body from the composite material that has the values of the elastic constants given by experimental realization and which for the time \(t>0\) is subjected to some external forces \(b(t,x)\). At the moment \(t=0\), the body is considered to occupy a domain \(\Omega\subset\mathbb{R}^{d}\) bordered by a smooth border \(\Gamma\), with d = 2 or 3. The border \(\Gamma\) is divided into two open and disjoint parts \(\Gamma_{D}\) and \(\Gamma_{N}\), so that \(\Gamma=\overline{\Gamma}_{D}\cup\overline{\Gamma}_{N}\), on which the boundary conditions will be imposed, and for \(t=0\), the initial conditions will be imposed. The dynamic equilibrium equation of linear viscoelasticity is of the following form:  
\[div\sigma(\dot{u}(t,x),u(t,x))+b(t,x)=\rho\ddot{u}(t,x), \tag{7}\]  
where \(\dot{u}(t,x)=\partial u(t,x)/\partial t\) and \(\ddot{u}(t,x)=\partial^{2}u(t,x)/\partial t^{2}\), \(t\in[0,T],T>0,x\in\Omega\), and the deformation tensor ([PERSON]'s principle) is as follows:  
\[\varepsilon(u(t,x))=(\
abla_{x}u(t,x)+(\
abla_{x}u(t,x))^{T})/2, \tag{8}\]  
The constitutive law of linear elasticity (Hooke's law) for viscoelastic material is the following:  
\[\sigma(\dot{u}(t,x),u(t,x))=A^{(1)}\varepsilon(\dot{u}(t,x))+A^{(2)}\varepsilon (u(t,x)). \tag{9}\]  
Boundary conditions:  
\[u(t,x)=0,\ x\in\Gamma_{D}, \tag{10}\]  
the homogeneous Dirichlet boundary condition and  
\[u(t,x)=f_{N}(x),\ x\in\Gamma_{N}, \tag{11}\]  
the Neumann boundary condition. The initial conditions are as follows:  
\[u(0,x)=\dot{u}(0,x)=0,\ \ x\in\Omega, \tag{12}\]  
where \(\sigma\) is the stress, \(u\) is the deformation, \(\
abla_{x}(\cdot)=grad(\cdot)\) is the gradient, and \(A^{(1)}\) and \(A^{(2)}\) are the fourth-order tensors of viscoelasticity which contains [PERSON]'s material constants. A first step towards approximation with the finite element method of the problem with boundary and initial values (7)-(12) is the transformation of the equation from the strong form to the weak form, a procedure based on the principle of virtual work, achieved by multiplying Equation (7) by admissible kinematic test functions \(\delta v\), also through member-by-member integration:  
\[\int\limits_{\Omega}div\sigma(\dot{u}(t,x),u(t,x))\delta vdV+\int\limits_{ \Omega}b(t,x)\delta vdV=\int\limits_{\Omega}\rho\ddot{u}(t,x)\delta vdV. \tag{13}\]  
Using the Green-Ostrogradski theorem for the first term in Equation (13), which allows us to include the boundary conditions in the weak form and to reduce the order of the differential operator to half, Equation (13) will become the following:  
\[\int\limits_{\Omega}\sigma(\dot{u}(t,x),u(t,x))\
abla_{x}(\delta v)dV-\int \limits_{\Omega}b(t,x)\delta vdV+\int\limits_{\Omega}\rho\ddot{u}(t,x)\delta vdV -\int\limits_{\Gamma_{N}}f_{N}(t,x)\delta vdS=0. \tag{14}\]The finite element formulation of the equilibrium equation for elastic bodies [19; 20; 25; 26] is standard, and we obtain the second-order differential system with constant coefficients with \(m\) equations, and after assembling the contributions of each finite element, we obtain the following:  
\[M\ddot{U}(t)+KU(t)=F(t) \tag{15}\]  
where with \(KU(t)\) was denoted the first term of Equation (14), with \(M\ddot{U}(t)\) the third term from (14), and with \(F(t)\), the second plus the fourth term from (14). The initial conditions will be as follows:  
\[U(0)=\dot{U}(0)=0 \tag{16}\]  
The size of the system is given by the total number of degrees of freedom of the partition with finite elements \(m\), the \(K\) matrix is the stiffness matrix, symmetric and positive definite, of size (\(mxm\)), \(M\) is the mass matrix, symmetric and semipositive definite, \(F(t)\) is the time-varying excitation force, and \(U(t)\) is the column vector of the nodal displacements.