## 2 Theory and Formulation
### Problem Description and Fundamental Equations  
The structural model relies on the following constitutive relations [35]  
\[\{\sigma\}=[Q]\{\epsilon\} \tag{7}\]  
where \(Q_{ij}=Q_{ij}(x)\), \(i,j=1,\dots,6\) stands for the stiffness properties, and \(\{\sigma\}^{T}=\left\{\sigma_{xx}\;\sigma_{yy}\;\sigma_{zz}\;\sigma_{xy}\; \sigma_{xz}\;\sigma_{yz}\right\}^{T}\), \(\{\epsilon\}^{T}=\left\{\epsilon_{xx}\;\epsilon_{yy}\;\epsilon_{zz}\; \epsilon_{xy}\;\epsilon_{xz}\;\epsilon_{yz}\right\}^{T}\) refer to the stress and strain vector, respectively. For orthotropic structures, the constitutive relations read as follows [36]  
\[\left\{ \begin{array}{c}\sigma_{xx}\\ \sigma_{yy}\\ \sigma_{zz}\\ \sigma_{xy}\\ \sigma_{xz}\\ \sigma_{yz}\end{array} \right\}=\left[ \begin{array}{cccccc}Q_{11}&Q_{12}&Q_{ 13}&0&0&0\\ Q_{21}&Q_{22}&Q_{23}&0&0&0\\ Q_{31}&Q_{32}&Q_{33}&0&0&0\\ 0&0&0&Q_{44}&0&0\\ 0&0&0&Q_{55}&0\\ 0&0&0&0&0&Q_{66}\end{array} \right]\left\{ \begin{array}{c}\epsilon_{xx}\\ \epsilon_{yy}\\ \epsilon_{zz}\\ \epsilon_{xy}\\ \epsilon_{xz}\\ \epsilon_{yz}\end{array} \right\} \tag{8}\]Based on a classic Euler-Bernoulli beam theory, the displacement field for a CNTs-strengthened composite beam is defined as  
\[u_{1}=\eta(x,t)-z\frac{\partial\xi(x,t)}{\partial x} \tag{9}\]  
\[u_{2}=\xi(x,t) \tag{10}\]  
accounting for the displacement components of the mid-plane of the beam in the length and thickness directions, \(\eta\) and \(\xi\), respectively. The nonzero strain component of the composite beam reads as follows  
\[\varepsilon_{xx}=\frac{\partial\eta(x,t)}{\partial x}-z\frac{\partial^{2}\xi (x,t)}{\partial x^{2}} \tag{11}\]  
Hereafter, we consider a beam embedded on an elastic medium, as depicted in Figure 1, accounting for a Pasternak coefficient, \(K_{\text{P}}\), as a shear layer and a Winkler coefficient, \(K_{\text{W}}\), as normal layer. The elastic foundation is generally defined by using the stiffness coefficient and power-law distribution function as [37]. \[K_{\text{P}}(x)=K_{\text{P0}}\Big{(}1-\beta\Big{(}\frac{x}{L}\Big{)}^{m}\Big{)} \tag{12}\]  
\[K_{\text{W}}(x)=K_{\text{W0}}\Big{(}1-\alpha\Big{(}\frac{x}{L}\Big{)}^{n}\Big{)} \tag{13}\]  
depending on the power term \(n\) and coefficient \(\alpha\). Moreover, \(K_{\text{W0}}\) stands for the initial value of the Winkler coefficient, whereas the Pasternak coefficient is likewise demonstrated using an initial value, power term, and coefficient of \(K_{\text{P0}}\), \(m\), and \(\beta\), respectively. The equations of motion for the CNTs-reinforced structure are determined based on the Hamiltonian principle, i.e.,  
\[\int\limits_{t_{0}}^{t_{1}}\delta(\Phi+W_{nc}-T)dt=0 \tag{14}\]  
where \(\Phi\) refers to the strain energy, \(T\) is the kinetic energy, \(W_{nc}\) is the work done by non-conservative external forces, and \(t\) is the time variable. Each energy quantity can be defined in variational form as  
\[\delta\,\Phi=\frac{1}{2}\!\int_{A}\int\limits_{0}^{L}\sigma_{xx}\delta \varepsilon_{xx}\text{d}x\text{d}A=\frac{1}{2}\!\int_{A}\int\limits_{0}^{L} \sigma_{xx}\bigg{(}\delta\frac{\partial\eta(x,t)}{\partial x}-z\delta\frac{ \partial^{2}\xi(x,t)}{\partial x^{2}}\bigg{)}\text{d}x\text{d}A \tag{15}\]  
\[\delta T=\frac{1}{2}\!\int_{0}^{L}\int_{A}\rho(x)\delta\bigg{(}\bigg{(}\frac{ \partial\eta(x,t)}{\partial t}-z\frac{\partial^{2}\xi(x,t)}{\partial x \partial t}\bigg{)}^{2}+\bigg{(}\frac{\partial\xi(x,t)}{\partial t}\bigg{)}^{ 2}\bigg{)}\text{d}A\text{d}x \tag{16}\]  
Figure 1: Geometrical model of composite beam on a Winkler-Pasternak foundation.