The rockeries in classical gardens are currently affected by various forms of deterioration, among which structural damage is of particular concern. Previous studies indicate that such damage mainly includes stress concentration, peak stone inclination, foundation settlement, adhesive detachment, and plant-induced root splitting1. In current research on the stability of traditional garden rockeries, some researchers have conducted structural risk assessments based on physical vulnerability indicators derived from finite element simulations2, while others have identified stress-concentrated areas in rockeries under self-weight and visitor loading through finite element analysis3.
In the field of historical monument studies, many structures exhibit significant deterioration, with cracks being one of the most distinct surface defects. Consequently, crack behavior has attracted considerable attention in the numerical analysis of historical monuments. Many researchers have investigated advanced computational methods for finite element modeling of crack propagation in brittle materials4,5. Detailed on-site visual inspections combined with finite element modeling are commonly used to evaluate the structural stability of such heritage objects. In our previous work, we investigated the crack expansion risk threshold of the Small Rock Mountain Adobe (SRMA) rockery in He Garden6.
It is worth noting that once cracks develop, they progressively accelerate the deterioration of cultural heritage, leading to further structural damage and even severe risks such as overall instability or collapse7. However, no existing literature has applied the finite element method to numerically simulate the staged evolution of overall instability in rockery structures or to determine the corresponding risk thresholds. Establishing a numerical simulation methodology for assessing the structural stability of cultural heritage is essential for scientifically predicting the damage process, implementing targeted monitoring and rescue measures, and prolonging the lifespan of heritage assets8. This approach also provides practical engineering guidance for heritage managers.
Finite element analysis offers a robust computational framework for evaluating structural risks in cultural heritage. Current research includes finite element modeling of heritage structures9, static and dynamic analyses of structural responses under various loading conditions10,11,12, and comprehensive assessments of structural stability risks13. These advancements facilitate sustainable conservation through deeper structural understanding and evidence-based restoration strategies14. For stone cultural relics such as masonry structures and grottoes, software like Abaqus and elastoplastic constitutive models are commonly used to evaluate structural stability, for example, employing plastic models to diagnose damage in historic masonry within large cultural heritage sites15. Some researchers also use the elastoplastic damage model and finite element method to analyze the structural safety of the hollow defense gun emplacements of the Great Wall16. In the finite element analysis of rock deformation, Drucker-Prager and MC constitutive models are often adopted.
In previous studies, the risk levels and corresponding risk thresholds for the expansion of the most structurally critical cracks in the central cave area of Small Rock Mountain Adobe (SRMA) in He Garden were evaluated6. Because the rockeries in Jiangnan gardens are predominantly built from limestone—a brittle rock material—crack propagation mainly occurs within the linear elastic stage, and different propagation risk levels can be classified using the key parameter SIF.
However, that study focused solely on the crack scale and did not address the overall instability risk threshold of the rockery when cracks advance into the high-risk stage. Moreover, it did not incorporate the deformation and failure behavior of the plastic stage17 into the numerical simulation framework for assessing structural stability. Therefore, this study attempts to incorporate an elastoplastic constitutive model to conduct a finite element numerical simulation of the entire process of the instability of the rockery structure, in order to determine the various stages of the instability risk of the rockery structure and the corresponding risk thresholds.
For brittle, hard limestone, a limit strain value exists under cyclic loading; as this limit is approached, the material experiences continuous deformation followed by rupture18. In current research on rock-structure stability, determining the overall failure behavior of rock materials through the L–D curve is a feasible approach19. Furthermore, based on the stress-strain curve, the failure process of brittle rocks can be divided into several stages20,21.
In this study, the assessment of the stability risk threshold of the rockery structure was extended from the crack-level stage to the structural-level stage. The crack-level stage refers to the stage from the initial formation of cracks to the stage of linear crack propagation in the rockery structure. The structural-level stage refers to the stage involving overall structural damage, instability, and eventual collapse of the rockery structure.
In previous studies, linear elastic fracture mechanics (LEFM)-based numerical simulations were used to model the cracking patterns and crack propagation paths of rockeries, and corresponding risk levels and risk thresholds were determined. LEFM theory (adopted to model crack propagation), which assumes a linear relationship between stress and strain until plastic deformation occurs, is effective for modeling brittle rock failure22. The stress-intensity factor solution obtained from the LEFM theory can be used to evaluate crack propagation23. The maximum principal stress criterion24,25 was adopted in this study, whereby damage initiation starts when the maximum principal stress at a point exceeds the tensile strength of the material. Crack propagation follows a linear evolution model based on the energy release rate, propagating in the direction that maximizes the total energy release.
The LEFM theory assumes a linear elastic response near crack tips and characterizes the stress field using the stress intensity factor (SIF). For Mode I loading, the SIF is defined as
$${K}_{I}=\sigma \sqrt{\pi a}$$
(1)
where \(\sigma\) is the applied stress, and \(a\) is the crack length.
The maximum principal stress criterion is used as the criterion for crack initiation, which holds that the crack will occur when the maximum principal stress reaches the ultimate strength of the material. Material failure is assumed to occur when the maximum principal stress \({\sigma }_{\max }\) reaches the ultimate strength of the material \({\sigma }_{c}\), i.e.,
$${\sigma }_{\max }={\sigma }_{c}$$
(2)
The energy release rate (\(G\)) describes the energy change during crack propagation. This represents the energy released by the system per unit increase in the crack area. For three-dimensional crack propagation, when the energy release rate \(G\) reaches the critical energy release rate \({G}_{c}\) of the material, the crack continues to propagate.
One of the main drawbacks of the aforementioned linear elastic fracture mechanics is that it neglects plasticity. Even when the overall stress is elastic, plastic deformation can still occur at the crack tip, even in materials that are generally considered brittle.
Therefore, an elastoplastic constitutive model is adopted.
Using the Drucker–Prager yield function, we get
$$F=\sqrt{{J}_{2}}+\alpha {I}_{1}-k\left(\kappa \right)$$
(4)
Where \({I}_{1}\) is the first invariant of the stress tensor, \({J}_{2}\) is the second invariant of the stress deviator, and \(\alpha\) is a material constant26. The internal variable \(\kappa ={\varepsilon }_{{eq}}^{p}\) represents the equivalent plastic strain (an output variable in Abaqus, PEEQ), and \(k\left(\kappa \right)\) denotes the strength parameter that governs the evolution of the yield surface during plastic deformation.
In this section, the ABAQUS numerical simulation adopts the uniaxial compression mode, that is, \({\sigma }_{1}\,=\,{\sigma }_{{\rm{c}}}\,\), and \({\sigma }_{2}\,=\,{\sigma }_{3}=0\). Therefore \(k\left(\kappa \right)\) is written as
$$\frac{{\sigma }_{{\rm{c}}}}{\sqrt{3}}+\alpha {\sigma }_{{\rm{c}}}=k\left(\kappa \right)$$
(5)
In this section, the hardening-softening function is adopted in the Drucker–Prager function setting of ABAQUS to observe the deformation behavior of the rockery material before and after the peak, and it is expressed as follows:
In the pre-peak stage, the material exhibits hardening behavior, which is expressed as:
$$\frac{{dk}}{d\kappa } > 0,0\le \kappa \le {\kappa }_{p}$$
(6)
In the post-peak stage, the material exhibits softening behavior, expressed as:
$$\frac{{dk}}{d\kappa } < 0,\kappa > {\kappa }_{p}$$
(7)
Among them, \({\kappa }_{p}={\varepsilon }_{{eq},p}^{p}\) denotes the equivalent plastic strain at the point of reaching the peak stress.
These key nodes include (1) crack initiation, (2) formation of a through-crack, (3) structural instability, and (4) structural collapse. Based on these risk thresholds, the instability risk of the rockery structure can be divided into five stages: Ⅰ Stable structure stage, Ⅱ Stage of linear crack propagation, Ⅲ Stage of structural damage, Ⅳ Stage of structural instability, and Ⅴ Stage of structural collapse.
Using the above formula, the key nodes of (1) crack initiation and (4) structural collapse can be obtained. Based on these results, and in combination with the L-D curve, the stability risk stages of the rockery structure can be further divided in a more detailed manner. The physical meanings and calculation bases of the key nodes are as follows:
The risk threshold during the crack-level stage includes the following two key points:
Crack initiation: The maximum principal stress criterion is used to determine crack initiation. According to this criterion, a crack occurs when the maximum principal stress reaches the ultimate strength of the material, i.e., \({\sigma }_{\max }={\sigma }_{c}\).
Formation of a through-crack: A significant turning point appears on the L-D curve.
The risk threshold during the structural-level stage includes the following two key points:
Structural instability: The trend of the L–D curve undergoes a significant change.
Structural collapse: The bearing strength of the rockery materials reaches its peak, i.e., \({\varepsilon }_{{eq},p}^{p}={\kappa }_{p}\).
The instability risk threshold proposed in this study is more innovative compared to previous work. In our previous work, we used LEFM combined with XFEM to predict the rockery cracking patterns and crack propagation paths.
In this study, it was observed that once the L-D curve reached the critical point corresponding to the “formation of a through crack”, its overall trend no longer changed significantly. This indicates that the overall structural response behavior of the rock formation can no longer be captured using LEFM alone. By incorporating the elastoplastic constitutive model, the key node corresponding to the peak point of stress-strain curve can be quantified to determine the critical node of “structural collapse.” Furthermore, by analyzing the evolution of the L-D curve, the risk stages of the rock formation structural instability can be further quantified. This approach enables the prediction of the entire process of structural stability risk in rock formations through numerical simulations, forming a comprehensive methodological framework with engineering practicality.
This addresses the existing gap in quantifying unstable risk stages in rockery structures and supports garden managers in implementing preventive and emergency protection measures, particularly in areas undergoing plastic deformation.
Comments are closed.