Multiscale Domain Decomposition Accurately Predicts Rolling Shear Failure in Timber Panels

A multiscale computational framework integrating domain decomposition with cohesive interface mechanics can accurately predict out-of-plane rolling shear failure and subsequent fracture progression in cross-laminated timber panels while substantially reducing computational overhead. Researchers J. Fernández, K. Saavedra, O. Allix, and P. Gosselet introduced the simulation approach in a research study published in Advanced Modeling and Simulation in Engineering Sciences. By aligning numerical substructure partitioning directly with the principal stress trajectories of timber lamellae and enhancing iterative convergence algorithms, the modeling strategy achieves up to a 55 percent reduction in degrees of freedom, a 70 percent decrease in required subdomains, and a 70.5 percent reduction in solver iterations compared to traditional partitioning configurations. Cross-laminated timber, widely recognized as a sustainable alternative to conventional building materials due to its favorable strength-to-weight ratio and stiffness, relies on multiple orthogonally bonded lumber layers whose structural capacity is frequently governed by rolling shear stresses within transverse plies.

Cross-laminated timber panels are composed of stacked lumber layers bonded with structural adhesives, configured in odd-numbered arrangements of three, five, or seven plies where each successive layer is oriented perpendicularly to adjacent boards. Although these engineered wood systems have been increasingly deployed in mid-rise and high-rise construction across Canada, Norway, France, and the United Kingdom, predicting their complex structural limitations remains a persistent engineering challenge. At the mesoscopic material scale, wood behaves as an orthotropic material with distinct elastic properties along its longitudinal, radial, and tangential coordinate axes. When subjected to out-of-plane bending loads, this structural arrangement induces critical stress states in transverse cross-layers characterized by low shear strength. These transverse shear stresses produce rolling shear failure, an internal mechanism where wood fibers twist and rupture, typically initiating within central lamellae in the absence of pre-existing cracks and complicating the application of classical fracture mechanics.

To evaluate the spatial distribution of internal forces, the research team applied the Shear Analogy Method following Mestek’s methodology to investigate four-point bending setups across panels featuring three, five, and seven layers of identical wood grade, layer width, and thickness. In these analytical evaluations, coniferous timber conforming to European Standard EN 338 was modeled as transversely isotropic in global coordinates, with transverse elastic moduli and shear stiffnesses expressed as proportional fractions of the longitudinal modulus, including a rolling shear modulus evaluated at one-hundred-and-sixtieth of the longitudinal stiffness. The analytical stress field analyses revealed that principal stresses in the transverse layers are dominated by rolling shear stresses, displaying characteristic inclination angles close to 45 degrees in three-layer and seven-layer panels and approximately 52 degrees in five-layer configurations. These calculated principal stress orientations align directly with experimental observations of diagonal crack paths, providing the fundamental geometric rationale for aligning numerical domain boundaries with physical fracture orientations.

Rather than relying solely on empirical strength failure criteria that compare local stress states against characteristic thresholds without capturing post-failure damage progression, the investigators adopted a damage mechanics framework employing cohesive zone models. The approach builds upon a mesomodel that represents laminated panels through volumetric layer entities coupled with two-dimensional interface elements. While earlier mesomodel formulations focused primarily on interlaminar interfaces to simulate delamination between plies, the current methodology introduces intralaminar cohesive interfaces to capture diverse failure modes within individual layers, including rolling shear, transverse fiber cracking, and normal tension failure under bending. Interface degradation is governed by thermodynamic forces reflecting dissipated surface energy, with the constitutive cohesive law derived from characteristic timber strengths, initial interface stiffness, and an exponent controlling the potential curve shape, thereby ensuring mathematical equivalence with Griffith fracture energy criteria across modes I, II, and III.

The governing equations are solved through an adapted multiscale Large Time Increment, or LaTIn, mixed domain decomposition strategy implemented within the parallel C++ finite element platform MULTI. To solve the structural boundary value problem, the spatial domain is partitioned into subdomains connected across interfaces, dividing the system into linear subdomain equations and local non-linear interface relations. Global kinematic and static admissibility conditions are maintained alongside linear orthotropic behavior within individual substructures, while scalability is preserved through a macroscopic admissibility condition formulated via an orthogonal projector. Iterative solutions cycle systematically between the linear admissibility manifold and the non-linear interface manifold through forward and backward search directions. Computations execute in parallel environments utilizing the METIS algorithm for automated domain partitioning and the Message Passing Interface library for inter-processor communication, allowing substructures to solve local problems independently before macroscopic reduction and assembly.

To establish baseline spatial accuracy, the authors examined an initial three-point bending benchmark on a three-layer beam measuring 420 millimeters in length, 50 millimeters in width, and 60 millimeters in total height, comprising outer MSG8 lumber layers and an MSG6 inner layer under an applied force of 1154 Newtons. The numerical representation discretized the specimen into 480 subdomains of 10-millimeter length using 131,040 linear three-dimensional wedge elements yielding 335,664 degrees of freedom. Comparisons between the finite element stress fields and the analytical Shear Analogy Method confirmed satisfactory agreement across normal and shear stress distributions at quarter-span. This controlled verification demonstrated that the multiscale LaTIn strategy accurately captures linear elastic stress distributions across laminated timber plies prior to the initiation of material non-linearities and localized cracking.

