# Dynamic Stability of Orthotropic Viscoelastic Rectangular Plate of an Arbitrarily Varying Thickness
## 2 Methods
### Computational Methods  
The solution of the system (1)-(3) is obtained by using the Bubnov-Galerkin method based on a polynomial approximation of the deflection and displacements; discretization is performed in spatial variables at each moment and is reduced to solving systems of nonlinear integrodifferential equations with weakly singular kernels with variable coefficients. The improved numerical method [48] based on quadrature formulas was applied to solve the obtained nondecaying system with a weakly singular Koltunov-Rzhanitsyn kernel [40; 41; 49]. A relaxation kernel \(\Gamma(t)=Ae^{-\beta t}t^{\alpha-1},\ (0<\alpha<1)\) with three rheological parameters was first proposed in [40]. Publications [40; 41] explain a method for determining the parameters \(A\), \(\alpha\), and \(\beta\). In [41], these parameters are tabulated. The parameters of specific materials are determined by combining the experimental curve and the theoretical curve plotted in logarithmic coordinates. The number of relaxation kernels for isotropic and orthotropic plates is different:  
1. The mechanical parameters of the viscoelastic isotropic are the same in all directions. Therefore, the viscoelastic properties of the plate material are described with one core with three rheological parameters. 2. Mechanical parameters of viscoelastic orthotropic plates are different in two mutually perpendicular directions, coinciding with the directions of the coordinate axes. In this case, there are five different relaxation nuclei. In addition, there are three rheological parameters in each relaxation core. Thus, there are 15 different rheological parameters in total. After applying the Bubnov-Galerkin method, most of the dynamic problems of viscoelastic thin-walled structures are reduced to solving nondecaying systems of integrodifferential equations of the following form:  
\[ \begin{array}{l}\sum\limits_{n=1}^{N}\sum\limits_{m=1}^{M}\left(a_{klmn}\ddot {u}_{nm}+\lambda_{nmkl}^{2}u_{nm}\right)=X_{nm}(t,u_{11},\dots,u_{NM},v_{11}, \dots,v_{NM},w_{11},\dots,w_{NM},\\ \int\limits_{0}^{t}\Phi_{1 nm}(t,\tau,u_{11}\dots,u_{NM},v_{11},\dots,v_{NM},w_ {11},\dots,w_{NM})\,d\tau),\\ \sum\limits_{n=1}^{N}\sum\limits_{m=1}^{M}\left(b_{klmn}\ddot{v}_{nm}+\chi_{ nmkl}^{2}v_{nm}\right)=Y_{nm}(t,u_{11},\dots,u_{NM},v_{11},\dots,v_{NM},w_ {11},\dots,w_{NM},\\ \int\limits_{0}^{t}\Phi_{1 nm}(t,\tau,u_{11}\dots,u_{NM},v_{11},\dots,v_{NM},w_ {11},\dots,w_{NM})\,d\tau),\\ \sum\limits_{n=1}^{N}\sum\limits_{m=1}^{M}\left(c_{klmn}\ddot{w}_{nm}+\omega_{ nmkl}^{2}w_{nm}\right)=Z_{nm}(t,u_{11},\dots,u_{NM},v_{11},\dots,v_{NM},w_ {11},\dots,w_{NM},\\ \int\limits_{0}^{t}\Psi_{1 nm}(t,\tau,u_{11}\dots,u_{NM},v_{11},\dots,v_{NM},w_ {11},\dots,w_{NM})\,d\tau)\\ u_{nm}(0)=u_{0 nm},\dot{u}_{nm}(0)=\dot{u}_{0 nm},\dot{v}_{nm}(0)=v_{0 nm},\dot{v }_{nm}(0)=\dot{v}_{0 nm},w_{nm}(0)=w_{0 nm},\dot{w}_{nm}(0)=\dot{w}_{0 nm},\\ n=1,2,\dots,N;m=1,2,\dots,M,\end{array}  \tag{4}\]  
where \(u_{nm}=u_{nm}(t)\), \(v_{nm}=v_{nm}(t)\), and \(w_{nm}=w_{nm}(t)\) are unknown functions of time; \(X_{nm},Y_{nm},Z_{nm},\Phi_{1 nm},\Phi_{1 nm}\) are continuous functions in the field of changing arguments; and \(a_{klmn}\), \(b_{klmn}\), \(c_{klmn}\), \(\lambda^{2}\),\(\lambda_{klmn}\), \(\lambda^{2}\),\(\lambda_{klmn}\), and \(\omega^{2}\),\(\lambda_{klmn}\) are given constant numbers. Integration of system (4) twice on \(t\) leads to its integral form. Then, \(t=t_{i}\), \(t_{i}=i\Delta t\), \(i=1,2,\dots\) (\(\Delta t=\text{const}\) is the interpolation step) is assumed.