Mechanical investigations of the peltate leaf of Stephania japonica (Menispermaceae): Experiments and a continuum mechanical material model

Stephania japonica is a slender climbing plant with peltate, triangular-ovate leaves. Not many research efforts have been devoted to investigate the anatomy and the mechanical properties of this type of leaf shape. In this study, displacement driven tensile tests with three cycles on different displacement levels are performed on petioles, venation and intercostal areas of the Stephania japonica leaves. Furthermore, compression tests in longitudinal direction are performed on petioles. The mechanical experiments are combined with light microscopy and X-ray tomography. The experiments show, that these plant organs and tissues behave in the finite strain range in a viscoelastic manner. Based on the results of the light microscopy and X-ray tomography, the plant tissue can be considered as a matrix material reinforced by fibers. Therefore, a continuum mechanical anisotropic viscoelastic material model at finite deformations is proposed to model such behavior. The anisotropy is specified as the so-called transverse isotropy, where the behavior in the plane perpendicular to the fibers is assumed to be isotropic. The model is obtained by postulating a Helmholtz free energy, which is split additively into an elastic and an inelastic part. Both parts of the energy depend on structural tensors to account for the transversely isotropic material behavior. The evolution equations for the internal variables, e.g. inelastic deformations, are chosen in a physically meaningful way that always fulfills the second law of thermodynamics. The proposed model is calibrated against experimental data, and the material parameters are identified. The model can be used for finite element simulations of this type of leaf shape, which is left open for the future work.

1 Introduction

The plant leaf is exposed to a variety of mechanical influences including environmental stresses and stresses from self-loading (Niklas, 1992; Niklas, 1999). The leaf must withstand these stresses and twist and bend without being damaged to secure its function as a platform for photosynthesis and production of organic compounds (Niklas, 1999; Adams and Terashima, 2018). Plants have evolved a multitude of different leaf shapes with many of them being composed of petiole and lamina. From a mechanical perspective, the petiole can be defined as a cantilevered beam fixed on one side and the free moving lamina attached to the other (Niklas, 1999). The lamina itself is composed of a scaffold of leaf veins supporting and supplying the intercostal areas in between that consist of mesophyll which is the main photosynthetic tissue of the leaf (Adams and Terashima, 2018). The internal composition and organisation of different tissues define the mechanical properties of the structural elements of the leaf and thus, its stability (Niklas, 1992). The thin-walled parenchyma and the collenchyma with partially thickened cell walls are hydrostatic tissues. Their mechanical behavior is greatly affected by turgor pressure (Falk et al., 1958; Niklas, 1992). In contrast, fibers and xylem elements are sclerenchymatous, non-hydrostatic tissues. They are characterized by thickened, lignified cell walls and elongated shape (Napp-Zinn, 1973). The polymer lignin within the cell walls of these tissue types acts as additional reinforcement resulting in significantly higher elastic moduli than parenchyma and collenchyma (Niklas, 1992; Gibson, 2012).

While anatomical and biomechanical properties of leaves are generally well studied, not many research efforts have been devoted to investigate the peltate leaf shape. Peltate leaves are defined by the petiole inserting on the abaxial side of the lamina resulting in a 3-dimensional spatial arrangement (Troll, 1932; Langer et al., 2021b; Wunnenberg et al., 2021). The peltate leaf shape is not very common in the plant kingdom. The most recent taxon list comprises approximately 350 peltate-leaved species representing 99 genera and 40 families (Wunnenberg et al., 2021). The 3D architecture of this leaf shape determines the stresses the leaf must withstand, especially in the petiole and the transition zone between petiole and lamina. Petioles of peltate leaves need to be resistant or flexible to bending and torsion in the same intensity in all directions due to the 3D spatial arrangement of the lamina. In comparison, leaves with a 2D configuration need to be resistant in bending in one direction while still being flexible in torsion (Langer et al., 2021a). The transition zone is characterized by a significant change in geometry from petiole to lamina, plays a crucial role in the dissipation of mechanical loads and is optimized to cope with different loads in comparison to the petiole (Sacher et al., 2019; Langer et al., 2021b; Langer et al., 2022).

An anatomical study on peltate leaves (Wunnenberg et al., 2021) already showed a complex fiber organisation in the petiole-lamina transition zone of two Stephania Lour. species. In Stephania delavyi Diels and Stephania venosa (Blume) Spreng., a ring-like structure was described for the transition zone making this genus especially interesting for further analysis.

Stephania japonica (Thunb.) Miers is a slender climbing plant species with peltate leaves that belongs to the Menispermaceae. It is native to tropical and subtropical Asia and Australia and grows in hedges, thickets, secondary growths, forests and along river banks at altitudes from 0 to 2000 m (Forman, 1986; Dataset POWO, 2022). It has herbaceous or woody stems that can grow to heights of 10 m. The leaves are petiolate. The petiole is 3 to 12 cm long and can be glabrous or pilose. The lamina is triangular-ovate to ovate in shape with acuminate apex and broadly rounded base and 4 to 17 × 4 to 14 cm in size. The adaxial side of the lamina is glabrous while the abaxial side can be glabrous or pilose (Forman, 1986).

Experiments have shown that the uniaxial tensile nonlinear stress-strain response of petiole, venation and intercostal area is rate dependent. Furthermore, anatomical investigations of the leaf have led to the assumption that the plant tissue is a fiber-reinforced material with one fiber direction. To model the mechanical behavior of such plant tissues, a number of complex features including finite deformation, anisotropy and viscoelasticity have to be taken into account. Although multi-scale methods can be helpful in describing such materials (see e.g. Stier et al., 2015), we focus here on the computationally cheaper macro-mechanical approach.

Several phenomenological approaches for fiber-reinforced composites have been proposed. For instance, Holzapfel and Gasser (2001) introduced a viscoelastic model for two groups of fibers. Nedjar (2007) and Nguyen et al. (2007) presented independently two models that treat the matrix and the fibers separately allowing as many bundles of fibers as desired. Latorre and Montáns (2015) presented an extension of the Reese and Govindjee (1998) framework to quasi-incompressible transversely isotropic and orthotropic materials, later they Latorre and Montáns (2016) presented a formulation based on a reverse (to that from Sidoroff) multiplicative decomposition of the deformation gradient. Furthermore, Liu et al. (2019) provided a model, based on the Holzapfel type anisotropic hyperelastic strain-energy function, with a new method to develop the evolution equations of the viscous internal variables.

Some works have explicitly focused on developing models for soft biological tissues. Bischoff et al. (2004) provided a rheological model for fibrous tissue, Nguyen et al. (2008) presented a constitutive model for the corneal stroma. Garcia-Gonzalez et al. (2018) developed a framework for the mechanical modelling of a wide variety of soft tissues that incorporates strain rate and temperature dependencies as well as the transverse isotropy. Moreover, Stumpf (2021) proposed a constitutive framework capable of modeling a wide variety of different soft tissues, such as skin, tendon, colorectal tissue and brain white matter.

However, according to the knowledge of the authors, an experimentally based sophisticated model for the physical mechanisms of the plant leaf as a whole and the peltate leaf of S. japonica cannot be found in the literature.

In perspective of possible engineering applications, in this study we aim to give a detailed anatomical and mechanical description of the peltate leaf of S. japonica. To simulate the mechanical behavior of the petiole, venation and intercostal areas, we present a constitutive model of transversely isotropic viscoelasticity at finite deformations. The model is based on a general phenomenological approach for fiber-reinforced composites, since leaf components are considered as a continuum mixture consisting of fibers orientated in one direction and embedded in a soft isotropic matrix. The model is based on a decomposition of the deformation gradient into elastic and viscous parts. Furthermore, we postulate the Helmholtz free energy, which is decomposed into isotropic and anisotropic parts. Both parts are additionally split into elastic and inelastic components. Here, we propose that the viscoelastic behaviors of the matrix and fibers mutually contribute to the nonlinear viscoelastic response of the leaf structures.

2 Materials and methods

2.1 Plant material and sampling

Leaf samples of S. japonica were provided by the Botanical Garden of Technische Universität Dresden, Germany. The plant (Figures 1A, B) is cultivated in a pot under open air conditions during summer and in the greenhouse during winter. Only intact, undamaged adult leaves were collected and transported in airtight containers. The fresh leaves were photographed (Lumix DMC-G81, Panasonic, Kodama, Osaka, Japan) from both sides and subsequently processed within a few hours or preserved in ethanol ( 70 % ) for later analysis.

2.2 Anatomy and morphology