The research team subsequently investigated crack initiation, propagation, and fracture process zone behavior using a standard double cantilever beam delamination benchmark. The modeled specimen measured 280 millimeters in length, 20 millimeters in width, and 20 millimeters in total height, with an initial crack length of 100 millimeters and a prescribed separation displacement applied in 32 incremental steps. Analytical formulations indicate that the characteristic length of the fracture process zone is governed by the layer elastic modulus, critical fracture energy, and directional strength, resulting in a theoretical cohesive process zone length spanning between 1.0 and 2.6 times the layer thickness. Numerical trials evaluating subdomain lengths of 5, 10, and 20 millimeters showed that matching substructure dimensions to the fracture process zone scale at 20 millimeters yielded superior surface stress representations without artificial discontinuities, while an element size of 2 millimeters satisfied the requirement of placing approximately ten finite elements across the active process zone.

Although the double cantilever beam benchmark reproduced the physical onset of delamination following initial linear elastic loading, initial tests exhibited pronounced iteration counts to achieve the convergence tolerance of 0.0001, prompting specific algorithmic optimizations. The authors implemented an open contact condition by setting the forward search direction to zero across delaminated interfaces, updated macroscopic projection operators whenever interface elements became completely detached, and introduced a linear extrapolation scheme to initialize interface traction fields across successive time increments during regular monotonic loading. The open contact update alone reduced total solver iterations by 52 percent, lowering computation time from 17.8 hours to 8.12 hours on a 64-node cluster equipped with 2.3-gigahertz AMD Opteron processors. Combining open contact updates with macro-operator refreshing and linear force extrapolation achieved a 70.5 percent reduction in overall iterations and reduced wall-clock solution time to 6.75 hours.

Building upon these numerical refinements, the framework was deployed to simulate a full four-point short-span bending experiment on a three-layer cross-laminated panel measuring 1260 millimeters in length, 600 millimeters in width, and 105 millimeters in total height, with individual 35-millimeter plies subjected to a prescribed downward displacement of 22 millimeters across 50 increments. Taking advantage of geometric symmetry, the simulation modeled one-quarter of the beam using 422,400 wedge elements and 166,464 quadrilateral interface elements across 1920 subdomains. The middle transverse layer incorporated triangular subdomains whose internal interfaces were oriented at 45 degrees, matching the analytically established rolling shear stress directions. Specific interface designations were incorporated across the domain: 45-degree rolling shear interfaces in the core layer, horizontal delamination and transverse strength interfaces along layer boundaries, and vertical bending strength interfaces in the bottom tensile ply.

Constitutive parameters for each interface category were assigned using the simplified wood-strength identification methodology, adopting an exponent of two to reflect the quasi-brittle character of timber. An initial interface stiffness of 300 Newtons per cubic millimeter was assigned to the damageable transverse interfaces to prevent overly compliant behavior, while a stiffer value of 800 Newtons per cubic millimeter was applied to vertical bottom interfaces to preserve initial flexural stiffness without distorting the global linear response. Rolling shear strengths of 1.6, 1.8, and 2.0 megapascals were systematically evaluated. When considered in isolation, rolling shear crack initiation produced an immediate 6.5 percent drop in overall panel stiffness. Once the critical shear stress threshold was exceeded, damage propagated rapidly along the 45-degree interfaces before triggering adjacent fracture mechanisms, executing in approximately 7.8 hours of total cluster computation time.

The multi-interface simulations revealed that rolling shear cracks do not trigger adhesive delamination at the ply boundary; instead, they transition into horizontal intralaminar cracks propagating through the wood fibers adjacent to the bond line, traveling outward along the lower boundary and inward along the upper interface. This progressive horizontal cracking degraded panel bending stiffness until tensile stresses in the lower longitudinal ply exceeded the material capacity, causing ultimate bending rupture beneath the loading point. The simulated mid-span moment-displacement curves, initial slope drop, horizontal crack paths, and ultimate fiber failure matched experimental test measurements reported in published literature. By tailoring subdomain geometries to expected failure paths, the approach achieved a 55 percent reduction in total degrees of freedom and a 70 percent decrease in required subdomains compared to previously published orthogonal partitioning strategies requiring 6528 subdomains and 2.9 million degrees of freedom.

The authors concluded that while numerical modeling does not replace physical laboratory experiments, the multiscale domain decomposition framework provides a robust and computationally tractable tool for simulating complex, interacting failure modes in cross-laminated timber structures. By capturing the progressive transition from localized rolling shear cracking to intralaminar fracture and catastrophic flexural rupture, the approach deepens understanding of timber mechanics and offers predictive support for structural design standards. Future research aims to extend the domain decomposition strategy to more complex architectural assemblies, including five-layer and seven-layer timber configurations as well as hybrid cross-laminated timber and concrete composite building systems.

Subject of Research: Multiscale domain decomposition modeling of rolling shear failure in cross-laminated timber

Article Title: Predicting rolling shear failure of cross-laminated panels using a multiscale domain decomposition approach

Article References: Fernández, J., Saavedra, K., Allix, O., & Gosselet, P. (2026). Predicting rolling shear failure of cross-laminated panels using a multiscale domain decomposition approach. Advanced Modeling and Simulation in Engineering Sciences, 13(1), Article 1. https://doi.org/10.1186/s40323-026-00318-9

Image Credits: AI Generated

DOI: 10.1186/s40323-026-00318-9

Keywords: cross-laminated timber, rolling shear failure, domain decomposition, LaTIn strategy, cohesive zone model, structural mechanics, finite element analysis, timber engineering, Predicting, rolling, shear, failure