Mathematical and Numerical Calculation of the Interlayer Slip of a Two-Layer Glued Beam

This study analyses the spatial variation of the deflection, rotation, and slip at the interface of two-layer bamboo scrimber-concrete composite beam simply supported under uniform transverse load on its entire length. The Timoshenko model is considered, and the equations are solved using the analytical methods and direct numerical simulation based on finite differences schemes and a Fortran code. The results obtained show that it is a good agreement between the results from the mathematical calculation and those obtained from the numerical simulation.


Introduction
Composite structures have excellent combined advantages that are increasingly important in different areas of engineering practice, but more widely used in civil engineering. e classic cases of such structures in civil engineering are composite steel-concrete beams, composite wood-concrete panels, reinforced concrete in buildings and bridges, woodconcrete floor systems, coupled bracing walls, concrete beams reinforced externally with laminates, sandwich beams, and many others [1][2][3][4][5][6][7][8]. Some advantages of this type of composite structure constituted, for example, by two layers can be seen as easy manufacture, taking advantage of the material properties. In fact, there are different grafting methods to make two-layer beams or plates either using glue or connector beads to link the two parts at the contact surface. Because of these advantages, many researchers have focused on the mechanical analysis of these structures. Existing analysis indicates that the interface connection property has a significant impact on displacements and strains in discrete layered structures. In engineering practice, epoxy adhesives are sometimes used to connect several structural components to form a composite structure. For viscoelastic materials, the dependence of the stress state on the stress history is intrinsic.
is study focuses on the mechanical behavior of two-layer beams bonded with adhesives. Adhesives have inherently viscoelastic characteristics, and thus, the mechanical properties of laminated beams vary over time and even under static loads. However, those assemblies which include preexisting slip planes or precracked interfaces may experience the interlayer slip due to either the flexibility of adhesive. erefore, the incorporation of the interlayer slip into beam theories becomes necessary in order to provide an increased capacity to simulate the degradation of the stiffness of the composite structure [9][10][11][12][13][14][15].
As all mechanical researchers know, the Euler-Bernoulli (EB) beam hypothesis assumes that the cross-sections of the beam remain plane and perpendicular to the neutral axis during deformation. Based on this assumption, several studies can be found for layered beams glued with the viscoelastic interlayer. e first theory of partial action, but certainly the most cited, was developed by Newmark [16]. So far, a number of elastic theories have been developed and presented in the professional literature [17][18][19][20]. Zhang and Wang [21] analyzed the viscoelastic property of concrete beams reinforced with fiber-reinforced polymers and gave the analytical solutions using the Laplace numerical transformation method. Lei et al. [22] studied the dynamic characteristics of nonlocal viscoelastic damped EB beams with different boundary conditions using the transfer function method. Wang et al. [23] presented the analytical solutions for the behavior in bending of viscoelastic beams under quasistatic loading. e bending and buckling behaviors of laminated beams glued by the viscoelastic interlayer have been studied by Galuppi and Royer-Carfagni [24,25]. In their study, two types of analytical models for the interlayer were proposed: the entirely viscoelastic model considering the memory effect of the strains and the quasielastic model ignoring the memory effect of the strains (the viscoelastic material is approximately modeled as a "aging" material). Li et al. [26] investigated the time varying behaviors of the functionally graded beams in layers with a viscoelastic interlayer using the analytical Laplace transformation.
e main drawback of the classical Euler-Bernoulli beam theory is that no transverse shear strain is allowed. is implies infinite shear strength of each layer. As there is in reality no material having such a property, the adequacy of the theory of the beams of Euler-Bernoulli to the composite beams with the intermediate slip can be questioned, in particular for the thick and short composite beams. In these cases, the application of the theory of Timoshenko beams is essential for an accurate prediction of the mechanical behavior of the aforementioned structures. A large number of homogeneous beam elements based on Timoshenko's theory of beams have been developed in the literature. e transverse shear strains in the laminated beams with the intermediate slip were taken into account for the first time by Murakami [27,28], where the nonlinear interface slip laws were taken into account. e approximate equations of the beams are derived from the principle of the virtual work. In this 21st century, several recent publications in structural mechanics have looked at the theories of the nonlinear behavior of two-layer and multilayer beams glued or in partial interaction of the Bernoulli and Timoshenko model with the interlayer slip by different methods, whether numerical or analytical [29][30][31][32][33][34][35][36][37][38]. But until now, the analytical method of the weighted residuals and the finite difference for numerical simulation of composite beams with two layers glued using deformable epoxy adhesive with a sliding at the interface have not been developed to the knowledge of authors. Most of the time, finite elements' codes are used by engineers as they have already been fixed in calculations codes.
e inspiration for this study comes from the work of Araujo published on 18 April 2019, who studied two-layer Timoshenko beams bonded with a viscoelastic intermediate layer. Pei et al. [39] emphasized the importance of applying an analytical stress and strain solution to two-layer Timoshenko beams glued with the viscoelastic interlayer under a uniform transverse load using the Laplace transformation method. erefore, an appropriate approximation of the inplane displacements over the thickness of the plate and the incorporation of transverse shear stresses that govern the sliding behavior at the interface are considered to be key elements for the development of accurate theories of beams which can simulate the degradation of the rigidity due to the interlaying slip.
In this present study, a new analytical method of the twolayer composite beam with two layers with an interlayer slip which takes into account deformations, slips, rotations, and displacements has been developed using the approximate analytical method of the weighted residuals known as the least squares method [40]. e governing equations and the coherent boundary conditions are then deduced by making appropriate displacement and transverse strain assumptions. e advantage of using the weighted residual method is that it automatically produces the appropriate shear correction factors for the constitutive transverse shear relationships. Comparisons with the approximate analytical solutions for the laws of deformations, linear slip, and rotation and with the numerical results of finite difference simulations [41] indicate that the current theory accurately simulates the important in-plane responses of composite beams. e organization of this study is as follows. Following this introduction (Section 1), the geometric and theoretical formulation of the static bending of two-layer composite beams based is presented in Section 2. Section 3 presents the mathematical and numerical methods. In Section 4, we consider some examples, in which one of the layers is made of bamboo while the other one is the concrete. Some general information on these materials will be given before the presentation of the main results showing the deflection, the rotation angle, and particularly the slip. We conclude the study in Section 5.