The software ImageJ (National Institutes of Health, Bethesda, Maryland, USA) was used for measurement of leaf and petiole dimensions for all samples. Leaf length was determined at the longest point along the midrib, leaf width perpendicular to the midrib at the widest part. MS Excel (Microsoft Corporation, Redmont, Washington, USA) was used for compilation of measured parameters.

Transverse sections of the petiole and transverse and longitudinal sections of the petiole-lamina transition zone were prepared with a razor blade or vibratome (Hyrax V50, Carl Zeiss AG, Jena, Germany). Sections were bleached in sodium hypochlorite solution (2.8 %), washed in water and stained with astrablue/safranin. Astrablue stains non-lignified tissue blue and safranin stains lignified tissue red. Subsequently, the sections were differentiated in ethanol (70 %). The stained sections were photographed using a light microscope (VHX-970F, Keyence AG, Osaka, Japan) and integrated camera (VHX-970F, CMOS-image sensor, Keyence AG, Osaka, Japan). The amount of lignified tissue in the petiole was measured for three large leaves. For each of these leaves, five sections from basal, central and apical part of the petiole were analyzed using the software ImageJ (National Institutes of Health, Bethesda, Maryland, USA). Following data analysis was performed using OriginLab 2021 Pro (OriginLab Corporation, Northhampton, Massachusetts, USA).

For X-ray tomography, the petiole-lamina transition zone was cut out of a fresh leaf using a razor blade (sample dimensions: 10 mm in height). The lamina was cut away from the transition zone without damaging the latter. Afterwards, the sample was wrapped in several layers of parafilm (Bemis Inc., Neenah, USA) and glued (UHU instant glue, UHU GmbH & Co. KG, Bühl, Germany) to a small wooden sample holder (dimensions: 5×5 mm ). Tomography measurements were done on the micro-CT ProCon CT-XPRESS (ProCon X-Ray GmbH, Sarstedt, Germany) and analyzed using the software Dragonfly 3.5 (Object Research Systems Inc., Montreal, Canada). The axial CT scans were performed with a source voltage of 30 kV and a source current of 70 µA (voxel size: 2.5 µm).

2.3 Mechanical tests

Individual samples of petiole, venation and intercostal areas were subjected to displacement driven tensile and compression tests with three cycles on different strain levels. The compression tests were only performed on petiole samples. Samples of 8 leaves in total were prepared and tested (for sample orientation, see Figure 1C). For tensile testing, petioles were cut into pieces of 50 mm . The venation was separated from intercostal areas using a razor blade and cut to a length of 50 mm . The intercostal areas were cut into 5×50 mm samples. For compression tests, samples of 2 mm length were cut from the petioles. Width of petiole and venation samples and thickness of intercostal area samples were measured using a digital caliper (Precise PS 7215, Burg-Wächter, Wetter-Volmarstein, Germany). Per testing speed, five samples of each structural element were tested. To minimize water loss, the samples were stored between wet paper sheets until testing and sprayed with water during testing.

FIGURE 1

Figure 1 Growth habit of Stephania japonica (A, B) and sampling method for mechanical testing (C). Petiolate leaves with a triangular-ovate to ovate lamina (A) and young climbing shoots (B). Examplary positions and orientations of leaf samples for mechanical tests (C), v - venation, ic - intercostal area, p - petiole, scale: 2 cm.

Two universal testing machines (zwickiLine and Zwick/Roell AllroundLine Z005, Zwick/Roell GmbH & Co. KG, Ulm, Germany) and corresponding software (testXpert II V3.5, Zwick/Roell GmbH & Co. KG, Ulm, Germany) were used. Clamping length was set to 20 mm; pre-load was set to 0.1 N for petiole samples and 0.05 N for venation and intercostal area samples. Resulting forces were recorded using a 50 N load cell for the zwickiLine (Type: KAP-Z, AST Angewandte Systemtechnik GmbH, Dresden, Germany) and a 5 kN load cell for the AllroundLine Z005 (Zwick/Roell xforce P, Zwick/Roell GmbH & Co. KG, Ulm, Germany). Tests were performed at three testing speeds: 1 % min⁻¹, 10 % min⁻¹ and 500 % min⁻¹. The engineering strain (ε: = u/l₀) levels for the three cycles were set to: 0.5 %, 1 % and 2 %. In each cycle the engineering strain was applied using a ramped constant rate load function. After the prescribed engineering strain value was achieved, the engineering strain was held for the some loading time, as presented in Figure 2, to show the relaxation process of the tissue. Next, the engineering strain was removed at a constant rate (see Figure 3). A small initial force was applied at the start in order to straighten the specimens. Tests were limited with the total testing time of 30 min, since the samples were drying out during the experiment.

FIGURE 2

Figure 2 Strain-driven cyclic experiments.

FIGURE 3

Figure 3 (A) Rheological model. (B) Initial, current and intermediate configurations of a body, and the multiplicative decomposition of the deformation gradient F.

Force-displacement measurements were obtained. For an objective comparison of experimental data, strain and stress measures have to be introduced. Since only uniaxial experiments were performed, stresses and strains can be computed in the following way

$\begin{array}{ll}{[\mathbf{\text{F}}]}_{11}:=\frac{u+l_0}{l_0},& (1)\end{array}$

$\begin{array}{ll}{[\mathbf{\text{P}}]}_{11}:=\frac{R}{A},& (2)\end{array}$

where R is the measured force, A is the initial cross sectional area of the specimen, l₀ is the initial length of the sample, and u is the measured displacement. F is the deformation gradient and P is the first Piola-Kirchhoff stress tensor (see Section 2.4 for other definitions).

2.4 Material model formulation

In the following, we develop a continuum mechanical transversely isotropic viscoelastic material model at finite deformations to model the behavior of the plant tissues. The three structural elements of the S. japonica leaf (petiole, venation, intercostal area) are represented as a continuum mixture consisting of one group of fibers bound by an isotropic matrix. To model the behavior of the tissue, we postulate a Helmholtz free-energy. The stress response of the tissue is schematically described by the rheological model shown in Figure 3A. In this rheological model, a single independent spring represents the time-independent equilibrium response which acts in parallel to two Maxwell elements representing two different time-evolving nonequilibrium processes. We assume that the time-dependent behavior of the tissue is dominated by both, the isotropic viscoelastic deformation of the matrix and the anisotropic viscoelastic deformation of the fibers. Thus, the equilibrium response of the tissue, represented in Figure 3A by the equilibrium spring, includes isotropic and anisotropic elastic behavior. The nonequilibrium response of the isotropic matrix is represented by the first Maxwell element in Figure 3A and the nonequilibrium response of the anisotropic part (fibers) is represented by the second Maxwell element.

2.4.1 Kinematic assumptions

The kinematics of the nonequilibrium process of the whole tissue is described by the multiplicative split of the deformation gradient F:=∂χ(X,t)/∂X (χ(X,t) motion, X position vector w.r.t. undeformed configuration) into elastic and inelastic (viscous) parts (see Reese and Govindjee, 1998; Nguyen et al., 2007),

$\begin{array}{ll}\mathbf{\text{F}}={\mathbf{\text{F}}}_e^{iso}{\mathbf{\text{F}}}_i^{iso}={\mathbf{\text{F}}}_e^{ani}{\mathbf{\text{F}}}_i^{ani}.& (3)\end{array}$

The viscous contributions, ${\mathbf{\text{F}}}_i^{iso}$ and ${\mathbf{\text{F}}}_i^{ani}$, define a mapping from the undeformed (reference) configuration to the intermediate configuration of the isotropic part ${\tilde{\Omega }}^{iso}$ for the matrix and the intermediate configuration of the anisotropic part ${\tilde{\mathrm{\Omega }}}^{ani}$ for the fibers, respectively (Figure 3B). The elastic components ${\mathbf{\text{F}}}_e^{iso}$ and ${\mathbf{\text{F}}}_e^{ani}$ are then the complementary mappings from ${\tilde{\Omega }}^{iso}$ and ${\tilde{\Omega }}^{ani}$, respectively, to the current configuration Ω. It is assumed that the fibers and matrix deform according to the right Cauchy–Green tensor:

$\begin{array}{ll}\mathbf{\text{C}}:={\mathbf{\text{F}}}^{\text{T}}\mathbf{\text{F}}.& (4)\end{array}$

Analogously we can define the elastic deformation tensors for the isotropic and anisotropic parts

$\begin{array}{ll}{\mathbf{\text{C}}}_e^{iso}:={\mathbf{\text{F}}}_e^{is{o}^T}{\mathbf{\text{F}}}_e^{iso}={\mathbf{\text{F}}}_i^{is{o}^{-T}}\mathbf{\text{C}}{\mathbf{\text{F}}}_i^{is{o}^{-1}},& (5)\end{array}$

