## 2 Basic Theory
### Thin-Layer Element
The strain component can be simplified to \(\varepsilon=[\varepsilon_{z}\ \gamma_{yz}\ \gamma_{zx}]^{T}\). Synthesizing the above analysis, it is assumed that the normal contact characteristics and tangential contact characteristics of the interface are independent of each other. The constitutive equation of the thin layer element which characterizes interface contact is:  
\[\left\{ \begin{array}{c}\sigma_{n}\\ \tau_{tx}\\ \tau_{ty}\end{array} \right\}=\left[ \begin{array}{ccc}E_{n}&0&0\\ 0&G_{t}&0\\ 0&0&G_{t}\end{array} \right]\left\{ \begin{array}{c}\varepsilon_{n}\\ \gamma_{tx}\\ \gamma_{ty}\end{array} \right\} \tag{5}\]  
\(E_{n}\) and \(G_{t}\) are the normal elastic modulus and tangential shear modulus, respectively. If the tangential and normal contact properties are coupled, the coupling term can be added to Equation (5). If the contact properties are coupled to normals and tangents, the coupling term can be added to the constitutive relation (5). The constitutive equation is used when using isotropic constitutive material to simulate the thin layer element. \[\left\{ \begin{array}{c}\sigma_{xx}\\ \sigma_{yy}\\ \sigma_{zz}\\ \tau_{xy}\\ \tau_{yz}\\ \tau_{zx}\end{array} \right\}=\left[ \begin{array}{ccccc}\lambda+2G&\lambda& \lambda\\ &\lambda+2G&\lambda\\ &&\lambda+2G&\\ &&&G&\\ &&&G&\\ sym&&&G&\\ &&&G\end{array} \right]\left\{ \begin{array}{c}\varepsilon_{xx}\\ \varepsilon_{yy}\\ \varepsilon_{zz}\\ \gamma_{xy}\\ \gamma_{yz}\\ \gamma_{zx}\end{array} \right\} \tag{6}\]  
\(\lambda\) is the Lame constant and \(G=E/2(1+\
u)\) is the shear modulus. As we know, the number of independent material parameters of the isotropic material is 2. When the  
\begin{table}
\begin{tabular}{c c c c c c c} \hline \multirow{2}{*}{**Gaussian**} & \multicolumn{3}{c}{\(\mathbf{I_{1}\times l_{2}\times e=1\times 1\times 0.1}\)} & \multicolumn{3}{c}{\(\mathbf{I_{1}\times l_{2}\times e=5\times 5\times 0.1}\)} \\ \cline{2-7}  & \(\partial\mathbf{N_{i}/\partial x}\) & \(\partial\mathbf{N_{i}/\partial y}\)  
thickness of the cell is much smaller than the feature size of the other two directions, let \(\varepsilon_{x}=\varepsilon_{y}=\varepsilon_{xy}\approx 0\). Finally, the constitutive equation of the material can be written as follows:  
\[\left\{ \begin{array}{c}\sigma_{z}\\ \tau_{yz}\\ \tau_{zx}\end{array} \right\}=\left[ \begin{array}{ccc}\lambda+2G&0&0\\ 0&G&0\\ 0&0&G\end{array} \right]\left\{ \begin{array}{c}\varepsilon_{z}\\ \gamma_{yz}\\ \gamma_{zx}\end{array} \right\} \tag{8}\]  
In this equation, the normal elastic constant and the tangent elastic constant are not independent. Combined with the constitutive relation from the aforementioned content, the element stiffness matrix of \(K\) can be obtained according to the principle of virtual work. \[\delta\mathbf{W}=\int_{0}^{l_{1}}\!\int_{0}^{l_{2}}\!\int_{0}^{d}\!\sigma^{T}\delta \varepsilon\mathrm{d}x\mathrm{d}y\mathrm{d}z=\delta\mathbf{u}_{nodal}^{T}\mathbf{K} \mathbf{u}_{nodal} \tag{9}\]  
Isoparametric transformation is used to calculate the stiffness matrix \(K\) of the thin layer element. \[K=\int_{V_{\varepsilon}}\mathbf{B}^{T}\mathbf{D}\mathbf{B}\mathrm{d}V=\int_{-1}^{1}\!\int_ {-1}^{1}\!\mathbf{B}^{T}\mathbf{D}\mathbf{B}det(\mathbf{I})\mathrm{d}\xi\mathrm{d}\eta\mathrm{ d}\zeta \tag{10}\]  
Figure 3 shows the isoparametric transformation of the thin layer element. \(\mathbf{\xi}\), \(\mathbf{\eta}\), and \(\mathbf{\zeta}\) are global coordinate symbols, and \(\mathbf{J}\) is the Jacobian matrix. When the global coordinate system is consistent with the local coordinate system, the following expression is given:  
\[\mathbf{J}=\left[ \begin{array}{ccc}\partial x/\partial\xi&\partial y/\partial\xi &\partial z/\partial\xi\\ \partial x/\partial\eta&\partial y/\partial\eta&\partial z/\partial\eta\\ \partial x/\partial\xi&\partial y/\partial\xi&\partial z/\partial\zeta\end{array}  \right]=\left[ \begin{array}{ccc}l_{1}/2&0&0\\ 0&l_{2}/2&0\\ 0&0&d/2\end{array}