Geometric Structure and Assumptions.
e geometric formulation of the problem is a two-layer composite beam simply supported with a length L and width b. It is comprised of two isotropic materials that have bonded together to the contact surface using an epoxy adhesive (DP3501 Gray), thin and relatively soft, assuming its thickness to be negligible, subjected to a uniform transverse load f z (x) applied along the length of the beam span, as shown in Figures 1 and 2. Two parts have different sizes in Cartesian (oxz) coordinates established with the origin of the lower beam surface with a thickness h a for the upper layer and h b for the lower layer. e elastic moduli of two outer layers are E a and E b , and the Poisson ratios are, respectively, μ a and μ b . e shear modulus of the epoxy interlayer depends on its characteristics and is denoted k sc [24]. To avoid repeating similar equations for the lower layer and the upper layer, we use the index a for the upper layer, the index b for the lower layer, and the sc index for the contact surface.
Before establishing the equations described, such a mechanical structure, we first indicate the assumptions considered in the study. e material in each layer is linear, elastic, and isotropic; the displacements are transverse; the slip and the rotation of the composite beam are small; there is no delamination phenomenon between two layers; the mass of the epoxy is much smaller than the mass of the beam (so the epoxy mass is neglected). Figures 1 and 2, the longitudinal displacement u(x, z) and the transverse displacement w(x) of the normal beam section along the direction of x can be written as follows:

Governing Equations. From
where u a (x) and u b (x) are, respectively, the displacements on the geometric axes of the upper and lower layers; θ(x) is the bent angle of rotation of the cross-section of the composite beam; w(x) represents the deflection of the composite beam; z a and z b represent the distances of the geometric axes of the upper and lower layers which take the point at the level of the lower surface of the beam; z represents the distance on the geometric axes of the upper and lower layers, indicating the passage points of the longitudinal displacements; h a and h b are, respectively, the upper and lower layer thicknesses of the composite beam. e interlayer slip g(x) is defined as the difference in axial displacements and the rotation angle of layers along the boundary surface, as follows: where H is the half-height sum of the upper and lower layers and H � (1/2)(h a + h b ). e interlayer shear force τ 0 can be written in the next form: where k sc is the shear modulus of epoxy interface which represents the stiffness of the connection of the upper and lower layers (the unit k sc is expressed in force by square length). e longitudinal normal strain ε xz and the transverse shear strain c xz in the upper and lower beam layers are given by Using the constitutive equations of elasticity and Hooke's law, we obtain the σ xz longitudinal stress and the transverse shear stress τ xz of the upper and lower beam layers as follows: Figure 2: Undeformed and deformed configuration of two-layer composite beam.
(E a , G a , k a , A a , I a ) Figure 1: Free body diagram of a nondeformable two-layer composite beam.
where G i is the shear modulus for isotropic materials and k i is the shear factor for rectangular cross-section of the beam component (i � a, b). e expressions of the forces (N i ), moments (M i ), and internal shear forces (T i ) of the component beam free body are defined as follows: where I i is the polar geometric moment with respect to the center of gravity and A i is the area of the upper and lower layers. By combining all the above expressions, the total energy of the composite beam can be written as where A is the area of the composite beam. By summing, we obtain knowing that w is only a function of x, the Lagrangian of the structure is e Lagrangian principles consider the following equations: By applying this formalism on the Lagrangian in equation (10), one obtains that the longitudinal deformations, the rotation angle, and the transverse deflection are governed by the following set of differential equations: Considering the simple supported ends at x � 0 and x � L, the boundary conditions can be written as follows: e set of equation (16) will be solved analytically using the weighted residuals method and numerically using the finite differences scheme implemented in a Fortran code. ese methods are presented in the next section.

e Weighted Residuals Approximation Method.
Analytically solving the differential equations obtained in (12)-(15), we will apply the weighted residual method of order four. e weighted residual methods seek to minimize the expression of the error in the differential equations or partial differential equations (and not on the unknown function itself ). For this purpose, approximate solutions are proposed with several unknown, whose number is reduced after applying the boundary conditions. en, one inserts the final expressions in each equations (12)- (15), each of which is called the residual R(x, a i ). e weighted residuals methods command that where domain corresponds to the physical domain on which the study is carried out. In this study, domain is the entire length of the beam. By the above minimization criteria, one obtains a set of algebraic equations whose solution gives the unknown coefficients.
Typically, in this study, we define u a (x), u b (x), w(x), and θ(x) as follows: where the coefficients a i , c i , and d i are unknown. By taking into account the boundary conditions posed in (17), we will have where the new unknown coefficients are a i , b i , c i , and d i . Substituting the above expressions (19) in the residuals, one obtains a set of 16 linear algebraic equations with the following unknown a i , b i , c i , and d i . e set of equations can thus be solved manually or using MATLAB to get the unknown values and consequently the mathematical expressions of the longitudinal displacement, rotation angle, and deflection. en, the slip can be calculated using equation (2).

e Finite Difference Method and the Numerical
Algorithm. To solve the set of differential equations (13) numerically, we have used the finite differences method with a second-order error. For any variable P, the following discretization schemes are used: where dx is the spatial step. Apart from the fact that we want to have a second order error, the discretization formula for the first-order derivative is also linked to the boundary conditions which are defined with the derivatives. Applying the discretization schemes, one arrives at the following set of coupled discrete equations: e set of discrete equations is then solved numerically using the iterative Gauss-Seidel method.

e Materials Used.
For the examples of applications, we have been interested by a composite beam made of a bamboo scrimber [42] and the concrete.
Our interest is based on bamboo. Indeed, bamboo scrimber has gained particular interest of researchers and manufacturers on account of its excellent mechanical characteristics and design. Bamboo scrimber also has the advantage of processing convenience. e manufacturer can utilize the wood-processing equipment and technology to process the bamboo scrimber. e structure of furniture made from bamboo scrimber could combine with each other by the traditional furniture tendon, and the modern connection structure could also be adopted. In addition to high mechanical performance, bamboo scrimber also has an excellent dimensional stability. e modulus of elasticity, short-beam strength, tensile strength, and compressive strength of bamboo also show overwhelmingly higher than those of other bamboo-based composites. e concrete is made of binders, sand, gravel, and mixing water. It is widely used in the field of construction in civil engineering, its compressive strength is very high, and its implementation is easy.