$\begin{array}{ll}{\mathbf{\text{C}}}_e^{ani}:={\mathbf{\text{F}}}_e^{an{i}^T}{\mathbf{\text{F}}}_e^{ani}={\mathbf{\text{F}}}_i^{an{i}^{-T}}\mathbf{\text{C}}{\mathbf{\text{F}}}_i^{an{i}^{-1}},& (6)\end{array}$

respectively. To model the transversely isotropic material behavior, we introduce the structure tensor M:= N⊗N in the reference configuration, where N is a unit orientation vector pointing in the direction of the fibers. We assume that the orientation vector N maps only into the anisotropic intermediate configuration. Thus, we can define the structure tensor (Reese et al., 2021)

$\begin{array}{ll}{\widetilde{\mathbf{\text{M}}}}^{ani}:={\widetilde{\mathbf{\text{N}}}}^{ani}\otimes {\widetilde{\mathbf{\text{N}}}}^{ani}=\frac{{\mathbf{\text{F}}}_i^{ani}{\mathbf{\text{MF}}}_i^{an{i}^T} }{{\mathbf{\text{C}}}_i^{ani}:\mathbf{\text{M}}}& (7)\end{array}$

in the intermediate configuration.

2.4.2 Helmholtz free energy

We assume the existence of a free Helmholtz energy with the property of being an isotropic function of the form $\text{Ψ}\left(\mathbf{\text{C}},\mathbf{\text{M}},{\mathbf{\text{C}}}_e^{iso},{\mathbf{\text{C}}}_e^{ani},{\widetilde{\mathbf{\text{M}}}}^{ani}\right)$ (see Boehler, 1977; Spencer, 1982; Reese et al., 2001; Reese, 2003). The proposed rheological model (Figure 3A) motivates the additive split of this scalar potential into three parts (see Reese and Govindjee, 1998; Nguyen et al., 2007):

$\begin{array}{ll}\text{Ψ}={\text{Ψ}}_{eq}\left(\mathbf{\text{C}},\mathbf{\text{M}}\right)+{\text{Ψ}}_{neq}^{iso}\left({\mathbf{\text{C}}}_e^{iso}\right)+{\text{Ψ}}_{neq}^{ani}\left({\mathbf{\text{C}}}_e^{ani},{\widetilde{\mathbf{\text{M}}}}^{ani}\right).& (8)\end{array}$

Here, Ψ_eq(C,M) represents the equilibrium (rate-independent) response of the matrix and the fibers together, the term ${\text{Ψ}}_{\mathrm{neq}}^{iso}\left({\mathbf{\text{C}}}_e^{iso}\right)$ the nonequilibrium (time-dependent) response of the matrix, and the term ${\text{Ψ}}_{neq}^{ani}\left({\mathbf{\text{C}}}_e^{ani},{\widetilde{\mathbf{\text{M}}}}^{ani}\right)$ the nonequilibrium response of the fibers. We additionally decompose the equilibrium part of the free energy further into isotropic (${\text{Ψ}}_{eq}^{iso}\left(\mathbf{\text{C}}\right)$) and anisotropic parts (${\text{Ψ}}_{eq}^{ani}(\mathbf{\text{C}},\mathbf{\text{M}})$). This leads to

$\begin{array}{ll}\text{Ψ}={\text{Ψ}}_{eq}^{iso}\left(\mathbf{\text{C}}\right)+{\text{Ψ}}_{eq}^{ani}\left(\mathbf{\text{C}},\mathbf{\text{M}}\right)+{\text{Ψ}}_{neq}^{iso}\left({\mathbf{\text{C}}}_e^{iso}\right)+{\text{Ψ}}_{neq}^{ani}\left({\mathbf{\text{C}}}_e^{ani},{\widetilde{\mathbf{\text{M}}}}^{ani}\right).& (9)\end{array}$

The stress response of the matrix is assumed to be nearly incompressible. To model such behavior, the following Neo-Hookean free energy density for the matrix is considered:

$\begin{array}{ll}{\text{Ψ}}_{eq}^{iso}\left(\mathbf{\text{C}}\right)=\frac{{\mu }_{eq}^{iso}}{2}\left(\text{tr} \mathbf{\text{C}}-3-\ln (\det \mathbf{\text{C}})\right)+\frac{{\Lambda }_{eq}^{iso}}{4}\left(\det \mathbf{\text{C}}-1-\ln (\det \mathbf{\text{C}})\right).& (10)\end{array}$

For the nonequilibrium contribution of the free energy for the matrix the energy

$\begin{array}{ll}{\text{Ψ}}_{neq}^{iso}\left({\mathbf{\text{C}}}_e^{iso}\right)=\frac{{\mu }_{neq}^{iso}}{2}\left(\text{tr }{\mathbf{\text{C}}}_e^{iso}-3-\ln (\det {\mathbf{\text{C}}}_e^{iso})\right)+\frac{{\Lambda }_{neq}^{iso}}{4}\left(\det {\mathbf{\text{C}}}_e^{iso}-1-\ln (\det {\mathbf{\text{C}}}_e^{iso})\right)& (11)\end{array}$

is assumed. For the equilibrium contribution of the fibers we choose (see Holzapfel et al., 2000)

$\begin{array}{ll}{\text{Ψ}}_{eq}^{ani}\left(\mathbf{\text{C}},\mathbf{\text{M}}\right)=\frac{K_{1,eq}^{ani}}{2K_{2,eq}^{ani}}\left(\exp \left[K_{2,eq}^{ani}{\left(\text{tr }(\mathbf{\text{CM}})-1\right)}^2\right]-1\right),& (12)\end{array}$

and for the nonequilibrium part we choose analogously

$\begin{array}{ll}{\text{Ψ}}_{neq}^{ani}\left({\mathbf{\text{C}}}_e^{ani},{\widetilde{\mathbf{\text{M}}}}^{ani}\right)=\frac{K_{1,neq}^{ani}}{2K_{2,neq}^{ani}}\left(\exp \left[K_{2,neq}^{ani}{\left(\text{tr }({\mathbf{\text{C}}}_e^{ani}{\widetilde{\mathbf{\text{M}}}}^{ani})-1\right)}^2\right]-1\right).& (13)\end{array}$

The quantities ${\mu }_{eq}^{iso}$, ${\text{Λ}}_{eq}^{iso}$, ${\mu }_{neq}^{iso}$, ${\text{Λ}}_{neq}^{iso}$, $K_{1,eq}^{ani}$, $K_{2,eq}^{ani}$, $K_{1,neq}^{ani}$, and $K_{2,neq}^{ani}$ are material constants which have to be fitted to the experimental results presented in section 3.2.

2.4.3 Dissipation inequality

The constitutive equations are needed to determine the global system of equations. For thermodynamic consistency, the internal dissipation inequality

$\begin{array}{ll}\frac{1}{2}\mathbf{\text{S}}:\dot{\mathbf{\text{C}}}-\dot{\text{Ψ}}\ge 0& (14)\end{array}$

has to be fulfilled. Here, S=F⁻¹P is is the second Piola-Kirchhoff stress tensor. Inserting the Helmholtz free energy (9) into the Clausius-Duhem inequality (14) and considering the fact that M is constant, yields

$\begin{array}{ll}\frac{1}{2}\mathbf{\text{S}}:\dot{\mathbf{\text{C}}}-\frac{\partial \text{Ψ}}{\partial \mathbf{\text{C}}}:\dot{\mathbf{\text{C}}}-\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{iso}}:{\dot{\mathbf{\text{C}}}}_e^{iso}-\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{ani}}:{\dot{\mathbf{\text{C}}}}_e^{ani}-\frac{\partial \text{Ψ}}{\partial {\widetilde{\mathbf{\text{M}}}}^{ani}}:{\dot{\widetilde{\mathbf{\text{M}}}}}^{ani}\ge 0.& (15)\end{array}$

Further, we define the rate tensors ${\dot{\mathbf{\text{C}}}}_e^{iso}$, ${\dot{\mathbf{\text{C}}}}_e^{ani}$ and the time derivative of ${\widetilde{\mathbf{\text{M}}}}^{ani}$ (Reese, 2003; Reese and Christ, 2008) by means of

$\begin{array}{ll}{\dot{\mathbf{\text{C}}}}_e^{iso}=-{\mathbf{l}}_i^{is{o}^T}{\mathbf{\text{C}}}_e^{iso}+{\mathbf{\text{F}}}_i^{is{o}^{-T}}\dot{\mathbf{\text{C}}}{\mathbf{\text{F}}}_i^{is{o}^{-1}}-{\mathbf{\text{C}}}_e^{iso}{\mathbf{l}}_i^{iso},& (16)\end{array}$

