2.1 Microstructural geometry of bamboo
a is the radial thickness of bamboo tube, and r is the radial coordinate.
As illustrated in Fig. 1(a), bamboo exhibits a mesoscale composition, consisting of vascular bundles and parenchyma cells. To create an simple and effective mesoscale geometry model, previous studies have often represented vascular bundles as cylindrical fibers and parenchyma cells as a continuous medium [14–17, 26]. From the cross-section of bamboo, it is evident that the fibers near the outer layer of the bamboo tube are thin and densely distributed. Moving from the outer to the inner layer, the distribution of fibers becomes sparser, accompanied by an increase in their radius. Therefore, our simulations incorporate a gradient distribution of vascular fibers in the cross-section [35–37]. The volume fraction of vascular bundles \({V}_{vb}\) can be defined as either a linear function Equ. (1) [38] or an exponential function Equ. (2) [25] of the dimensionless radius \(\stackrel{-}{r}\).
$${V}_{vb}={C}_{1}\stackrel{-}{r}+{C}_{2}$$
1
$${V}_{vb}={C}_{3}\text{exp}\left({C}_{4}\stackrel{-}{r}\right)$$
2
where \(\stackrel{-}{r}=r/a\), \(r\) is the radial coordinate, and \(a\) is the thickness of the bamboo tube. That is, \(\stackrel{-}{r}=0\) corresponds to the inner surface of the bamboo tube, and \(\stackrel{-}{r}=1\) corresponds to the outer surface of the bamboo tube. \({C}_{1}\), \({C}_{2}\), \({C}_{3}\) and \({C}_{4}\) are material parameters.
In order to accurately capture the mechanical behavior of the mesostructure while saving computational time, a parametric reconstruction methodology was utilized in this study. This approach simplified the representation of the vascular bundles and vessels within the bamboo structure, reducing them to cylindrical fibers. The linear Equ. (1) was employed to represent the volume fraction of the vascular bundles within the bamboo mesostructure, allowing for the modeling of its mechanical properties. This expression was subsequently compared to the Equ. (2). Figure 1(b) displays the reconstructed mesostructure of bamboo, illustrating the simplified representation of the vascular bundles and parenchyma cells.
2.2 Numerical and experimental specimens under different load conditions
The figures depicting the experimental and simulation specimens of bamboo are presented in Fig. 2. Mottled bamboo tubes from Chongqing City in China were chosen for the experiments. These bamboo tubes were approximately 5 years old and had an outer diameter of around 100 mm, a wall thickness ranging from 8 mm to 12 mm, and no visible defects. To control the water content, the tubes were dried in a blast oven, resulting in a moisture level of 8%-15%. Each group of mechanical experiments consisted of four to five specimens. This study specifically investigates the mechanical behavior of bamboo under various load conditions, including longitudinal tensile and compressive, transverse tensile and compressive, and longitudinal shear. The experimental specimens' geometries and the loading method used in both the experiments and numerical simulations were based on ISO 22157 [39] and JG/T 199 [40]. In Fig. 2, the numerical specimen (4 mm × 8 mm × t mm) chosen for longitudinal tension corresponds to the effective portion of the experimental specimen. By controlling the boundary conditions, the same loading effect can be achieved. For longitudinal compression, the specimen size is a regular rectangular (square) shape measuring 15 mm × 15 mm × t mm, while transverse compression is also conducted using regular rectangular specimens of the same dimensions. The specific size for the transverse tension specimen was not specified in the aforementioned experimental standards. Therefore, after conducting repeated tests to ensure stable responses, a size of 35 mm × 35 mm × t mm was selected. Here, t = 8 mm corresponds to the thickness of the bamboo tube in the simulation specimens, aligning with the experimental conditions. The simulated and experimental shear specimens possess identical shapes, and their sizes are determined according to the aforementioned standards.
The specimens subjected to longitudinal tension were clamped at both ends, ensuring that the narrow side of each specimen made contact with the jaws of the testing machine. Approximately 20 mm of each end was exposed near the curved section, and the specimens were mounted vertically on the testing machine. To determine their breaking point, the specimens were subjected to a uniform loading rate of 200 N/mm2 per minute. The load measurement accuracy was 10 N. Following the breaking test, the longitudinal tensile strength was converted to the strength at a water content of 12% using Equ. (3), with an accuracy of 0.1 N/mm2.
$$\left\{\begin{array}{c}{f}_{t,12}={K}_{t,w}{f}_{t,w} \\ {K}_{t,w}=\frac{1}{1.1-0.57{e}^{-0.15w}}\end{array}\right.$$
3
where \({f}_{t,12}\) is the longitudinal tensile strength at 12% water content, N/mm2. \({K}_{t,w}\) is the correction coefficient of the water content for the longitudinal tensile strength. \(w\)% is the water content of the specimen. In order to achieve uniform deformation, the transverse tension specimens were affixed to the chuck's end face instead of being clamped due to the absence of a clamping mechanism. Once fixed in place, the specimens were subjected to a uniform loading rate of 200 N/mm2 per minute until they broke, with a precision of 10 N for load measurements. The conversion equation for determining the transverse tensile strength at a water content of 12% is the same as Equ. (3).
For the longitudinal and transverse compression experiments, the specimens were centrally positioned on the spherical support of the testing machine. A consistent force of 80 N/mm2 per minute was then gradually applied until the specimen was completely crushed, with an accuracy of 10 N. Equ. (4) was used to convert the longitudinal and transverse compressive strength, based on a water content of w%, to a water content of 12%, with a precision of 0.1 N/mm2.
$$\left\{\begin{array}{c}{f}_{c,12}={K}_{c,w}{f}_{c,w} \\ {K}_{c,w}=\frac{1}{0.79+1.5{e}^{-0.16w}}\end{array}\right.$$
4
where \({f}_{c,12}\) is the longitudinal or transverse compressive strength at 12% water content, N/mm2. \({K}_{c,w}\) is the correction coefficient of the water content for the longitudinal compressive or transverse strength.
The shear experiments involve placing the specimen on the testing machine and applying a uniform load of 10 N/mm2 per minute until it reaches a state of damage. The load at which damage occurs is measured with an accuracy of 10 N. Equ. (5) is then utilized to calculate the shear strength at a moisture content of 12%.
$$\left\{\begin{array}{c}{f}_{v,12}={K}_{v,w}{f}_{v,w} \\ {K}_{v,w}=\frac{1}{0.67+0.77{e}^{-0.77w}}\end{array}\right.$$
5
where \({f}_{v,12}\) is the shear strength at 12% water content, N/mm2. \({K}_{v,w}\) is the correction coefficient of the water content for the shear strength.
2.3 Constitutive models for vascular bundles and parenchyma cells
Bamboo is a natural composite material that exhibits distinct characteristics from glass-reinforced and fiber-reinforced composites due to the anisotropy present in both its vascular bundles and parenchyma cells. Since the distribution of vascular bundles in bamboo differs from conventional composite materials, the direct application of a macroscopic constitutive model [41] is not appropriate. Consequently, two constitutive models must be defined, one for the fibers (vascular bundles) and the other for the matrix (parenchyma cells). At the mesoscale, both the fibers and matrix of bamboo are regarded as anisotropic materials, and their elastic components are determined using Hooke's law.
$$\left(\begin{array}{c}\begin{array}{c}{\sigma }_{1}\\ {\sigma }_{2}\\ {\sigma }_{3}\end{array}\\ \begin{array}{c}{\tau }_{12}\\ {\tau }_{23}\\ {\tau }_{31}\end{array}\end{array}\right)=\left[\begin{array}{cc}\begin{array}{ccc}{E}_{1}& -{E}_{1}/{\nu }_{12}& -{E}_{1}/{\nu }_{13}\\ & {E}_{2}& -{E}_{2}/{\nu }_{23}\\ & & {E}_{3}\end{array}& \begin{array}{ccc}0 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{array}\\ \begin{array}{ccc} & & \\ sym& & \\ & & \end{array}& \begin{array}{ccc}{G}_{12}& & \\ & {G}_{23}& \\ & & {G}_{31}\end{array}\end{array}\right]\left(\begin{array}{c}\begin{array}{c}{\epsilon }_{1}\\ {\epsilon }_{2}\\ {\epsilon }_{3}\end{array}\\ \begin{array}{c}{\gamma }_{12}\\ {\gamma }_{23}\\ {\gamma }_{31}\end{array}\end{array}\right)$$
6
where \({E}_{1}\), \({E}_{2}\), and \({E}_{3}\) are elastic modules for each direction, \({G}_{12}\), \({G}_{23}\), and \({G}_{31}\) are shear modules for the x-y, y-z and z-x planes, respectively.
Typically, the failure behavior of bamboo can be categorized into two groups: fiber failure and matrix failure. The constitutive models for the mesostructure of bamboo must satisfy two specific criteria. Firstly, the maximum stress failure criterion for fibers is derived from Rankine's theory [42] of maximum principal stress, as applied to isotropic materials. In this criterion, the stress acting on the cross-section is decomposed into normal and shear stresses. If the normal stress in the longitudinal direction exceeds the allowable strength of the fiber, the fiber is deemed to have failed. This approach considers the stress components separately and assumes that failure occurs when any one of the stress components reaches the allowable strength in the corresponding direction. Using this criterion, the following equalities are evaluated for fiber failure,
$$\left\{\begin{array}{c}\begin{array}{c}{\text{F}\text{I}}_{1}=\frac{{\sigma }_{1}}{{S}_{1}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{1}=-\frac{{\sigma }_{1}}{{S}_{1}^{C}}\\ {\text{F}\text{I}}_{2}=\frac{{\sigma }_{2}}{{S}_{2}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{2}=-\frac{{\sigma }_{2}}{{S}_{2}^{C}}\end{array}\\ \begin{array}{c}{\text{F}\text{I}}_{3}=\frac{{\sigma }_{3}}{{S}_{3}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{3}=-\frac{{\sigma }_{3}}{{S}_{3}^{C}}\\ {\text{F}\text{I}}_{12}=\left|\frac{{\tau }_{12}}{{S}_{12}^{T}}\right|\end{array}\\ \begin{array}{c}{\text{F}\text{I}}_{23}=\left|\frac{{\tau }_{23}}{{S}_{23}^{T}}\right|\\ {\text{F}\text{I}}_{31}=\left|\frac{{\tau }_{31}}{{S}_{31}^{T}}\right|\end{array}\end{array}\right.$$
7
where FI (Failure indices) indicates the occurrence of a failure, so that a value over 1 indicates that a failure has occurred. \({S}_{i}^{m}\)(\(i=1, \text{2,3}\); \(m=T,C\)) are the allowable stresses, and superscripts \(T\) and \(C\) stand for tension and compression, while subscripts 1, 2, and 3 correspond to the X, Y and Z directions, respectively; \({S}_{ij}\) is the allowable shear stress of the \(i-j\) plane. Failure would occur if any of these equalities were satisfied.
The failure criterion for the matrix is derived from St. Venant's theory of maximal normal strain and Tresca's theory of maximal shear stress [43]. The maximum tension strain, as proposed by these theories, significantly contributes to the occurrence of brittle fractures in bamboo. When subjected to longitudinal compression, transverse tension, or compression loads, the matrix fractures when the maximum tension strain at a specific point in the material exceeds its limit. Similar to the maximal stress theory, this criterion also adheres to the allowable strain theory. The failure criterion is expressed through the following equalities,
$$\left\{\begin{array}{c}\begin{array}{c}{\text{F}\text{I}}_{1}=\frac{{\epsilon }_{1}}{{e}_{1}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{1}=-\frac{{\epsilon }_{1}}{{e}_{1}^{C}}\\ {\text{F}\text{I}}_{2}=\frac{{\epsilon }_{2}}{{e}_{2}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{2}=-\frac{{\epsilon }_{2}}{{e}_{2}^{C}}\end{array}\\ \begin{array}{c}{\text{F}\text{I}}_{3}=\frac{{\epsilon }_{3}}{{e}_{3}^{T}} \text{o}\text{r} {\text{F}\text{I}}_{3}=-\frac{{\epsilon }_{3}}{{e}_{3}^{C}}\\ {\text{F}\text{I}}_{12}=\left|\frac{{\gamma }_{12}}{{e}_{12}}\right|\end{array}\\ \begin{array}{c}{\text{F}\text{I}}_{23}=\left|\frac{{\gamma }_{23}}{{e}_{23}}\right|\\ {\text{F}\text{I}}_{31}=\left|\frac{{\gamma }_{31}}{{e}_{31}}\right|\end{array}\end{array}\right.$$
8
where \({e}_{i}^{m}\)(\(i=1, \text{2,3}\); \(m=T,C\)) is the allowable strain in the corresponding direction, The superscripts and subscripts are the same as the maximum stress failure criteria. \({e}_{ij}\) is the allowable engineering shear strain of the \(i-j\) plane.
Following the definition of failure criteria, the selection of a damage progression model becomes the subsequent crucial step. In this study, the continuous damage mechanics (CDM) model is adopted to characterize the internal damage of the material by introducing six internal state variables, corresponding to each direction [44]. Although an instantaneous damage model may initially come to mind due to the immediacy of material failure, it tends to lack convergence. Therefore, this study employs a recursive model [45] to progressively diminish the material properties over multiple solution steps. This approach involves applying damage variables, ranging from 0 to 1, to the current constitutive matrix during each step following the initiation of damage, until the stiffness is significantly reduced (i.e., to 5% of the initial value). Consequently, damage evolution laws are associated with the internal damage variables within each anisotropic direction (i = 1, 2, …, 6).
2.4 The subroutine of constitutive models
The incremental approach is employed to numerically integrate the constitutive models for the fiber (vascular bundles) and matrix (parenchyma cells) into the user material subroutine UMAT of the nonlinear finite element software CalculiX [46]. This subroutine, written in FORTRAN, incorporates the constitutive models and their implementation process, as illustrated in Fig. 3. In the bamboo mesostructure, two materials, namely "fiber" and "matrix," are present, each with separate sub-UMATs (MAT-1 for fiber and MAT-2 for matrix) embedded in the main UMAT to handle the correct sub-UMAT utilization. The subroutine follows the steps outlined below,
1st, the local coordination of fiber and matrix, as well as the load and boundary conditions, are read. The FIs (Failure Index) and state variables are initialized, with the incremental index set to i = 1.
2nd, the process starts from the i-th increment, wherein the material properties of the fiber and matrix are read individually. Consequently, the local stress and strain tensors for both materials are calculated using linear equations. Subsequently, the corresponding FI values are determined.
3rd, the subroutine then checks whether the FI for the fiber or matrix is greater than 1. If the FI for the fiber exceeds 1, the program proceeds to the "subroutine for fiber damage" in order to execute the recursive damage procedure. This procedure involves calculating the damage variables, degrading the stiffness, and updating the stress and strain tensors correspondingly. Additionally, if the FI value for the matrix is less than 1, the state variable is updated. In cases where the FI for the matrix is greater than 1, the program advances to the "subroutine for matrix damage," wherein the calculation process is similar to fiber failure.
4th, the program then checks if the current increment is the last one. If not, the incremental index is updated by i = i + 1, and the computation skips back to step 2 for the calculation of a new increment. Once the last increment is reached, the UMAT subroutine is exited, and the stress, strain, and state variable results are outputted.
of the developed constitutive models.
2.5 Material parameters and boundary conditions
Based on the constitutive model of bamboo developed in Section 2.2, Table 1 gives the required material parameters for the fiber (vascular bundles) and matrix (parenchyma cells) of a type of mottled bamboo. The number of material parameters for both the fiber and matrix is 18, and the physical significance of each parameter is described in detail in the table.
Table 1
Material parameters for fiber and matrix
Matrix | Fiber |
No. | Name | Description | Value | Name | Description | Value |
1 | E11 | Young’s modulus, X direction | 160 | E11 | Young’s modulus, X direction | 3600 |
2 | E22 | Young’s modulus, Y direction | 160 | E22 | Young’s modulus, Y direction | 11300 |
3 | E33 | Young’s modulus, Z direction | 160 | E33 | Young’s modulus, Z direction | 3600 |
4 | µ12 | Poisson rations, 1–2 plane | 0.34 | µ12 | Poisson rations, 1–2 plane | 0.3 |
5 | µ13 | Poisson rations, 1–3 plane | 0.34 | µ13 | Poisson rations, 1–3 plane | 0.3 |
6 | µ23 | Poisson rations, 2–3 plane | 0.34 | µ23 | Poisson rations, 2–3 plane | 0.3 |
7 | G12 | Shear modulus, 1–2 plane | 59.7 | G12 | Shear modulus, 1–2 plane | 1384 |
8 | G13 | Shear modulus, 1–3 plane | 59.7 | G13 | Shear modulus, 1–3 plane | 4346 |
9 | G23 | Shear modulus, 2–3 plane | 59.7 | G23 | Shear modulus, 2–3 plane | 1384 |
10 | \({e}_{1}^{T}\) | Allowable tension strains, X direction | 0.035 | \({S}_{1}^{T}\) | allowable tension strengths, X direction | 248 |
11 | \({e}_{1}^{C}\) | Allowable compression strains, Y direction | 0.53 | \({S}_{1}^{C}\) | Allowable compression strengths, X direction | 120 |
12 | \({e}_{2}^{T}\) | Allowable tension strains, Y direction | 0.035 | \({S}_{2}^{T}\) | allowable tension strengths, Y direction | 248 |
13 | \({e}_{2}^{C}\) | Allowable compression strains, Y direction | 0.53 | \({S}_{2}^{C}\) | allowable compression strengths, Y direction | 120 |
14 | \({e}_{2}^{T}\) | Allowable tension strains, Z direction | 0.035 | \({S}_{3}^{T}\) | allowable tension strengths, Z direction | 248 |
15 | \({e}_{2}^{C}\) | Allowable compression strains, Z direction | 0.53 | \({S}_{3}^{C}\) | allowable compression strengths, Z direction | 120 |
16 | \({e}_{12}\) | allowable engineering shear strains, 1–2 plane | 0.12 | \({S}_{12}\) | allowable shear strengths, 1–2 plane | 124 |
17 | \({e}_{23}\) | allowable engineering shear strains, 2–3 plane | 0.12 | \({S}_{23}\) | allowable shear strengths, 2–3 plane | 124 |
18 | \({e}_{31}\) | allowable engineering shear strains, 3 − 1 plane | 0.12 | \({S}_{31}\) | allowable shear strengths, 3 − 1 plane | 60 |
*Note: The Young's modulus and strength values in the table are in MPa, and strain is a dimensionless quantity.
In the mesoscale simulations, all numerical specimens are modelled using hexahedral meshes for accurate calculation. For longitudinal and transverse tensions, one end of the specimen was fixed and a tensile force was applied to the other end at a rate of 200 N/mm2 per minute. The Y direction is longitudinal and the X direction is transverse, see Fig. 2 (a)(b). For longitudinal and transverse compression, one end was fixed and the compression force was loaded at 80 N/mm2 per minute onto the other end. The Y direction is longitudinal and the X direction is transverse, see Fig. 2 (c)(d). For shear simulations, the bottom of the specimen is fixed and then a downward load of 10 N/mm2 per minute is applied at the notched position, see Fig. 2 (e).