Results.
We consider a simply supported bamboo scrimber-concrete composite beam with a span length of 4 m loaded transversely (Figure 3 ). e bamboo scrimber crosssection is rectangular (300 mm width and 150 mm height). e concrete is 300 mm wide and 50 mm thick. e material parameters used for calculation are as follows: E a � 12 · 10 9 Pa, E b � 14 · 10 9 Pa, and μ a � μ b � 0.3. e shear stiffness k sc is taken equal to 2.43 · 10 6 Pa. A load of 1 kN/m is applied on the structure. Figures 4-7 present the spatial variations of the deflection, the slip, the rotation, and the bending moments. Both the results from the weighted residuals methods and that of the direct numerical methods are plotted. One finds an interesting agreement between the results presented by both methods. e profiles obtained are similar to the ones presented in [43]. Figure 4 shows the deflection distribution along the span, while Figure 5 shows the variation of the bending moment.
We can see that the maximum deflection (same for the bending moment) obtained numerically is very close to that obtained analytically. e small quantitative difference in         x (m) x (m)  x (m)  shows the variation of the rotation angle. At the center of the beam, the slip and rotation angle are equal to zero. en, they increase when one moves the center to the ends of the beams. ere was also not a slight difference between the results obtained from the mathematical development and those obtained from the direct numerical simulation. One notices that the observed smooth slip curve is effectively related to the smooth rotation curve. e comparison of the values of the deflection and that of the slips at important locations are given in Table 1.
An interesting question is to know what are the effects of the load and how the glue shear modulus affect the deflection and the slips. Let us mention that g(L) and θ(L) are, respectively, equal to −g(0) and −θ(0). To have an idea on such effects, we have first varied the values of the transversal load with k sc � 2.43 × 10 6 Pa. e purpose of this manipulation under the action of various transverse loads is to see if the composite beam made of concrete and bamboo can withstand the surface weight of a partitioning of building wall which varies between 40 daN/m 2 -100 daN/m 2 depending on its height. Figures 8-11 show the variation of the deflection, bending moment, slip, and rotation for three values of the load (these are the results of the weight residuals calculation). As one expects, the magnitudes of the variation of different mechanical variables increase with the load intensity. e threshold of the rupture of this composite beam is equal to 8 mm (L/500).
us, when a load of 3 kN/m is applied on the structure, one reaches the threshold of the rupture of this composite beam. e transverse load value at failure is 3.2 kN/m. Now, considering the effect of the glue type or the glue shear modulus, this is important since there are various types in the market. We have thus also varied the glue shear modulus, and Figure 12 shows the slip variations for three values of k sc . It is observed that the slip increases when the glue shear modulus decreases. is means that when the glue shear modulus is low, then that can cause large displacement between the two materials at their interface.

Conclusion
e objective of this study was to evaluate the spatial variation of the deflection, rotation, slip, and bending moment in a two-layer composite beam simply supported at its ends and under the action of a uniform transverse load on its entire length. e specific case used is a composite beam where one layer is made of bamboo and the other layer is the concrete. Timoshenko's beam model has been used. We have demonstrated that one can use the weighted residuals method to derive the mathematical expressions of the longitudinal deformation, the deflection, the rotation, and thus the slip and bending moment. Moreover, we have also demonstrated that one can use the finite differences method coupled to an iterative method for a set of algebraic equations to obtain through numerical simulation the deflection, rotation, and longitudinal displacement in the twolayer composite beam. An interesting agreement has been observed between the results obtained from the mathematical method and those from the direct numerical simulation using the finite differences method.
Supported by these mathematical and numerical methods, we have analyzed the effects of the load values and that of the glue shear modulus on the characteristic variables. It has been found that when the load increases, the deflection increases and the slip increases as well. When the glue shear modulus increases, it has been found that the slip becomes small. ese interesting results, which are qualitatively similar to those obtained by other researchers on different two-layer composite beams, require experiment for more confirmation of the developed mathematical and numerical simulation methods.
As it is known, most of the time, a problem such as the one solved here calls for the use of the finite element methods which already have black-boxes in which one just needs to introduce some specific parameters. e use of the finite difference method is thus another interesting way of solving the problem numerically, since the structure geometry is regular. Although, we have limited the study on a beam which is a one-dimensional system, extension can be undertaken in the case of plates.

Data Availability
No data were used to support the findings of this study.

Conflicts of Interest
e authors declare that they have no conflicts of interest.