$\begin{array}{ll}{\dot{\mathbf{\text{C}}}}_e^{ani}=-{\mathbf{l}}_i^{an{i}^T}{\mathbf{\text{C}}}_e^{ani}+{\mathbf{\text{F}}}_i^{an{i}^{-T}}\dot{\mathbf{\text{C}}}{\mathbf{\text{F}}}_i^{an{i}^{-1}}-{\mathbf{\text{C}}}_e^{ani}{\mathbf{l}}_i^{ani},& (17)\end{array}$

$\begin{array}{ll}{\dot{\widetilde{\mathbf{\text{M}}}}}^{ani}={\mathbf{l}}_i^{ani}{\widetilde{\mathbf{\text{M}}}}^{ani}+{\widetilde{\mathbf{\text{M}}}}^{ani}{\mathbf{l}}_i^{an{i}^T}-2\left({\mathbf{\text{d}}}_i^{ani}:{\widetilde{\mathbf{\text{M}}}}^{ani}\right){\widetilde{\mathbf{\text{M}}}}^{ani}.& (18)\end{array}$

In the latter equations, the definitions ${\mathbf{l}}_i^{iso}:={\dot{\mathbf{\text{F}}}}_i^{iso}{\mathbf{\text{F}}}_i^{is{o}^{-1}}$, ${\mathbf{l}}_i^{ani}:={\dot{\mathbf{\text{F}}}}_i^{ani}{\mathbf{\text{F}}}_i^{an{i}^{-1}}$ have been used. The tensor ${\mathbf{\text{d}}}_i^{ani}$ is the symmetric part of ${\mathbf{l}}_i^{ani}$. Exploiting the symmetry of $\partial \text{Ψ}/\partial {\mathbf{\text{C}}}_e^{iso}$ and $\partial \text{Ψ}/\partial {\mathbf{\text{C}}}_e^{ani}$, the Clausius-Duhem inequality (15) can be rewritten to read

$\begin{array}{ll}\left(\mathbf{\text{S}}-2\frac{\partial \mathrm{\Psi }}{\partial \mathbf{\text{C}}}-2{\mathbf{\text{F}}}_i^{is{o}^{-1}}\frac{\partial \mathrm{\Psi }}{\partial {\mathbf{\text{C}}}_e^{iso}}{\mathbf{\text{F}}}_i^{is{o}^{-T}}-2{\mathbf{\text{F}}}_i^{an{i}^{-1}}\frac{\partial \mathrm{\Psi }}{\partial {\mathbf{\text{C}}}_e^{ani}}{\mathbf{\text{F}}}_i^{an{i}^{-T}}\right):\frac{1}{2}\dot{\mathbf{\text{C}}}+{\text{Σ}}^{iso}:{\mathbf{l}}_i^{iso}+{\mathbf{\Gamma }}^{ani}:{\mathbf{l}}_i^{ani}\ge 0.& (19)\end{array}$

Here, Σ^iso and Γ^ani are symmetric tensors (Svendsen, 2001) defined by

$\begin{array}{ll}{\text{Σ}}^{\mathrm{iso}}:=2{\mathbf{\text{C}}}_e^{iso}\frac{\partial \mathrm{\Psi }}{\partial {\mathbf{\text{C}}}_e^{iso}},& (20)\end{array}$

$\begin{array}{ll}{\mathbf{\Gamma }}^{ani}:=2{\mathbf{\text{C}}}_e^{ani}\frac{\partial \mathrm{\Psi }}{\partial {\mathbf{\text{C}}}_e^{ani}}-2\frac{\partial \mathrm{\Psi }}{\partial {\widetilde{\mathbf{\text{M}}}}^{ani}}{\widetilde{\mathbf{\text{M}}}}^{ani}+2\frac{\partial \mathrm{\Psi }}{\partial {\widetilde{\mathbf{\text{M}}}}^{ani}}:\left({\widetilde{\mathbf{\text{M}}}}^{ani}\otimes {\widetilde{\mathbf{\text{M}}}}^{ani}\right).& (21)\end{array}$

2.4.4 Constitutive equations

Following the Coleman-Noll procedure (see Coleman and Noll, 1961) the inequality (19) will be sufficiently satisfied if the second Piola-Kirchhoff stress tensor S is given by the relation

$\begin{array}{ll}\mathbf{\text{S}}=2\frac{\partial \text{Ψ}}{\partial \mathbf{\text{C}}}+2{\mathbf{\text{F}}}_i^{is{o}^{-1}}\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{iso}}{\mathbf{\text{F}}}_i^{is{o}^{-T}}+2{\mathbf{\text{F}}}_i^{an{i}^{-1}}\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{ani}}{\mathbf{\text{F}}}_i^{an{i}^{-T}}.& (22)\end{array}$

Further, the inequality ${\text{Σ}}^{iso}:{\mathbf{\text{d}}}_i^{iso}+{\mathbf{\Gamma }}^{ani}:{\mathbf{\text{d}}}_i^{ani}\ge 0$ has to hold. This can be achieved by making the following choices for the evolution equations:

$\begin{array}{ll}{\mathbf{\text{d}}}_i^{iso}=\frac{1}{{\tau }^{iso}}\frac{\partial g^{iso}\left({\Sigma }^{iso}\right)}{\partial {\Sigma }^{iso}},& (23)\end{array}$

$\begin{array}{ll}{\mathbf{\text{d}}}_i^{ani}=\frac{1}{{\tau }^{ani}}\frac{\partial g^{ani}\left({\mathbf{\Gamma }}^{ani}\right)}{\partial {\mathbf{\Gamma }}^{ani}}.& (24)\end{array}$

Here, τ^iso and τ^ani are relaxation times, and the quantities g^iso and g^ani are relaxation potentials. The potential of the isotropic part g^iso can be additively decomposed into isochoric and the volumetric parts according to

$\begin{array}{ll}{g}^{iso}\left({\text{Σ}}^{iso}\right):=\frac{1}{2{\mu }_{neq}^{iso}}\text{tr }\left({\left({\text{dev Σ}}^{iso}\right)}^2\right)+\frac{1}{9{\kappa }_{neq}^{iso}}\text{tr }{({\text{Σ}}^{iso})}^2,& (25)\end{array}$

where ${\kappa }_{neq}^{iso}:={\text{Λ}}_{neq}^{iso}+2/3{\mu }_{neq}^{iso}$ represents the bulk modulus. The relaxation potential of the anisotropic part g^ani describes the stresses in the direction of the fibers and is defined by

$\begin{array}{ll}{g}^{ani}\left({\mathbf{\Gamma }}^{ani}\right):=\frac{1}{2K_{1,neq}^{ani}}\text{tr }{\left({\mathbf{\Gamma }}^{ani}\text{sym }({\mathbf{\text{C}}}_e^{ani}{\widetilde{\mathbf{\text{M}}}}^{ani})\right)}^2.& (26)\end{array}$

To represent all the constitutive equations in the reference configuration, the evolution equations (23) and (24) need to be pulled back to the reference configuration. A detailed description of this procedure is beyond the scope of this work – the interested reader is referred to e.g. (Marsden and Hughes, 1994; Reese and Christ, 2008). The pull-back leads to the evolution equations

$\begin{array}{ll}\begin{array}{cccc}{\dot{\mathbf{\text{C}}}}_i^{iso}=2& {\mathbf{\text{F}}}_i^{is{o}^T}& {\mathbf{\text{d}}}_i^{iso}& {\mathbf{\text{F}}}_i^{iso}\end{array}=\frac{1}{{\tau }^{iso}}{\mathbf{\text{f}}}^{iso},& (27)\end{array}$

$\begin{array}{ll}\begin{array}{cccc}{\dot{\mathbf{\text{C}}}}_i^{ani}=2& {\mathbf{\text{F}}}_i^{an{i}^T}& {\mathbf{\text{d}}}_i^{ani}& {\mathbf{\text{F}}}_i^{ani}\end{array}=\frac{1}{{\tau }^{ani}}{\mathbf{\text{f}}}^{ani},& (28)\end{array}$

where f^iso and f^ani are defined as

$\begin{array}{ll}{\mathbf{\text{f}}}^{iso}:=2{\mathbf{\text{F}}}_i^{is{o}^T}\frac{\partial g^{iso}\left({\Sigma }^{iso}\right)}{\partial {\Sigma }^{iso}}{\mathbf{\text{F}}}_i^{iso},& (29)\end{array}$

