## 2 Linear Elastic Traction-Separation Constitutive Behavior of the Interface Model  
In this approach, a zero-thickness cohesive interface element was introduced into the finite element computations for composite materials. The zero-thickness means that the geometric thickness of cohesive elements is equal to zero, and not the constitutive thickness. Since the initial geometric thickness of the interface element is zero, the deformation state of the zero-thickness element cannot be described by the classical definition of strain. Instead, the constitutive response of the zero-thickness interface element is established in terms of a linear elastic traction-separation description. In general, the available traction-separation model consists of three parts: linear elastic traction-separation; damage initiation criteria; and damage evolution laws [34; 35]. However, only initially linear elastic behavior is considered due to stress distribution being the only concern in this paper.  
The constitutive response of this interface was obtained using a procedure introduced by [PERSON] and [PERSON] [36]. Consider a typical interface that uses a traction-separation description between two plies as shown in Figure 1. The two coincident nodes of the interface element are a pair of nodes, and each interface element has four pairs of nodes, and each node has three degrees of freedom (DOFs); thus, there are 24 DOFs in total for each interface element. The nodes global displacement vector \(u_{N}\) is defined as follows:  
\[u_{N}=\left\{u_{ix}^{+},u_{iy^{\prime}}^{+},u_{iz}^{-},u_{ix^{\prime}}^{-},u_{ iz}^{-}\right\}^{\mathrm{T}} \tag{1}\]  
where \(u_{ix^{\prime}}^{+},u_{iy^{\prime}}^{+},u_{iz}^{-},u_{iz^{\prime}}^{-},u_{iz}^ {-}\) are the displacements of the \(i\) nodes of the interface element in the \(x,y\), and \(z\) directions.  
The relative displacement \(\Delta u\) for each pair of nodes can be described as  
\[\Delta u=u_{k}^{+}-u_{k}^{-}=(-I_{12\times 12}\big{|}I_{12\times 12})u_{N} \tag{2}\]  
where \(I_{12\times 12}\) is the identity diagonal matrix.  
Converting \(\Delta u\) to the relative displacement at the interface local coordinate, defining the relative displacements tensor \(\delta_{k}=\{\delta_{n},\delta_{s},\delta_{t}\}^{\mathrm{T}}\) in local coordinates. Subscript \(n\), \(s\) and \(t\) represent the local 3-direction, the 1- and 2-direction, respectively, in which \(\delta_{k}\) can be written as  
\[\delta_{k}=B\Delta u_{k} \tag{3}\]  
where \(B\) is the geometric matrix.  
Unlike continuum elements, the interface element stiffness matrix before softening onset requires the penalty stiffness \(K\) of the interface material [36]. The penalty stiffness \(K\) relates to the element traction stresses \(\sigma_{k}=\{\sigma_{n},\sigma_{s},\sigma_{t}\}^{\mathrm{T}}\) to the element relative displacements tensor \(\delta_{k}\):  
Figure 1: Zero-thickness cohesive interface element.  
\[\sigma_{k}=K\delta_{k}. \tag{4}\]  
From Equation (4), it can be seen that traction-separation constitutive behavior is not stress-strain constitutive behavior. The interface element can be visualized as a spring between the initially coincident nodes of the element. After softening onset, the initially coincident nodes will open or slide relative to each other [37]. Due to the zero-thickness interface element describing interlaminar behavior, the three traction stresses can be considered as corresponding interlaminar stresses. In addition, the penalty stiffness \(K\) is the interface element constitutive parameter, which has been introduced by some researchers [38]. The value of \(K\) is usually expected to be higher than the elastic modulus of the material, but it should also be noted that a larger value for \(K\) in FEM might lead to a bad convergence [31].