A continuum-based model for a layer jamming beam

Layer jamming structures (LJSs) are widely used as variable-stiffness components in collaborative robots to ensure safe human–robot interaction. However, existing analytical models often fail to adequately describe the mechanical behavior of LJSs with a large number of layers, particularly in capturing detailed stress distributions and deformations. This paper introduces a continuum-based layer jamming model (CLJM), which treats the LJS as a continuous medium under the assumption of infinitely many thin layers. The CLJM comprehensively analyzes internal stress distribution – including both shear and normal stresses – and deformation across different mechanical states (full jamming, half slipping, and full slipping). The model is validated through finite-element analysis (FEA) and experimental tests, showing strong agreement in both deformation response and stress profiles. The results demonstrate that the CLJM provides an effective and accurate theoretical tool for designing LJSs in variable-stiffness applications.

1 Introduction

For many structures, stiffness – an essential mechanical characteristic – remains constant during operation. However, as machines increasingly operate in unstructured environments, variable-stiffness structures are more effective than constant-stiffness ones in interacting with non-cooperative targets (Morrison and Su, 2020). The pursuit of variable stiffness is a pervasive theme in modern mechanical and robotic design, extending beyond layer jamming to include innovative mechanisms in various applications. For instance, in vehicle suspension systems, Anubi et al. (2013) demonstrated a passive variable-stiffness suspension that outperforms constant-stiffness designs by optimally controlling energy dissipation, highlighting the broad performance benefits of adaptive stiffness. In the realm of robotic joints, Yigit and Boyraz (2017) designed a cable-driven parallel-series hybrid variable-stiffness joint, addressing critical needs for compliance and workspace adaptability in humanoid robotics through a sophisticated mechanism incorporating elastic elements. These diverse approaches underscore the field's active exploration of variable-stiffness principles. Among various variable-stiffness techniques, jamming structures have garnered significant attention from researchers and engineers and have been applied in numerous fields (Blanc et al., 2017; Caruso et al., 2022a).

Jamming refers to the physical process in which discrete particles, layers, or fibers confined within a limited volume transition from liquid-like to solid-like behavior under increased stress (Brancadoro et al., 2020; Liu and Nagel, 1998). As a variable-stiffness technique, layer jamming beams offer features such as rapid and reversible transitions between different stiffness states and ease of fabrication (Blanc et al., 2017; Caruso et al., 2022a). Kawamura et al. (2002, 2003) and Tabata et al. (2000) were the first to introduce layer jamming structures (LJSs) into the engineering field as a variable impedance mechanism for wearable devices. Subsequently, LJSs have been employed in diverse applications, including medical endoscopic robots (Kim et al., 2012, 2013), soft actuators (Wall et al., 2015), grippers (Elgeneidy et al., 2019; Hou et al., 2019; Zhu et al., 2019), and cooperative robot arms (Zhou et al., 2020).

To further expand the engineering applications of LJSs, their mechanical behavior must be well understood to guide the design process. This area has been the focus of several studies. However, due to the complex interactions within jamming structures, effective numerical or analytical methods for modeling these structures are still needed. Specifically for LJSs, a widely accepted model focuses on the second moment of the cross-section, yielding $I=b\left(n\mathit{\delta }\right)\mathrm{3}/\mathrm{12}$ for a solid-like state and $I=bn\mathit{\delta }\mathrm{3}/\mathrm{12}$ for a liquid-like state (Henke and Gerlach, 2014; Ibrahimi et al., 2021; Kawamura et al., 2002), where n is the number of layers, δ is the thickness of the layers, and b is the width of the layers. Narang and colleagues (Narang, 2018; Narang et al., 2018a, b, 2020) were the first to present a comprehensive model of two-layer jamming cantilevers using Euler–Bernoulli beam theory. Their model not only predicts the stiffness of LJSs but also provides transition criteria, and it has been extended to multi-layer beams using finite-element methods. Furthermore, Caruso et al. (2022b) developed a model for multi-layer jamming cantilevers based on stress analysis, resulting in a summation form. They also introduced an analytical model to describe the three-point bending of multi-layer jamming beams (Caruso et al., 2022a, 2023). LJS beams bear similarity to composite beams with interlayer slip. Newmark (1951) studied composite beams with incomplete interaction, while Schnabl et al. (2007) provided an analytical solution for two-layer beams with interlayer slip and shear deformation. Tarsi and Afshin (2022) studied the delamination propagation in composite beams. However, these studies cannot adequately describe the deformation behavior of LJSs with significant slip.

Figure 1 Continuum-based layer jamming model. (a) LJM based on continuum assumption. (b) The stress state and its deforming path in the τ–σ stress space of the CLJM. (c) The elastoplastic deformation process of the LJM.