$\begin{array}{ll}{\mathbf{\text{f}}}^{ani}:=2{\mathbf{\text{F}}}_i^{an{i}^T}\frac{\partial g^{ani}\left({\mathbf{\Gamma }}^{ani}\right)}{\partial {\mathbf{\Gamma }}^{ani}}{\mathbf{\text{F}}}_i^{ani}.& (30)\end{array}$

2.4.5 Numerical implementation

The presented model has been implemented into the software Matlab (The MathWorks, Inc., Massachusetts, USA). Software package AceGen has been employed to automatically generate source code for the computation of the tangent operators (see e.g. Korelc, 2002; Korelc, 2009; Korelc and Stupkiewicz, 2014).

In summary, at the Gauss point level, the following equations have to be solved:

$\begin{array}{ll}\mathbf{\text{S}}=2\frac{\partial \text{Ψ}}{\partial \mathbf{\text{C}}}+2{\mathbf{\text{F}}}_i^{is{o}^{-1}}\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{iso}}{\mathbf{\text{F}}}_i^{is{o}^{-T}}+2{\mathbf{\text{F}}}_i^{an{i}^{-1}}\frac{\partial \text{Ψ}}{\partial {\mathbf{\text{C}}}_e^{ani}}{\mathbf{\text{F}}}_i^{an{i}^{-T}},& (31)\end{array}$

$\begin{array}{ll}{\dot{\mathbf{\text{C}}}}_i^{iso}=\frac{1}{{\tau }^{iso}}{\mathbf{\text{f}}}^{iso},& (32)\end{array}$

$\begin{array}{ll}{\dot{\mathbf{\text{C}}}}_i^{ani}=\frac{1}{{\tau }^{ani}}{\mathbf{\text{f}}}^{ani}.& (33)\end{array}$

For the implicit time integration of the evolution equations (32) and (33), the exponential mapping algorithm has been used (see e.g. Reese and Govindjee, 1998; Vladimirov et al., 2008; Lamm et al., 2022). The final form of the update formulas read:

$\begin{array}{ll}{({\mathbf{\text{C}}}_i^{is{o}^{-1}})}_n={\mathbf{\text{U}}}_i^{is{o}^{-1}} \text{exp}\left(\begin{array}{cccc}\frac{\Delta t}{{\tau }^{iso}}& {\mathbf{\text{U}}}_i^{is{o}^{-1}}& {\mathbf{\text{f}}}^{iso}& {\mathbf{\text{U}}}_i^{is{o}^{-1}}\end{array}\right){\mathbf{\text{U}}}_i^{is{o}^{-1}}& (34)\end{array}$

$\begin{array}{ll}{({\mathbf{\text{C}}}_i^{an{i}^{-1}})}_n={\mathbf{\text{U}}}_i^{an{i}^{-1}} \text{exp}\left(\begin{array}{cccc}\frac{\Delta t}{{\tau }^{ani}}& {\mathbf{\text{U}}}_i^{an{i}^{-1}}& {\mathbf{\text{f}}}^{ani}& {\mathbf{\text{U}}}_i^{an{i}^{-1}}\end{array}\right){\mathbf{\text{U}}}_i^{an{i}^{-1}}& (35)\end{array}$

where ${\mathbf{\text{U}}}_i^{iso}:=\sqrt{{\mathbf{\text{C}}}_i^{iso}}$ and ${\mathbf{\text{U}}}_i^{ani}:=\sqrt{{\mathbf{\text{C}}}_i^{ani}}$ are the inelastic right stretch tensors, and Δt: =t_n+1−t_n is the time interval.

3 Results

3.1 Anatomy and morphology

The analyzed leaves of S. japonica (Figure 4A) show the triangular-ovate shape described for this species and an average lamina length of 14.9±1.6 cm and width of 12.3±1.3 cm. The average lamina area is 120.2±23 cm². The petiole is 2.5±0.2 mm in diameter and 11.8±2.7 cm long. The petiole is inserted eccentric with an average distance of the petiole insertion point from the basal leaf margin of 3.0±0.5 cm. The lamina is palmately veined. From the petiole insertion point, the main vein runs towards the leaf apex. Larger veins run towards the leaf sides, while smaller veins supply the base of the leaf, together building the scaffolding of the lamina and surrounding the intercostal areas.

FIGURE 4

Figure 4 Morphology and anatomy of the Stephania japonica leaf. General leaf shape (A), transverse (B) and longitudinal section (C) of the petiole-lamina transition zone, transverse section of the petiole (D) and topview (E) and sideview (F) of the fiber segmentation from CT data in the transition zone. (for B, D, E: abaxial side at the bottom, staining: astrablue/safranin, scales: 2 cm (A) and 100 µm (B-D).

The petiole shows a more or less circular shape in cross-section (Figure 4D). The vascular bundles are arranged in a circle with the xylem positioned towards the center and the phloem towards the periphery. The bundles on the abaxial side are larger than those on the adaxial side. Sclerenchymatous tissue accompanies the vascular bundles on the phloem side. This tissue can be arranged as sclerenchymatous caps or in a closed ring. The percentage of lignified tissue in the petiole averages at 7.2±3.4 % . Statistical analysis (Mann-Whitney U test) shows, that the basal part of the petiole has significant higher percentages than the apical part ( 10.6±2.4 % vs. 3.5±1.0 %).

The petiole-lamina transition zone shows a ring-like structure (Figures 4B, C). The xylem bundles run parallel in the petiole and from there into the transition zone where each bundle cross-links and merges with at least the neighboring bundle to form a new bundle that runs into a vein of the lamina. The three larger bundles from the abaxial side of the petiole run into the transition zone and subsequently, form the three largest veins in the lamina. The cross-linking in the transition zone is achieved by single xylem elements from one bundle connecting to xylem elements of the neighboring bundles (Figures 4B, E). The additional sclerenchymatous fibers associated with the vascular bundles run into the transition zone and into the veins of the lamina accompanying the vascular bundles (Figure 4F).

3.2 Mechanical tests

The stress-time curves of cyclic uniaxial tensile experiments are plotted in Figures 5, 7, 8 for the petiole, venation and intercostal area, respectively. Furthermore, stress-time curves from uniaxial compression tests of the petiole are plotted in Figure 6.

FIGURE 5

Figure 5 Uniaxial cyclic tensile stress-relaxation test of petiole.

FIGURE 6

Figure 6 Uniaxial cyclic compression stress-relaxation test of petiole.

FIGURE 7

Figure 7 Uniaxial cyclic tensile stress-relaxation test of venation.

FIGURE 8

Figure 8 Uniaxial cyclic tensile stress-relaxation test of intercostal areas.

Results show a clear viscoelastic behavior of the petiole, venation and intercostal area. At higher engineering strain rates tissues respond with increased stiffness (Figures 5–8). After the loading the stress relaxes to the equilibrium state. In the results we also see, that equilibrium states are independent of the engineering strain rate. The equilibrium states at fixed engineering strains 0.5 %, 1.0 % and 2.0 % are listed in Table 1. Petiole and venation show a similarly stiff response. Intercostal area responds to tension with a ten times lower stiffness than the petiole and venation (see Table 1). The equilibrium strain-stress response is nonlinear (see Table 1).

TABLE 1

Table 1 Equilibrium state of each fixed strain.

All tests (Figures 5–8) show, that in the second and third cycle when the engineering strain has been removed, the material does not recover to the initial state.

3.3 Simulation

The simulation of the uniaxial cyclic tensile tests of the petiole, venation, intercostal area, and the uniaxial compression test of the petiole of S. japonica were performed at a material point level. The implemented material model includes 10 parameters introduced in Section 2.4, which have to be identified beforehand. This can be achieved in different ways. We decided for the curve fitting approach via optimization, where we used an algorithm, which finds the parameters by minimizing the squared difference between the experimental data and the simulation results.

The parameters of the model were fitted only to the averaged stress-time curves of the tests where the loading with strain rate of 10 % min^-1 was applied (Figure 2). Furthermore, the fits are restricted to uniaxial tension, and by petiole also to uniaxial compression. All the identified parameters are listed in Table 2.

TABLE 2

Table 2 Fitted parameters of the transversely isotropic viscoelastic material model.

The Lamé parameters ${\mu }_{eq}^{iso}$ and ${\text{Λ}}_{eq}^{iso}$ can be expressed as

$\begin{array}{ll}{\mu }_{eq}^{iso}=\frac{E_{eq}^{iso}}{2\left(1+{\nu }_{eq}^{iso}\right)},& (36)\end{array}$

$\begin{array}{ll}{\text{Λ}}_{eq}^{iso}=\frac{{\nu }_{eq}^{iso}{E}_{eq}^{iso}}{\left(1-2{\nu }_{eq}^{iso}\right)\left(1+{\nu }_{eq}^{iso}\right)},& (37)\end{array}$

where $E_{eq}^{iso}$ is a Young’s module of the isotropic part and ${\nu }_{eq}^{iso}$ is the corresponding Poisson’s ratio. If we assume, that the Poisson’s ratios of the isotropic part (matrix) of all three structural components of the leaf of S. japonica (petiole, venation and intercostal area) are the same, we can compute its average value from the fitted parameters to read ${\overline{v}}_{iso}$ = 0.480. This result can be used in further computations.

The fibers in soft tissues usually cannot carry compressive loadings. Due to this fact, a clear tension-compression asymmetry is observed in the results. The obtained petiole parameters for tension differ from those for compression. If we compute the Young’s modulus of the matrix part $E_{eq}^{iso}$ for the petiole using eq. (36) and (37), we obtain an averaged value $E_{eq}^{iso}=16.029 \text{MPa}$ for tension and $E_{eq}^{iso}=13.496 \text{MPa}$ for compression. The $E_{eq}^{iso}$ values for tension for venation and intercostal areas averaged 10.248 MPa and 4.088 MPa, respectively. Hence, petiole possesses the highest and intercostal areas the lowest tensile stiffness of the isotropic part of the tissue (matrix).

The anisotropic constants $K_{1,eq}^{ani}$ determine the equilibrium stiffness of the anisotropic part of the tissues (fibers). Our results show that venation possesses a higher $K_{1,eq}^{ani}$ value than petiole (see Table 2). This corresponds to the fact that venation consists mainly of xylem fibers (see Sec. 3.1). The $K_{1,eq}^{ani}$ value of intercostal area is very low, which also corresponds to its anatomy (see Sec. 3.1). Thus, the anisotropic part of the intercostal area can be neglected in future computations. $K_{1,eq}^{ani}$ constants determine the nonequilibrium stiffness of the anisotropic part. The higher they are, the higher the extra (nonequilibrium) stress at some constant strain rate. Accordingly, the petiole and the venation show a similar amount of relaxed stress.

Viscous effects in the material are distributed between the isotropic and the anisotropic parts of the tissue. The parameters 1/τ^iso and 1/τ^ani are directly related to the viscosity of the material which influences the nonequilibrium stress. From our results, it is evident that the anisotropic parts (fibers) have a lower relaxation time than the isotropic parts. This means that fibers relax faster than the matrix.

We simulated four different cyclic stress-relaxation experiments, where the 10 % min⁻¹ strain rate has been applied for loading and unloading. The strain was applied in the same way as it was applied in the experiments - three cycles of loading, holding and unloading (see Figure 2). The simulated results are plotted in the Figures 9–12. Black curves represent averaged experimental results, red curves represent simulation results and the green regions indicate a standard deviation from the averaged experimental results.

FIGURE 9

Figure 9 Simulation of the cyclic uniaxial tensile stress-relaxation test of the petiole (red curve). Green region indicates a standard deviation from the averaged experimental result (black curve).

FIGURE 10

Figure 10 Simulation of the cyclic uniaxial compression stress-relaxation test of the petiole (red curve). Green region indicates a standard deviation from the averaged experimental result (black curve).

FIGURE 11

Figure 11 Simulation of the cyclic uniaxial tensile stress-relaxation test of the venation (red curve). Green region indicates a standard deviation from the averaged experimental result (black curve).

FIGURE 12

Figure 12 Simulation of the cyclic uniaxial tensile stress-relaxation test of the intecostal area (red curve). Green region indicates a standard deviation from the averaged experimental result (black curve).

The simulated results of the tensile test of petiole stays well within the experimental error. Simulations of tensile tests of venation and intercostal area still show a good match but deviate slightly from the experimental data, especially after the unloading. The reason for this is that after each unloading, viscoelastic material should theoretically relax to a zero line. This different behavior of the samples has been already observed in Section 3.2. The reason for this is explained in Section 4.

To compare the experimental and simulation results, we define a normalized root-mean-square error (NRMSE )

$\begin{array}{ll}NRMSE:=\frac{1}{P_{11}^{exp}{ }_{max}-P_{11}^{exp}{ }_{min}}\sqrt{\frac{1}{n}{\displaystyle \sum _{i=1}^n{\left([P_{11}^{sim}{]}_i-[P_{11}^{\exp }{]}_i\right)}^2}}& (38)\end{array}$

where $P_{11}^{exp}{ }_{max}-P_{11}^{exp}{ }_{min}$ is the range of the measured data ($P_{11}^{exp}{ }_{max}$ is the maximum value, $P_{11}^{exp} min$ is the minimum value), n is the number of data points, ${[P_{11}^{exp}]}_i$ is the averaged experimental data point for the first component of the first Piola-Kirchhoff stress tensor at i-th time, and ${[P_{11}^{sim}]}_i$ is the simulation result at the corresponding time. In Table 3, the NRMS errors of the simulated experiments (Figures 9–12) are listed. In the computation of the error only the result from the beginning of the first to the end of the third cycle were considered. The error stays below 5 % for all four simulations. This result indicates, that the model provides a good fit of the cyclic uniaxial stress-relaxation tests.

TABLE 3

Table 3 Normalized root-mean-square errors comparing the simulation and experimental data for cyclic uniaxial stress-relaxation tests.

4 Discussion

The main goal of this study was to investigate the mechanical behavior of peltate leaf of Stephania japonica. First, the mechanical tests on petiole, venation and intercostal areas were performed. To arrive at a better understanding of the experimental results a transversely isotropic viscoelastic model for finite deformation was developed and tested. Additionally, the anatomy of petiole and petiole-lamina transition zone was analyzed.

4.1 Mechanical tests

The available measurement technique is unfortunately limited in the way that it is not possible to measure Poisson’s ratio. The material model proposed in this work includes ten parameters. The identification of those requires several independent experiments and a challenging testing procedure, which is a research topic in its own. In the future, digital image correlation may help to calibrate the model in a more accurate way, since the whole deformation field can be obtained.

Since only one individuum of Stephania japonica was available, the total number of leaves available for testing was limited. This resulted into a relatively small number of samples tested. Furthermore, the samples came from different leaves causing a higher deviation in the experimental results which is not unusual for plant samples. Nevertheless, these results do show the fundamental mechanical behavior of the plant parts. This is sufficient for the calibration of the proposed model. To obtain statistically more accurate results, the sample size should be increased which is left open for future work.

The differences in structure and composition of samples, such as amount of lignin and the degree of lignification, should also be considered in experiments. Since the tests were performed in the first place with the aim to identify the major mechanical mechanisms of plants, only the averaged behavior of plant leaves was taken into account.

4.2 Mechanical behavior of the structural elements of the leaf

Mechanical experiments reveal that all parts of the leaf show a rate-dependent behavior. This behavior is typical for plant tissues, especially parenchymatous tissue (Pitt and Chen, 1983; Niklas, 1992) which is, in different proportions, present in all structural elements of the leaf. When parenchyma tissue is stressed, the stresses are redistributed to the cell walls throughout the tissue. The cell walls flatten in reaction to the applied load and fluids from the cell cross the plasma membrane. With higher strain-rates, less time is available for the fluids to diffuse which results in a higher elastic modulus (Niklas, 1992).

Additionally, experimental data (see Figures 5–8) show that all samples in the second and third cycle, after the strain is removed, do not recover to the initial state (zero stress). This means that all leaf components change their mechanical properties around 1% of strain. This phenomenon can be a result of different processes within the samples. One possibility is occurring damage in the tissues, which is a result of the deformed cell walls due to fluid loss in hydrostatic tissues when loads are applied (as mentioned above), cell–cell debonding or the rupture of cell membrane and walls (Pitt and Chen, 1983; Niklas, 1992; Zdunek and Umeda, 2006). Another reason for this change in mechanical properties could be the change of turgor in the cells as a consequence of water loss during testing since testing times exceeded 15 min and especially in venation and intercostal areas, larger surfaces from cutting were exposed to air. Water loss has a significant impact on the properties of plant tissues including hydrostatic and fibrous tissues. While the stiffness of hydrostatic tissues decreases with water loss, plant fibers usually exhibit a higher stiffness when dried (Niklas, 1992; Niklas, 1999; Faisal et al., 2010). For a clearer understanding of this phenomenon more tests have to be conducted.

Assuming plastic deformation, we can still formulate some statements. Petiole and venation can resist quite high tensile loads – they can withstand around 1 MPa of inner stress at 0.5 % strain (see Table 1). As the petiole shows several large xylem fiber bundles, as well as additional sclerenchyma fibers (Figure 4) resulting in an average fiber percentage of 7.2 % and the venation consists mostly of vascular tissue, including the lignified xylem fibers, resistance to higher tensile loads is expected. The intercostal areas have a more than ten times lower resistance to tensile stress. This can be attributed to their anatomical features. The intercostal areas of the leaf are mainly composed of spongy parenchymatous mesophyll. Since the stiffness of parenchymatous tissue also depends on the packaging density of the cells and cell-to-cell connections, tissues with a high amount of air-filled spaces like the mesophyll, show only low resistance to mechanical stresses (Gao and Pitt, 1991; Niklas, 1992). The elastic range of the leaf components is below 0.5 % of tensile strain (see Table 1). This rather low elastic range, which is typical for solid materials, is also a known phenomenon in hydrostatic plant tissues (Niklas, 1992).

The petiole can withstand around 0.2 MPa pressure at 0.5 % compression strain. In plants, compression strains are mostly carried by tissues that are stronger in compression (e.g., parenchyma) than fibers. Thus, in the petiole of S. japonica the tissues bearing these loads are the parenchyma and epidermis.

4.3 The petiole-lamina transition zone of the peltate leaf of Stephania japonica

The petiole and petiole-lamina transition zone of S. japonica showed the same anatomy and fiber arrangement like S. venosa and S. delavayi (Wunnenberg et al., 2021). Especially the ring-like arrangement of the vascular bundles in the transition zone could be shown in high resolution in this study (Figure 4). The segmented CT data allows for understanding of the transition zone on a single-fiber level. Studies have shown the importance of the transition zone in load dissipation and bearing part of the torsional load affecting the leaf (Sacher et al., 2019; Langer et al., 2022). The stability of the transition zone and its ability to dissipate load from lamina and petiole could be achieved by the cross-linking of several xylem elements that branch off of the petiole vascular bundles and merge with elements from the neighboring bundle forming the vascular bundles of the venation. The additional sclerenchyma fibers that run from the periphery of the petiole into the periphery of the leaf veins allow for extra resistance to mechanical loads in the transition zone. Further experiments and modelling could shed further light on the mechanical properties of the petiole-lamina transition zone.

4.4 The material model and simulations

We developed a constitutive framework that is capable of predicting an anisotropic viscous behavior (strain-stress curves) of the leaf of S. japonica. The fitted parameters are physically meaningful. They give a first insight into plant tissue properties which correspond with experimental results. The simulated results agree well with experimental data. The NRMS errors comparing the experimental and simulation results were less than 5 % for all of the tests.

One of the main goals of this work was to develop and test the proposed continuum mechanical material model. Since the calibration of the model is a highly time consuming process, we identified the material parameters only for the strain rate 10 % min⁻¹. This set of parameters is sufficient to test the model. However, it may not be the unique parameter set, which would allow the model to accurately predict any other possible loading rate. To achieve this goal, further calibrations, also in combination with further experiments (e.g. torsion, bending), need to be performed. However, this is beyond the scope of the present contribution.

With this model more complex stress states, experimentally not achievable, can be simulated. Such simulations are very suitable to provide an even better understanding of the behavior of peltate leaves in general and the peltate leaf of S. japonica. The model can be easily implemented into a finite element environment. There are some features that can be considered to improve this model, for instance non-constant fiber orientations or temperature dependence of parameters. With a possible extension of this model, not only the response to mechanical loads on petiole, venation and intercostal areas of the leaf can be simulated, but also the complex petiole-lamina transition zone and its mechanical significance in load dissipation. This stays open for future research.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author contributions

DM and AR reviewed the relevant existing literature, planned the experiments and wrote this article. DM performed all simulations and interpreted the experimental and simulation results. AR performed all experiments. HH worked out the theoretical material model and implemented it into the software Matlab (The MathWorks, Inc., Massachusetts, USA). HH, StR, JWS, TL, CN and SR gave conceptual advice, contributed in the discussion of the results, read the article and gave valuable suggestions for improvement. All authors contributed to the article and approved the submitted version.

Funding

This research was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)-SFB/TRR 280. Project-ID: 417002380.