The analytical models of LJSs proposed in previous research often lack sufficient detail in stress distribution and model multi-layer mechanical behavior on a layer-by-layer basis. This approach becomes impractical when the number of layers is large. When a sufficient number of thin layers is present, it is unnecessary to distinguish between individual layers, and LJSs can be treated as a continuum. The continuum assumption is rigorously justified when the number of layers is sufficiently large, rendering the discrete effects negligible. This approach is analogous to methods employed in modeling laminated composites or granular materials, where the macroscopic behavior emerges from the collective interaction of a large number of constituents. In this work, we focus on the regime where the layer thickness δ is much smaller than the beam height h (i.e., $\mathit{\delta }/h\ll \mathrm{1}$), which means the number of layers n is large, ensuring the validity of treating stress and strain as continuous fields across the beam's height.

This paper is organized as follows. Section 2 introduces the basic concept and properties of the continuum-based layer jamming model (CLJM). Section 3 provides a detailed analysis of the internal forces and stresses in LJSs. Section 4 describes the deformation behavior of LJSs. Section 5 compares the CLJM results with experimental and FEA data. Finally, Sect. 6 presents the conclusions.

2 Continuum-based layer jamming model

The transition from a discrete layered system to a continuum is grounded in the limit where the layer thickness δ approaches zero, while the total height h remains constant, equivalent to an infinite number of layers. Under this assumption, which holds effectively for $n>\sim \mathrm{10}$ layers in practical applications, as shown in the FEA comparison (Sect. 5), mechanical quantities such as stress and strain become continuous functions of the spatial coordinate y. When the sheet is exceptionally thin, the layer jamming structure (LJS) can be considered to be a beam composed of a special type of continuous medium, based on the assumption that the layer thickness approaches zero, while the number of layers becomes infinite. In this case, all physical quantities can be treated as continuous functions along the beam height, enabling the use of continuum mechanics in the modeling process. The mechanical behavior of the LJS is described in this section using fundamental concepts from continuum mechanics.

From a continuous media perspective, the LJS may be modeled as a curved beam with a rectangular cross-section. To facilitate the description, the coordinate system in Fig. 1a is introduced: the s axis represents the central surface of the beam, the y axis represents the beam height, and θ is the angle between the y axis at s and at s=0. The radius of the curvature at s is $r=\mathrm{d}s/\mathrm{d}\mathit{\theta }$. To illustrate the relationship between the LJS and continuum mechanics, the stress state of the differential volume element at (s, y) in Fig. 1a is analyzed.

Similarly, Euler–Bernoulli beam theory assumes that all material points in the beam are in a plane stress state, neglecting interlaminar normal strain. Under this assumption, the stress state at (s, y) can be fully described by σ and τ, where σ represents the tensile force within a sheet, and τ represents the friction between sheets. The complete stress state can be represented as a point in the τ–σ diagram, as shown in Fig. 1b. Under an applied load, the differential volume element at (s, y) deforms, as illustrated in Fig. 1a. This deformation is characterized by normal strain ε_s and shear strain γ_sy. It should be noted that ε_s arises from the extension of sheets, while γ_sy is caused by both the extension of sheets and the slippage between them.

From the perspective of continuous media, LJSs can be viewed as curved beams made of specific elastoplastic materials, and their similarities are discussed later. Since interlayer normal strain is neglected, no elastic deformation-induced pressure exists between sheets. Consequently, the maximum static friction per unit area between sheets is given by μp, where p is the vacuum pressure acting on the LJS, and μ is the friction coefficient between layers. When $\left|\mathit{\tau }\right|<\mathit{\mu }p$, there is no slipping between sheets, only extension, and this stress state is referred to as the jamming state. Conversely, when $\left|\mathit{\tau }\right|\ge \mathit{\mu }p$, slip occurs, and the stress state is referred to as the slipping state. In the τ–σ diagram, the jamming state region appears as a gray strip around the σ axis, while the slipping state region is represented by two yellow strips above and below. Deformation in the jamming state is entirely reversible upon unloading, indicating elastic behavior. However, sliding in the slipping state is irreversible and can be regarded as plastic deformation. The stress state (σ, τ) exists only on the boundary or within the jamming state region.

The elastoplastic deformation behavior of the LJS medium is illustrated in Fig. 1c and can be divided into three stages: loading, yielding, and unloading. Each stage follows distinct principles, necessitating clear criteria for stage differentiation. These criteria, as illustrated in Fig. 1b, are as follows:

$$\left\{\begin{array}{l}\mathrm{Loading}\mathrm{stage}:\mathit{\tau }<\mathit{\mu }p\mathrm{and}\mathit{\tau }\mathrm{d}\mathit{\tau }>\mathrm{0}\\ \mathrm{Yielding}\mathrm{stage}:\mathit{\tau }=\mathit{\mu }p\mathrm{and}\mathit{\tau }\mathrm{d}\mathit{\tau }\ge \mathrm{0}\\ \mathrm{Unloading}\mathrm{stage}:\mathit{\tau }\le \mathit{\mu }p\mathrm{and}\mathit{\tau }\mathrm{d}\mathit{\tau }<\mathrm{0}.\end{array}\right.$$

It is important to note that shear strain is independent of shear stress and can take any value in the yielding stage. Aside from that, there is no one–one mapping between shear strain and shear stress. All of this suggests that the deformation of the LJS medium must be solved incrementally.

3 Internal force and stress analysis

3.1 Internal force analysis

From the perspective of the continuum, the mechanical behavior of LJS beams can be studied by utilizing the equivalent concentrated internal force on a section, that is, the axial force N, shear force Q, and bending moment M, as in Euler–Bernoulli beam theory.

Figure 2 The internal forces of a differential segment of a curved beam.

By concentrating on the differential beam segment ds in Fig. 2, the relationship between internal force and the loads on LJM beams is explored. As shown in Fig. 2, the center angle corresponding to ds is dθ, and the radius of curvature of the axis s is r. N and N+ dN are the axial forces on the left- and right-end sections, Q and Q+ dQ are the shear forces on the left- and right-end sections, and M and M+ dM are the bending moments on the left- and right-end sections of the differential beam segment, respectively. q is the distributed load on the curved beam. The arrows in Fig. 2 indicate the positive direction of each mechanical quantity. By ignoring the second-order small quantities and considering the approximation of cos dθ≈1 and sin dθ≈ dθ, the force and moment balance equations on the left-end section can be expressed as follows:

$$\begin{array}{}\text{(1)}& \left\{\begin{array}{l}\frac{\mathrm{d}N}{\mathrm{d}s}=-Q\mathit{\theta }{}^{\prime }-q_{\mathrm{s}}\\ \frac{\mathrm{d}Q}{\mathrm{d}s}=N\mathit{\theta }{}^{\prime }+q_y\\ \frac{\mathrm{d}M}{\mathrm{d}s}=Q\end{array}\right.\iff \left\{\begin{array}{l}\frac{\mathrm{d}N}{\mathrm{d}s}=-\frac{Q}{r}-q_{\mathrm{s}}\\ \frac{dQ}{\mathrm{d}s}=\frac{N}{r}+q_y\\ \frac{\mathrm{d}M}{\mathrm{d}s}=Q\end{array}\right.,\end{array}$$

where q_s and q_y represent the components of q along the s and y axes, respectively. The first two equations of Eq. (1) can give the solution written as

$$\begin{array}{}\text{(2)}& \begin{array}{l}\left[\begin{array}{c}N\\ Q\end{array}\right]=\left[\begin{array}{cc}\cos \mathit{\theta }& -\sin \mathit{\theta }\\ \sin \mathit{\theta }& \cos \mathit{\theta }\end{array}\right]\left(\left[\begin{array}{c}{N}_{\mathrm{0}}\\ Q_{\mathrm{0}}\end{array}\right]\right.\\ \left.+\underset{s_{\mathrm{0}}}{\overset{s}{\int }}\left[\begin{array}{cc}\cos \mathit{\theta }& \sin \mathit{\theta }\\ -\sin \mathit{\theta }& \cos \mathit{\theta }\end{array}\right]\left[\begin{array}{c}-q_{\mathrm{s}}\left(\mathit{\xi }\right)\\ q_y\left(\mathit{\xi }\right)\end{array}\right]\mathrm{d}\mathit{\xi }\right),\end{array}\end{array}$$

where N₀ and Q₀ are the internal forces at s₀ that can be obtained by solving the LJS beams' reaction force. Following that, M(s) can be calculated by integrating from Q(s).

It should be noted that the conclusion of this section holds for the region of both jamming and slipping LJS beams since its deduction depends only on the continuum assumption.

Figure 3 The section deformation state of the LJM beam. (a) Full-jamming state. (b) Half-slipping state. (c) Full-slipping state.

3.2 Stress analysis

Next, we focus on the shear stress in LJS beams. It is evident that when the shear force Q on a section is sufficiently small, there will be no sliding region, and the section is considered to be in the full-jamming state, which is shown in Fig. 3a. For a section in the full-jamming state, the shear stress distribution of a rectangular straight beam can approximate that of a curved beam because the radius of curvature r of LJS beams is much larger than their height h. This implies that for a section with a free surface, the following equation holds:

$$\begin{array}{}\text{(3)}& \mathit{\tau }\left(y,Q\right)={\displaystyle \frac{\mathrm{3}}{\mathrm{2}}}{\displaystyle \frac{Q}{A}}\left(\mathrm{1}-{\displaystyle \frac{\mathrm{4}{y}^{\mathrm{2}}}{h^{\mathrm{2}}}}\right),\end{array}$$

where A is the area of the rectangular cross-section, A=bh, and b is the width of the beam.

As Q on a section increases, the central region of the section transitions from the jamming state to the slipping state because τ follows a quadratic distribution along the y axis. Such a section is referred to as being in the half-slipping state, as shown in Fig. 3b. However, in general, sections with a free surface in an LJS beam cannot transition entirely into the slipping state. This is because the outermost layers will always remain in the elastic (jamming) state. As Q continues to increase, the limiting case is reached when no jamming region remains for $\left|y\right|\le h/\mathrm{2}$ – δ, where δ is the thickness of a single layer. This state is referred to as the full-slipping state, as depicted in Fig. 3c.

The distribution of τ_sy on a section in the half-slipping state can then be derived. Within the jamming region, τ follows a form similar to Eq. (3), but Q is replaced by Q_eq, which is determined by solving τ(y_s, Q_eq)=μp, where y_s is the boundary between the jamming and sliding regions. The result is

$$\begin{array}{}\text{(4)}& Q_{\mathrm{eq}}=\sin gQ{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}{\displaystyle \frac{\mathit{\mu }pA}{\mathrm{1}-\mathrm{4}{y}_{\mathrm{s}}^{\mathrm{2}}/h^{\mathrm{2}}}}.\end{array}$$

The distribution of τ in the jamming region can then be obtained by substituting Q in Eq. (3) with Q_eq as given in Eq. (4).

The critical shear force Q_slip, which represents the threshold between the full-jamming and half-slipping states, can be determined by solving max_y [τ(y, Q)] =μp. From Eq. (3), it is evident that τ reaches its maximum value at y=0. When the section transitions from the half-slipping state to the full-slipping state, the boundary y_s satisfies $y_{\mathrm{s}}=h/\mathrm{2}-\mathit{\delta }$. In this case, substituting Q in Eq. (3) with Q_eq in Eq. (4) and integrating τ across the section yield another critical shear force, Q_max:

$$\begin{array}{}\text{(5)}& {\displaystyle }{Q}_{\mathrm{slip}}={\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\mathit{\mu }pbh,\text{(6)}& {\displaystyle }{Q}_{\max }=\mathit{\mu }pb{\displaystyle \frac{\mathrm{3}{h}^{\mathrm{2}}-\mathrm{6}\mathit{\delta }h+\mathrm{4}{\mathit{\delta }}^{\mathrm{2}}}{\mathrm{3}\left(h-\mathit{\delta }\right)}}.\end{array}$$

These expressions are derived based on the fact that τ equals sign(Q)μp in the sliding region $\left|y\right|\le y_{\mathrm{s}}$.

To summarize, the shear stress τ(y, s) in an LJS beam can be expressed as follows:

$$\begin{array}{}\text{(7)}& \mathit{\tau }\left(y,s\right)=\left\{\begin{array}{ll}\frac{\mathrm{3}}{\mathrm{2}}\frac{Q}{A}\left(\mathrm{1}-\frac{\mathrm{4}{y}^{\mathrm{2}}}{h^{\mathrm{2}}}\right),& \left|Q\right|\le Q_{\mathrm{slip}}\\ \mathrm{sign}\left(Q\right)\mathit{\mu }p,& \left|y\right|\le y_{\mathrm{s}}\mathrm{and}{Q}_{\mathrm{slip}}\le \left|Q\right|\le Q_{\max }\\ \mathrm{sign}\left(Q\right)\mathit{\mu }p\frac{\mathrm{1}-\mathrm{4}{y}^{\mathrm{2}}/h^{\mathrm{2}}}{\mathrm{1}-\mathrm{4}{y}_{\mathrm{s}}^{\mathrm{2}}/h^{\mathrm{2}}},& \left|y\right|>y_{\mathrm{s}}\mathrm{and}{Q}_{\mathrm{slip}}\le \left|Q\right|\le Q_{\max }\\ c_{\mathrm{2}}{y}^{\mathrm{2}}\pm c_{\mathrm{1}}y+c_{\mathrm{0}},& \left|y\right|>y_{\mathrm{s}}\mathrm{and}\left|Q\right|>Q_{\max }.\end{array}\right.\end{array}$$

Here, the coefficients c₀, c₁, and c₂ can be determined by the following conditions: $\mathit{\tau }(h/\mathrm{2},s)=\mathrm{0}$, τ(y_s, s)= sign(Q)μp, and ${\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\mathit{\tau }\mathrm{d}y=Q$. The resulting expressions for c₀, c₁, and c₂ are as follows:

$$\begin{array}{}\text{(8)}& \begin{array}{l}\left\{\begin{array}{l}{c}_{\mathrm{2}}=-\frac{\mathrm{12}\left(\mathit{\mu }pbh+\mathrm{2}\mathit{\mu }pb{y}_{\mathrm{s}}-\mathrm{2}Q\right)}{b{\left(h-\mathrm{2}{y}_{\mathrm{s}}\right)}^{\mathrm{3}}}\\ c_{\mathrm{1}}=\frac{\mathrm{4}\left(\mathrm{3}Qh+\mathrm{6}Q{y}_{\mathrm{s}}-\mathrm{2}\mathit{\mu }pb{h}^{\mathrm{2}}-\mathrm{4}\mathit{\mu }pbh{y}_{\mathrm{s}}-\mathrm{8}\mathit{\mu }pb{y}_{\mathrm{s}}^{\mathrm{2}}\right)}{b{\left(h-\mathrm{2}{y}_{\mathrm{s}}\right)}^{\mathrm{3}}}\\ c_{\mathrm{0}}=\frac{h\left(\mathit{\mu }pb{h}^{\mathrm{2}}+\mathrm{2}\mathit{\mu }pbh{y}_{\mathrm{s}}+\mathrm{16}\mathit{\mu }pb{y}_{\mathrm{s}}^{\mathrm{2}}-\mathrm{12}Q{y}_{\mathrm{s}}\right)}{b{\left(h-\mathrm{2}{y}_{\mathrm{s}}\right)}^{\mathrm{3}}}.\end{array}\right.\end{array}\end{array}$$

Finally, the most critical parameter, y_s, can be obtained by solving ${\int }_{-h/\mathrm{2}}^{h/\mathrm{2}}\mathit{\tau }\mathrm{d}y=Q$. Setting y_s=0 in the full-jamming state and $y_{\mathrm{s}}=h/\mathrm{2}-\mathit{\delta }$ in the full-slipping state, y_s can then be expressed as

$$\begin{array}{}\text{(9)}& \begin{array}{l}{y}_{\mathrm{s}}=\left\{\begin{array}{ll}\mathrm{0},& \left|Q\right|<Q_{\mathrm{slip}}\\ \frac{\mathrm{3}Q-\mathrm{2}\mathit{\mu }pA+\sqrt{\mathrm{3}\left(\mathrm{3}Q-\mathrm{2}\mathit{\mu }pA\right)\left(Q+\mathrm{2}\mathit{\mu }pA\right)}}{\mathrm{8}\mathit{\mu }pb},& Q_{\mathrm{slip}}\le Q<Q_{\max }\\ \frac{-\mathrm{3}Q+\mathrm{2}\mathit{\mu }pA+\sqrt{\mathrm{3}\left(\mathrm{3}Q-\mathrm{2}\mathit{\mu }pA\right)\left(Q+\mathrm{2}\mathit{\mu }pA\right)}}{\mathrm{8}\mathit{\mu }pb},& -Q_{\max }<Q\le -Q_{\mathrm{slip}}\\ h/\mathrm{2}-\mathit{\delta },& \left|Q\right|\ge Q_{\max }.\end{array}\right.\end{array}\end{array}$$

Next, we turn to the distribution of normal stress in LJS beams. The distribution of normal stress σ varies depending on the deformation state. For the full-jamming state, σ distribution is equivalent to that of an elastic beam. In the jamming regions of the half-slipping and full-slipping states, the normal stress follows a linear distribution for reasons previously discussed. Considering an LJS beam subject to a bending moment only, each layer can be treated as a beam. Consequently, the function σ(s, y) in the sliding region can be represented as a series of polylines, as shown in Fig. 4. The amplitude of these polylines decreases to zero as the layer thickness approaches zero, indicating an infinite number of layers. Combining these observations, the normal stress distribution can be expressed as

$$\begin{array}{}\text{(10)}& \mathit{\sigma }=\left\{\begin{array}{ll}{N}_J/A_J+K_J\left(y\mp y_{\mathrm{s}}\right),& \left|y\right|>y_{\mathrm{s}}\\ {\mathit{\sigma }}_{\mathrm{S}}\left(y\right),& \left|y\right|<y_{\mathrm{s}}.\end{array}\right.\end{array}$$

Here, the subscript J denotes quantities defined in the jamming region; that is, N_J represents the axial force in the jamming region, A_J is the cross-sectional area of the jamming region ($A_J=b(h-\mathrm{2}{y}_{\mathrm{s}})$), and K_J is a proportional coefficient. In addition, the symbol ∓ in Eq. (10) takes the minus sign when y≥y_s and the plus when $y\le -y_{\mathrm{s}}$. σ_S denotes the normal stress distributing on the sliding regions, which is discussed in detail later.

Figure 4 The normal stress analysis in the jamming region.

The relationships between NJ, MJ (for momentum on the jamming region) and N, M are

$$\begin{array}{}\text{(11)}& \begin{array}{r}\left\{\begin{array}{l}N=\underset{\mathrm{Jamming}}{\int }{\mathit{\sigma }}_{J}b\mathrm{d}y+\underset{\mathrm{Sliding}}{\int }{\mathit{\sigma }}_{\mathrm{S}}b\mathrm{d}y\\ M=\underset{\mathrm{Jamming}}{\int }{\mathit{\sigma }}_{J}yb\mathrm{d}y+\underset{\mathrm{Sliding}}{\int }{\mathit{\sigma }}_{\mathrm{S}}yb\mathrm{d}y\end{array}\right.\Rightarrow \\ \left\{\begin{array}{l}N=N_J+\underset{-y_{\mathrm{s}}}{\overset{y_{\mathrm{s}}}{\int }}{\mathit{\sigma }}_{S}b\mathrm{d}y\\ M=M_J+\underset{-y_{\mathrm{s}}}{\overset{y_{\mathrm{s}}}{\int }}{\mathit{\sigma }}_{S}yb\mathrm{d}y,\end{array}\right.\end{array}\end{array}$$

where σ_S denotes the normal stress in sliding regions. And M_J can be determined by integrating Eq. (10) over the jamming region:

$$\begin{array}{}\text{(12)}& M_J={\displaystyle \frac{b}{\mathrm{12}}}\left(h+y_{\mathrm{s}}\right){\left(h-\mathrm{2}{y}_{\mathrm{s}}\right)}^{\mathrm{2}}{K}_J.\end{array}$$

Figure 5 The normal stress analysis in the sliding region. (a) The normal stress distribution in the sliding region. (b) Neither of the sliding region's end sections is in the half-slipping state. (c) At least one of the end sections is in the half-slipping state.

To determine the distribution of σ_S in the sliding region, we apply the differential equation of stress balance:

$$\begin{array}{}\text{(13)}& {\displaystyle \frac{\partial {\mathit{\sigma }}_{\mathrm{S}}}{\partial s}}+{\displaystyle \frac{\partial {\mathit{\tau }}_{\mathrm{S}}}{\partial y}}=\mathrm{0}.\end{array}$$

Since τ_S in the sliding region is a constant sign (Q)μp, we have $\partial {\mathit{\tau }}_{\mathrm{S}}/\partial y=\mathrm{0}$. Solving Eq. (13) yields σ_S=f(y), where f is an unknown function, indicating that σ_S depends only on s in the sliding regions. This implies that σ_S takes the same value at all points with the same y in the sliding region (e.g., all points on the red line in Fig. 5a). Consequently, σ_S(y) dy=σ_S(ys) dy_s and σ_S(y)y dy=σ_S(y_s)y_s dy_s.

Using this property and the relation dy_s=y_s^′ ds, the integrals over y in Eq. (11) can be transformed to integrals over s for sections in the half-slipping state. There are two cases. (1) Neither end of the sliding region is in the half-slipping state, as shown in Fig. 5b. (2) At least one end of the sliding region is in the half-slipping state, as shown in Fig. 5c. For the first case, the original integral over y (domain indicated by the red line in Fig. 5b) can be divided into two integrals over s (domain indicated by the blue curves), as mathematically described:

$$\begin{array}{}\text{(14)}& \left\{\begin{array}{ll}\underset{-y_{\mathrm{s}}}{\overset{y_{\mathrm{s}}}{\int }}{\mathit{\sigma }}_{\mathrm{S}}b\mathrm{d}y& =\underset{s}{\overset{s_{\mathrm{0}}}{\int }}{\mathit{\sigma }}_{\mathrm{S}}\left(-y_{\mathrm{s}}\right)b{y}_{\mathrm{s}}{}^{\prime }\left(\mathrm{d}s\right)+\underset{s_{\mathrm{0}}}{\overset{s}{\int }}{\mathit{\sigma }}_{\mathrm{S}}\left(y_{\mathrm{s}}\right)b{y}_{\mathrm{s}}{}^{\prime }\mathrm{d}s\\ & =\underset{s_{\mathrm{0}}}{\overset{s}{\int }}\left[{\mathit{\sigma }}_{\mathrm{S}}\left(y_{\mathrm{s}}\right)+{\mathit{\sigma }}_{\mathrm{S}}\left(-y_{\mathrm{s}}\right)\right]b{y}_{\mathrm{s}}{}^{\prime }\mathrm{d}s\\ \underset{-y_{\mathrm{s}}}{\overset{y_{\mathrm{s}}}{\int }}{\mathit{\sigma }}_{\mathrm{S}}yb\mathrm{d}y& =\underset{s}{\overset{s_{\mathrm{0}}}{\int }}{\mathit{\sigma }}_{\mathrm{S}}\left(-y_{\mathrm{s}}\right)\left(-y_{\mathrm{s}}\right)b{y}_{\mathrm{s}}{}^{\prime }\left(\mathrm{d}s\right)+\underset{s_{\mathrm{0}}}{\overset{s}{\int }}{\mathit{\sigma }}_{\mathrm{S}}\left(y_{\mathrm{s}}\right)y_{\mathrm{s}}b{y}_{\mathrm{s}}{}^{\prime }\mathrm{d}s\\ & =\underset{s_{\mathrm{0}}}{\overset{s}{\int }}\left[{\mathit{\sigma }}_{\mathrm{S}}\left(y_{\mathrm{s}}\right)+{\mathit{\sigma }}_{\mathrm{S}}\left(-y_{\mathrm{s}}\right)\right]y_{\mathrm{s}}b{y}_{\mathrm{s}}{}^{\prime }\mathrm{d}s,\end{array}\right.\end{array}$$

where s₀ is one end of the half-slipping state segment.

Based on the continuum assumption, the normal stress should remain continuous across the boundary, which means that there are σ_S(y_s)=σ_J(y_s) and ${\mathit{\sigma }}_{\mathrm{S}}(-y_{\mathrm{s}})={\mathit{\sigma }}_J(-y_{\mathrm{s}})$. Substituting Eqs. (14) and (10) into Eq. (11) and noticing ${\mathit{\sigma }}_J\left(y_{\mathrm{s}}\right)-{\mathit{\sigma }}_J(-y_{\mathrm{s}})=\mathrm{0}$ can give

$$\begin{array}{}\text{(15)}& \left\{\begin{array}{l}N=N_J+\underset{s}{\overset{s}{\int }}\mathrm{2}\frac{N_J}{A_J}b{y}_{\mathrm{s}}{}^{\prime }\mathrm{d}\mathit{\xi }\\ M=M_J.\end{array}\right.\end{array}$$

The derivative of Eq. (15) gives the differential equation about N_J:

$$\begin{array}{}\text{(16)}& N{}^{\prime }=N_J{}^{\prime }+N_J{\displaystyle \frac{y_{\mathrm{s}}{}^{\prime }}{h/\mathrm{2}-y_{\mathrm{s}}}}=-Q\mathit{\theta }{}^{\prime }.\end{array}$$

Considering that at s₀ there is N_J=N, the solution of Eq. (16) is

$$\begin{array}{}\text{(17)}& N_J=\left(\mathrm{2}{y}_{\mathrm{s}}/h-\mathrm{1}\right)\left[\underset{s}{\overset{s}{\int }}-{\displaystyle \frac{Q\mathit{\theta }{}^{\prime }}{\mathrm{2}{y}_{\mathrm{s}}/h-\mathrm{1}}}\mathrm{d}\mathit{\xi }+N\left(s_{\mathrm{0}}\right)/\left(\mathrm{2}{y}_{\mathrm{s}}\left(s_{\mathrm{0}}\right)/h-\mathrm{1}\right)\right].\end{array}$$

For the second case, the original integral over y, whose domain is indicated by the red line in Fig. 5c, can be divided into three parts: two integrals over s, whose domain is indicated by the blue curves, and an integral over y, whose domain is indicated by the purple line. Applying the same transformation technique used in Eq. (14), we can obtain a differential equation that is the same as Eq. (15), and it can also give a similar solution to the first case, but the values of N_J at s₀ should correspond to the distribution of σ_J at s₀. Supposing that the axial force at s₀ is evenly distributed on the section gives the normal stress distribution at s₀:

$$\begin{array}{}\text{(18)}& \left\{\begin{array}{ll}{\mathit{\sigma }}_J\left(s_{\mathrm{0}}\right)& =\frac{N\left(s_{\mathrm{0}}\right)}{A}+K_J\left(s_{\mathrm{0}}\right)\left(y\mp y_{\mathrm{s}}\right)\\ & =\frac{N\left(s_{\mathrm{0}}\right)}{A}+\frac{\mathrm{12}M\left(s_{\mathrm{0}}\right)}{b\left(h+y_{\mathrm{s}}\right){\left(h-\mathrm{2}{y}_{\mathrm{s}}\right)}^{\mathrm{2}}}\left(y\mp y_{\mathrm{s}}\right)\\ {\mathit{\sigma }}_{\mathrm{S}}\left(s_{\mathrm{0}}\right)& =N\left(s_{\mathrm{0}}\right)/A,\end{array}\right.\end{array}$$

where K_J is solved by combining Eqs. (12) and (15). It is then easy to compute that N_J(s₀)=N(s₀)2y_s $\left(s_{\mathrm{0}}\right)/h$.

4 Deformation analysis

The goal of the deformation analysis is to obtain Δθ(s) under a given load q. Using q and θ(s), the internal forces can be obtained from Eq. (2). Notably, the deformation of an LJS beam derived from the jamming region and sliding region is consistent. This consistency allows the deformation of an LJS beam to be calculated solely based on the internal moment within the jamming region.

For small deformations of an LJS beam, the following relationship holds in the jamming region, based on Euler–Bernoulli beam theory:

$$\begin{array}{}\text{(19)}& \mathrm{\Delta }\mathit{\theta }=\underset{\mathrm{0}}{\overset{s}{\int }}{\displaystyle \frac{M}{E{I}_J}}\mathrm{d}\mathit{\xi },\end{array}$$

where I_J represents the moment of inertia of the jamming region in LJS beams. Its value can be calculated as follows:

$$\begin{array}{}\text{(20)}& I_J=\left\{\begin{array}{ll}\frac{b{h}^{\mathrm{3}}}{\mathrm{12}},& y_{\mathrm{s}}=\mathrm{0}\\ \frac{\mathrm{2}}{\mathrm{3}}\left[{\left(\frac{h}{\mathrm{2}}\right)}^{\mathrm{3}}-y_{\mathrm{s}}^{\mathrm{3}}\right],& y_{\mathrm{s}}\ne \mathrm{0}.\end{array}\right.\end{array}$$

The mechanical response of the LJS beam under a large load q must be solved incrementally, based on the small deformation assumption. A detailed algorithm for this step-by-step procedure is presented in the Algorithm.

Figure 6 FEA analysis of straight LJM cantilever's stress distribution. (a) Distribution of σ and τ under 3N load. (b) Distribution of σ and τ under 4N load. (c) Distribution of σ and τ under 5.5N load.

5 Results

To further validate the mechanical behavior predicted by the proposed method, stress distributions obtained using finite-element analysis (FEA) provide a more detailed comparison. The LJS cantilevers analyzed in FEA measure 40 mm in length, 20 mm in thickness, and 5 mm in height. These cantilevers consist of 10 and 25 layers, with assumed material parameters of a friction coefficient μ=0.5 and Young's modulus E=0.3 GPa. The finite-element models were constructed in Abaqus, utilizing CPS4R elements with a side length of 0.1 mm. A pressure of 0.1 MPa was applied on the upper and lower surfaces, and a concentrated force was evenly distributed along the free ends of each layer to facilitate convergence.

The shear and normal stress distributions predicted by the continuum-based model and FEA are shown in Fig. 6. For full-jamming sections, Fig. 6a illustrates the stress distribution at s=30 mm under a 3N load. Here, the stress distributions predicted by FEA and the continuum-based model align closely, supporting the assumption that the jamming region behaves like an elastic material. However, near the ends of the cantilever, such as at s=1 mm under a 3N load, the shear stress distribution changes significantly depending on the number of layers. As the number of layers increases, the amplitude of τ(y) predicted by FEA decreases and becomes more continuous, though it deviates from the parabolic shape due to boundary effects.

For half-slipping sections, Fig. 6b shows the stress distribution at s=10, s=10, and s=30 mm. In the middle segments of the cantilever, the stress distributions predicted by FEA and the continuum-based model show minimal differences. However, the sliding region predicted by FEA is narrower than that predicted by the continuum-based layer jamming model (CLJM). This discrepancy arises from the CLJM's neglect of normal stress along the y axis. Ignoring this stress leads to an overestimation of the sliding region, as the normal stress along the y axis can raise the critical shear stress, reducing the extent of sliding. Consequently, the slope of the normal stress in the jamming region predicted by FEA is smaller than that predicted by CLJM.

Figure 7 Experimental analysis of straight LJM cantilever's deformation. (a) The measuring device. (b) Experimental result.

For full-slipping sections, Fig. 6c highlights the differences in stress distributions. Similar to the half-slipping sections, the sliding region predicted by FEA is narrower than that predicted by CLJM, and the shear stress in the sliding region obtained by FEA fluctuates around the values predicted by CLJM. However, a significant discrepancy exists in the jamming region's shear stress between FEA and CLJM. This error likely occurs because LJS beam segments in the full-slipping state behave more like curved beams, while CLJM assumes that the shear stress distribution in the jamming region follows the pattern of elastic straight beams, introducing notable inaccuracies.

The relationship between load and deflection in a straight LJS cantilever can be easily measured using a tensile machine. The experimental setup used in this study is shown in Fig. 7a. In this setup, the straight LJS cantilever is fixed to a frame that slides freely along a rail to eliminate horizontal restraints. A flexible PVC membrane encasing 20 PVC sheets forms the test specimen, which measures 100 mm in length, 20 mm in thickness, and 5 mm in height.

The load-deflection results under 60 kPa pressure are represented as gray dots in Fig. 7b. The orange curve denotes the experimental mean, while the error bars spanning ±1 standard deviation are shown as the yellow-shaded area. The analytical results, depicted by the red curve in Fig. 7b, were obtained using the previously discussed algorithm and the material parameters: Young's modulus E=755 MPa and friction coefficient μ=0.55. The analytical curve lies entirely within the 1σ interval, demonstrating the validity of the proposed analysis method for describing the mechanical behavior of the LJS cantilever.

These discrepancies, particularly near the boundaries and in the full-slipping state, highlight the inherent limitations of the current continuum assumption, which does not explicitly account for end effects or the curved beam behavior under large slip. Nevertheless, the proposed CLJM demonstrates high accuracy in predicting global deformation and stress distribution in the central regions of the beam under small-to-moderate slip conditions, which are of primary interest in most practical applications. Future refinements could incorporate boundary correction factors or a more detailed treatment of curvature in the jamming region.

6 Conclusions

In this study, an analytical model based on continuum mechanics is presented to describe the mechanical behavior of LJS beams. Previous works typically treated each layer in the LJS beam as an independent beam, which is inadequate for handling LJSs with numerous layers and describing stress distribution in detail. In contrast, the continuum layer jamming model presented here treats the LJS as a continuous medium, assuming that all mechanical quantities are continuous across the material. This model enables the description of various mechanical behaviors of the LJS beam, such as stress distribution, deformation, and other related phenomena. The effectiveness of CLJM is validated through experiments and FEA. The deformation predicted by CLJM aligns well with experimental measurements, and the stress distribution predicted by CLJM closely matches the results obtained from FEA. Ultimately, it can be concluded that the assumption of treating LJSs as a continuum is reasonable, providing engineers with a valuable tool for designing materials to meet the specific demands of various engineering applications.

The CLJM provides a robust analytical foundation for designing LJSs with a large number of layers. It is most accurate for predicting global mechanical behavior away from constraint boundaries and under partial-slip conditions. The model's simplicity and computational efficiency make it particularly suitable for preliminary design and optimization, while FEA remains recommended for final verification, especially near free ends or under extreme loading leading to full slipping.