Acknowledgments

The authors would like to thank the staff of the Botanical Garden Dresden for providing the plant material and Dr. Frank Liebold and Franz Wagner from the Institute of Photogrammetry and Remote Sensing of the Technische Universität Dresden for performing and reconstructing the X-ray tomography scans. All authors gratefully acknowledge the financial support of the research work by the German Research Foundation (DFG, Deutsche Forschungsgemeinschaft) within the transregional Collaborative Research Center SFB/TRR 280, project number 417002380, subprojects A01 and A03. Furthermore, HH and SR acknowledge the financial support through the DFG project PAK 961 “Towards a model based control of biohybrid implant maturation” (RE 1057/45-1, Project number 403471716). SR acknowledges the financial support through DFG projects StentsX “Medikamente freisetzende Koronarstents in stenotischen Arterien: medizinische Untersuchung und numerische Modellierung” (RE 1057/44-1, Project number 395712048) and SPP 2311 “Robuste Kopplung kontinuumsbiomechanischer in silico Modelle für aktive biologische Systeme als Vorstufe klinischer Applikationen - Co-Design von Modellierung, Numerik und Nutzbarkeit” (RE 1057/53-1, Project number 441884911).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Adams, W. W., III, Terashima, I. (2018). The leaf: a platform for performing photosynthesis Vol. 44 (Cham, Switzerland:Springer).

Google Scholar

Bischoff, J. E., Arruda, E. M., Grosh, K. (2004). A rheological network model for the continuum anisotropic and viscoelastic behavior of soft tissue. Biomechanics modeling mechanobiology 3, 56–65. doi: 10.1007/s10237-004-0049-4

CrossRef Full Text | Google Scholar

Boehler, J. (1977). On irreducible representations for isotropic scalar functions. ZAMM-Journal Appl. Mathematics Mechanics/Zeitschrift für Angewandte Mathematik und Mechanik 57, 323–327. doi: 10.1002/zamm.19770570608

CrossRef Full Text | Google Scholar

Coleman, B. D., Noll, W. (1961). Foundations of linear viscoelasticity. Rev. modern Phys. 33, 239. doi: 10.1103/RevModPhys.33.239

CrossRef Full Text | Google Scholar

Dataset POWO (2022) Plants of the World Online, facilitated by the Royal Botanical Gardens, Kew. Available at: https://powo.science.kew.org/.

Google Scholar

Faisal, T. R., Abad, E. M. K., Hristozov, N., Pasini, D. (2010). The impact of tissue morphology, cross-section and turgor pressure on the mechanical properties of the leaf petiole in plants. J. Bionic Eng. 7, S11–S23. doi: 10.1016/S1672-6529(09)60212-2

CrossRef Full Text | Google Scholar

Falk, S., Hertz, C. H., Virgin, H. I. (1958). On the relation between turgor pressure and tissue rigidity. I: Experiments on resonance frequency and tissue rigidity. Physiologia Plantarum 11, 802–817. doi: 10.1111/j.1399-3054.1958.tb08274.x

CrossRef Full Text | Google Scholar

Forman, L. L. (1986). “Menispermaceae,” in Flora Malesiana, vol. 10. (Dodrecht, Boston, Lancaster: Martinus Nijhoff Publishers), 157–253.

Google Scholar

Gao, Q., Pitt, R. E. (1991). Mechanics of parenchyma tissue based on cell orientation and microstructure. Trans. ASAE 34, 232–0238. doi: 10.13031/2013.31651

CrossRef Full Text | Google Scholar

Garcia-Gonzalez, D., Jérusalem, A., Garzon-Hernandez, S., Zaera, R., Arias, A. (2018). A continuum mechanics constitutive framework for transverse isotropic soft tissues. J. Mechanics Phys. Solids 112, 209–224. doi: 10.1016/j.jmps.2017.12.001

CrossRef Full Text | Google Scholar

Gibson, L. J. (2012). The hierarchical structure and mechanics of plant materials. J. R. Soc. Interface 9, 2749–2766. doi: 10.1098/rsif.2012.0341

PubMed Abstract | CrossRef Full Text | Google Scholar

Holzapfel, G. A., Gasser, T. C. (2001). A viscoelastic model for fiber-reinforced composites at finite strains: Continuum basis, computational aspects and applications. Comput. Methods Appl. mechanics Eng. 190, 4379–4403. doi: 10.1016/S0045-7825(00)00323-6

CrossRef Full Text | Google Scholar

Holzapfel, G. A., Gasser, T. C., Ogden, R. W. (2000). A new constitutive framework for arterial wall mechanics and a comparative study of material models. J. elasticity Phys. Sci. solids 61, 1–48. doi: 10.1023/A:1010835316564

CrossRef Full Text | Google Scholar

Korelc, J. (2002). Multi-language and multi-environment generation of nonlinear finite element codes. Eng. Comput. 18, 312–327. doi: 10.1007/s003660200028

CrossRef Full Text | Google Scholar

Korelc, J. (2009). “Automation of primal and sensitivity analysis of transient coupled problems,” in Computational mechanics (Springer), vol. 44, 631–649.doi: 10.1007/s00466-009-0395-2

CrossRef Full Text | Google Scholar

Korelc, J., Stupkiewicz, S. (2014). Closed-form matrix exponential and its application in finite-strain plasticity. Int. J. Numerical Methods Eng. 98, 960–987. doi: 10.1002/nme.4653

CrossRef Full Text | Google Scholar

Lamm, L., Holthusen, H., Brepols, T., Jockenhövel, S., Reese, S. (2022). A macroscopic approach for stress-driven anisotropic growth in bioengineered soft tissues. Biomechanics Modeling Mechanobiology 21, 627–645. doi: 10.1007/s10237-021-01554-1

CrossRef Full Text | Google Scholar

Langer, M., Hegge, E., Speck, T., Speck, O. (2022). Acclimation to wind loads and/or contact stimuli? A biomechanical study of peltate leaves of Pilea peperomioides. J. Exp. Bot. 73, 1236–1252. doi: 10.1093/jxb/erab541

PubMed Abstract | CrossRef Full Text | Google Scholar

Langer, M., Kelbel, M. C., Speck, T., Müller, C., Speck, O. (2021a). Twist-to-Bend ratios and safety factors of petioles having various geometries, sizes and shapes. Front. Plant Sci. 12, 765605–765605. doi: 10.3389/fpls.2021.765605

PubMed Abstract | CrossRef Full Text | Google Scholar

Langer, M., Speck, T., Speck, O. (2021b). Petiole-lamina transition zone: A functionally crucial but often overlooked leaf trait. Plants 10, 774. doi: 10.3390/plants10040774

PubMed Abstract | CrossRef Full Text | Google Scholar

Latorre, M., Montáns, F. J. (2015). Anisotropic finite strain viscoelasticity based on the sidoroff multiplicative decomposition and logarithmic strains. Comput. Mechanics 56, 503–531. doi: 10.1007/s00466-015-1184-8

CrossRef Full Text | Google Scholar

Latorre, M., Montáns, F. J. (2016). Fully anisotropic finite strain viscoelasticity based on a reverse multiplicative decomposition and logarithmic strains. Comput. Structures 163, 56–70. doi: 10.1016/j.compstruc.2015.09.001

CrossRef Full Text | Google Scholar

Liu, H., Holzapfel, G. A., Skallerud, B. H., Prot, V. (2019). Anisotropic finite strain viscoelasticity: Constitutive modeling and finite element implementation. J. Mechanics Phys. Solids 124, 172–188. doi: 10.1016/j.jmps.2018.09.014

CrossRef Full Text | Google Scholar

Marsden, J. E., Hughes, T. J. (1994). Mathematical foundations of elasticity (New York, USA:Courier Corporation).

Google Scholar

Napp-Zinn, K. (1973). Anatomie des Blattes II, Blattanatomie der Angiospermen. 2 edn (Berlin, Stuttgart: Gebrüder Bornträger).

Google Scholar

Nedjar, B. (2007). An anisotropic viscoelastic fibre–matrix model at finite strains: continuum formulation and computational aspects. Comput. Methods Appl. Mechanics Eng. 196, 1745–1756. doi: 10.1016/j.cma.2006.09.009

CrossRef Full Text | Google Scholar

Nguyen, T., Jones, R., Boyce, B. (2007). Modeling the anisotropic finite-deformation viscoelastic behavior of soft fiber-reinforced composites. Int. J. Solids Structures 44, 8366–8389. doi: 10.1016/j.ijsolstr.2007.06.020

CrossRef Full Text | Google Scholar

Nguyen, T., Jones, R., Boyce, B. (2008). A nonlinear anisotropic viscoelastic model for the tensile behavior of the corneal stroma. J. Biomechanical Eng. 130 (4), 041020. doi: 10.1115/1.2947399

CrossRef Full Text | Google Scholar

Niklas, K. J. (1992). Plant biomechanics: an engineering approach to plant form and function (Chicago, USA:University of Chicago).

Google Scholar

Niklas, K. J. (1999). A mechanical perspective on foliage leaf form and function. New Phytol. 143, 19–31. doi: 10.1046/j.1469-8137.1999.00441.x

CrossRef Full Text | Google Scholar

Pitt, R. E., Chen, H. L. (1983). Time-dependent aspects of the strength and rheology of vegetative tissue. Trans. ASAE 26, 1275–1280. doi: 10.13031/2013.34116

CrossRef Full Text | Google Scholar

Reese, S. (2003). Meso-macro modelling of fibre-reinforced rubber-like composites exhibiting large elastoplastic deformation. Int. J. Solids Structures 40, 951–980. doi: 10.1016/S0020-7683(02)00602-9

CrossRef Full Text | Google Scholar

Reese, S., Brepols, T., Fassin, M., Poggenpohl, L., Wulfinghoff, S. (2021). Using structural tensors for inelastic material modeling in the finite strain regime–a novel approach to anisotropic damage. J. Mechanics Phys. Solids 146, 104174. doi: 10.1016/j.jmps.2020.104174

CrossRef Full Text | Google Scholar

Reese, S., Christ, D. (2008). Finite deformation pseudo-elasticity of shape memory alloys–constitutive modelling and finite element implementation. Int. J. Plasticity 24, 455–482. doi: 10.1016/j.ijplas.2007.05.005

CrossRef Full Text | Google Scholar

Reese, S., Govindjee, S. (1998). A theory of finite viscoelasticity and numerical aspects. Int. J. solids structures 35, 3455–3482. doi: 10.1016/S0020-7683(97)00217-5

CrossRef Full Text | Google Scholar

Reese, S., Raible, T., Wriggers, P. (2001). Finite element modelling of orthotropic material behaviour in pneumatic membranes. Int. J. solids structures 38, 9525–9544. doi: 10.1016/S0020-7683(01)00137-8

CrossRef Full Text | Google Scholar

Sacher, M., Lautenschläger, T., Kempe, A., Neinhuis, C. (2019). Umbrella leaves–biomechanics of transition zone from lamina to petiole of peltate leaves. Bioinspiration Biomimetics 14, 046011. doi: 10.1088/1748-3190/ab2411

PubMed Abstract | CrossRef Full Text | Google Scholar

Spencer, A. (1982). “The formulation of constitutive equation for anisotropic solids,” in Mechanical behavior of anisotropic Solids/Comportment méchanique des solides anisotropes (Springer), 3–26. doi: 10.1007/978-94-009-6827-1_1

CrossRef Full Text | Google Scholar

Stier, B., Simon, J.-W., Reese, S. (2015). Comparing experimental results to a numerical meso-scale approach for woven fiber reinforced plastics. Composite structures 122, 553–560. doi: 10.1016/j.compstruct.2014.12.015

CrossRef Full Text | Google Scholar

Stumpf, F. (2021). An accurate and efficient constitutive framework for finite strain viscoelasticity applied to anisotropic soft tissues. Mechanics Materials 161, 104007. doi: 10.1016/j.mechmat.2021.104007

CrossRef Full Text | Google Scholar

Svendsen, B. (2001). On the modelling of anisotropic elastic and inelastic material behaviour at large deformation. Int. J. Solids Structures 38, 9579–9599. doi: 10.1016/S0020-7683(01)00140-8

CrossRef Full Text | Google Scholar

Troll, W. (1932). Morphologie der schildförmigen Blätter. Planta 17, 153–230. doi: 10.1007/BF01909626

CrossRef Full Text | Google Scholar

Vladimirov, I. N., Pietryga, M. P., Reese, S. (2008). On the modelling of non-linear kinematic hardening at finite strains with application to springback–comparison of time integration algorithms. Int. J. Numerical Methods Eng. 75, 1–28. doi: 10.1002/nme.2234

CrossRef Full Text | Google Scholar

Wunnenberg, J., Rjosk, A., Neinhuis, C., Lautenschläger, T. (2021). Strengthening structures in the petiole-lamina junction of peltate leaves. Biomimetics 6, 25. doi: 10.3390/biomimetics6020025

PubMed Abstract | CrossRef Full Text | Google Scholar

Zdunek, A., Umeda, M. (2006). Extension and fracture of cell walls after parenchyma tissue deformation. Biosyst. Eng. 93, 269–278. doi: 10.1016/j.biosystemseng.2005.12.007

CrossRef Full Text | Google Scholar