Growth spin in plant cell wall growth and remodelling

Abstract

A growth spin model is proposed for the modelling of the expansive growth and fibre remodelling of cell wall. The introduction of the growth spin allows to relax the perfectly bonding assumption in the kinematical growth, which provides a more sophisticated and flexible kinematical description for modelling the selective control and flexible regulation of anisotropic growth. Extended hardening laws are proposed for the growth and remodelling, respectively, aiming at further clarification of the dynamic balance between hardening and softening in a growing cell wall. The proposed model may shed a new light into the micro-structural interpretation of the “fictitious” intermediate (growth) configuration in the kinematical growth modelling of soft matter. A case study of the cell wall as a growing cylindrical wall is presented to demonstrate the proposed model. Spencer’s deviatoric stress tensor and its rate format are shown to play a key role in the modelling of the cell wall as a fibre-reinforced soft tissue.

1. Introduction

How cells expand despite the presence of cell walls is a central issue in plant biology [1–6]. The plant cell wall is a composite in which cellulose microfibrils (CMFs) are embedded into an amorphous matrix composed of hemicellulose, pectins, and proteins (see Figure 1). Microfibrils, several nanometres in width and many micrometres long, have a tensile strength similar to that of steel, and are relatively inert and inextensible in growth [7,8]. Neighbouring microfibrils tend to be roughly parallel, giving the cell wall a mat-like appearance and a distinct structural anisotropy similar to a laminated composite.

Figure 1.

A schematic REV of cell wall showing distribution of cellulose microfibrils (CMFs), hemicelluloses, and pectins in a laminated material.

The cell wall is synthesised in an integrative effort: mobile cellulose synthesis complexes produce long, thin, strong, stiff CMFs at the cell surface, while matrix polysaccharides and glycoproteins are deposited to the cell surface via the secretory apparatus [4]. The newly deposited material may be transported through the wall matrix along the wall thickness direction and integrated into the existing microstructure [9,10]. Together, such a process of deposition, transport, and integration creates the volumetric growth of the wall.

However, wall synthesis without mechanical relaxation would only cause wall thickening [11,12]. To grow, plant cells must physically expand their restraining cell walls with the help of turgor pressure (about 5 bar) applied on its inner surface while satisfying wall yielding (or called wall loosening) threshold condition [3,13]. The wall yielding/loosening refers to a molecular modification of the wall network that results in relaxation of wall stress. Wall stress relaxation might result from scission of a stress-bearing crosslink or from sliding of such a crosslink along a scaffold; in both cases, it results in a reduction in wall stress without a substantial change in wall dimensions [2]. Therefore, in addition to the contribution from the volumetric growth, cell wall expansive growth has the contribution from the passively stress-driven irreversible isochoric in-plane expansion, which is analogues to the conventional viscoplasticity and may be considered playing a dominant and leading role in cell wall growth. Nevertheless, wall expansion without the matching volumetric growth inevitably results in wall thinning and eventually ruptures of material. Thus, the matching volumetric growth may be considered playing a supporting role to the wall expansion for maintaining stable wall geometry and strength.

Therefore, the cell wall must satisfy two contradictory requirements. It must be strong enough to resist mechanical stress (which may exceed ${10}^8\mathrm{N}{\mathrm{m}}^{-2}$ (1000 bar)) generated by cell turgor pressure, and at the same time, it must be sufficiently compliant to permit irreversible wall expansion. Although it is widely accepted that the cell wall accommodates these requirements through its composite structure, the earliest molecular depictions of growing cell walls were offered as concept maps of how the composite components might be spatially arranged and connected to one another as hypothetical summaries rather than quantitative or predictive models [5].

By contrast, macroscopic phenomenological models provided both quantitative prediction and qualitative interpretation from the perspective of continuum theory. As a milestone, Lockhart [14] summarised a wide range of experimental data on wall extensibility in a formula, i.e., the well-known “Lockhart equation” which is a one-dimensional (1D) viscoplasticity-like model with a yield surface from the point of view of solid mechanics.

In principle, a three-dimensional (3D) growth model may be proposed using an elasto-viscoplasticity-like theoretical framework taking the Lockhart equation as its consistent volumetric/uniaxial format. Boudaoud [15] investigated the growth of isolated walled cells, where cell growth is similar to a perfectly plastic deformation. Dumais et al. [16] developed an anisotropic viscoplastic thin shell model of cell walls which was the first attempt at integrating mechanical deformation driven by the turgor pressure with new material deposition and transport. Huang et al. [17,18] proposed the finite strain fibre-reinforced models of cell wall growth which took into account the non-linear elasticity, the mechanical inelastic behaviour (e.g., viscoelastic behaviour which performs differently in loading and unloading), and the temperature-dependent growth driven by stress and sustained with new material deposition and transport. In consistency with the aforementioned biological postulates of cell wall growth [2–4,9,10], those phenomenological models [16,17] proposed a key assumption that cell wall expansive growth may be decomposed into (i) a directly stress-driven irreversible isochoric in-plane expansion and (ii) the associated out-of-plane volumetric growth achieved by material deposition, transport, and integration (see Figure 2(a)).

Figure 2.

Cell wall growth compared with additive manufacturing (AM): (a) Cell wall material is schematically represented by a laminated composite material with one family of CMFs in each layer pointing to different directions. Cell wall growth includes the in-plane isochoric expansion (surface enlargement) and the matching volumetric growth maintaining a stable wall thickness. The volumetric growth results from the material deposition on the wall surface and the subsequent material transport in the same direction from the wall surface through the wall thickness. CMFs may (or may not) rotate depending on the stress distribution. (b) AM shows out-of-plane volumetric increase in a direction opposite to the direction of material deposition.

During the last decade, finite element method has been increasingly applied to cell wall growth modelling (see, e.g., Yanagisawa et al. [19] and the reviews [4,5]), which provided a computational tool to support the new experimental development using atomic force microscopy and X-ray scattering [20,21]. A number of publications [4,5,6,22–24] provided the extensive reviews about the recent progress on the multiscale modelling of cell wall and its growth from different points of view. The more general theories about the mathematical modelling of biological growth and remodelling could be found in the previous studies [25–28].

One of the most noticeable progresses in recent years is the continuously deepening understanding of the cell wall structure, the integration mechanisms of the wall components, and the corresponding mechanical regulations of growth at nano- and larger scale [4,5]. The knowledge of the geometry, structure, properties, and roles of the individual component of the cell wall (i.e., CMF, hemicellulose, pectin) have been improved so extensively as not only to challenge some dominated concepts in textbook [5,29] but also to stimulate new developments of modelling at different scales. From the point of view of continuum modelling, those new developments at nano- or larger scale may serve a similar role as that of crystal plasticity providing a microscopic interpretation to the phenomenological plasticity theory. Thus, this study aims to connect those recent progresses to the macroscopic continuum theory. For this purpose, we pay special attention to the following observation.

Among the new postulates proposed recently, the one suggested that “… selective control of side-by-side slippage and lateral separation of CMFs may provide a flexible means to modulate the directionality of cell growth …” [5] actually is consistent with the phenomenological model using Spencer’s deviatoric stress tensor [30] of fibre-reinforced composites as the driving force of the wall expansive growth [17]. However, in more detail, the other new postulate suggested that, at a closer look, instead of a tethered network in which the mechanical connections between CMFs were made by extended xyloglucan chains (part of hemicellulose) that adhered non-covalently to CMF surfaces, CMFs make a load-bearing network via close physical contacts with one another in bundled regions that potentially function as sites of cell wall loosening and creep [5,31]. Moreover, the third key component of cell wall, i.e., pectin, is supposed mainly to play a role of space-filling, not being tethered to CMFs. The detailed study outlined that the various wall properties are not tightly coupled, and therefore reflect distinctive aspects of wall structure [32].

This new concept about the load-bearing network is actually relevant to a prevailing assumption in the fibre-reinforced material modelling, i.e., the perfectly bonding constraint between reinforcement fibres and a matrix [33,34] so that “the fibres convect with the (matrix) material” [30]. This perfectly bonding assumption, although being convenient and successful to certain extent till now, restrains the modelling from capturing more complicated micromechanical behaviours of cell wall growth where wall properties are not tightly coupled [32], and breakage and reconnection of the crosslinks at nano-scale must be taken into account [9,10,35,36].

Even at the macroscopic scale, the perfectly bonding assumption also causes some conceptional difficulties. Different from additive manufacturing (AM) where the deposition of new materials constantly increases the total number of layers and total size of structure in an out-of-plane direction opposite to the deposition direction (Figure 2(b)), the deposition of newly synthesised material on the innermost surface of the cell wall and the following transport of the deposited material through the wall matrix for the assembly of components, which is in the same direction of the deposition, does not substantially change the thickness of the cell wall throughout the growth process [16] (see Figure 2(a)). This indicates that the new and existing CMFs need to be integrated together to maintain stable geometry and strength of the cell wall, which leads to postulates of remodelling of CMFs beyond the perfectly bonding assumption.

Remodelling may refer to changes of material microstructures which can be characterised by constituents, density, porosity, tortuosity, higher order correlation functions, and directionality [37]. In this study, the remodelling is focused on and confined to re-orientation of CMFs, which, in this study, is defined as the extra rotation to deviate from a direction complying with the perfectly bonding assumption. This naturally relaxes the perfectly bonding assumption. Moreover, antisymmetric part of the growth deformation gradient, which has been missing from most, if not all, of the existing growth models, is considered here to enrich the description of growth. We propose a theory framework of grow spin and demonstrate that both the antisymmetric part of the growth deformation gradient of the cell wall and the remodelling of CMFs can be covered in this framework and denoted as growth spin of the matrix and remodelling spin of CMFs, respectively. Special attention is paid to hardening laws of both growth and remodelling as one of the key parts to characterise living cell walls.

Noting that Dafalias has proposed a non-affine rotation theory based on the concept of plastic spin theory [38] which aimed at relaxing perfectly bonding assumption in a variety of engineering materials. The orientation distribution function (ODF) was used to represent the collective rotations of a constituent in the microstructure of a continuum. As the motivation of this study is to understand the remodelling of individual family of CMFs and to connect to nano-scale findings, the growth spin of the matrix and the remodelling spin of each and every family of CMFs are assumed independent and having their own constitutive laws. This may give the fidelity for understanding the details in the microstructure and provide the micromechanical information for the more general and more computationally efficient method proposed by Dafalias.

Through such modelling of individual family of CMFs, it is noted that there may be some similarities between cell wall growth and crystal plasticity even though seemingly being totally different subjects: (1) anisotropy; (2) multiple slip/sliding systems; (3) evidence of the existence of an anisotropic yield surface as a continuum; (4) deviatoric stresses as the driving forces; (5) isochoric deformation as the leading/dominant part of inelastic deformation; (6) inelastic volumetric deformation playing a supporting/ignorable role; (7) existence of hardening and softening.

Therefore, the growth spin model proposed in the present work may be capable to take into account the new experimental observation and new postulates. The model may also serve to bridge the gap between continuum theory and atom-/nano-scale computational methods [20,39,40]. As such, the proposed model may shed a new light into micro-structural interpretation of the “fictitious”local intermediate (growth) configuration in the kinematical growth [41] in a way similar to the crystal plasticity interpretation (e.g., slip/twining systems) of the intermediate configuration in the macroscopic plasticity theory [42–44]. This study shows broad connections between cell wall growth and the conventional viscoplasticity. It is also interesting to understand the connection between cell wall and the Cosserat models of fibre-reinforced materials [45–47] which is out of the scope of this study but may be worthy of further exploration.

This study is outlined as follows. First, a kinematical model including the growth spin and the remodelling spin are proposed. Second, a framework of constitutive model is discussed by adopting and extending Dafalias’ viscoplasticity model reported in Dafalias [48]. Third, a discussion about (i) the invariants and (ii) the spin representation in terms of Spencer’s deviatoric stress tensor [30] is presented to provide a technical foundation of modelling. A case study of cell wall as a cylindrical wall is presented to demonstrate the proposed theory, where Spencer’s deviatoric stress tensor plays a key role for the modelling of both growth and remodelling.

2. Kinematics of the model

2.1. Definition of the different configurations

We consider a spatial representative volume (REV) composed of a representative domain of a deformable solid continuum in the current configuration, ${\mathcal{B}}_t$ (see Figure 3). At the current time, $t\in [0,+\infty )$ , the REV of a spatial point $\mathit{x}\in {\mathcal{B}}_t\subset {\mathbb{R}}^3$ is a geometric domain, $\mathcal{U}(\mathit{x},t)\subset {\mathcal{B}}_t\subset {\mathbb{R}}^3$ , with a closed surface, $\mathcal{S}(\mathit{x},t)$ , embedded into the Euclidean space ${\mathbb{R}}^3$ . Due to growth, this domain represents the material domain of an open system. For an infinitesimal increment of time, a solid phase of this domain may have incremental growth and microstructure remodelling (G & R) imposed upon it. So, it is reasonable to establish an incremental growth theory on $\mathcal{U}(\mathit{x},t)$ . It is assumed that through such an incremental theory, we may trace back to a reference configuration, ${\mathcal{B}}_0$ , at a certain initial time, $t_0$ , with a material point labelled by $\mathit{X}\in {\mathcal{U}}_0(\mathit{X},t_0)\subset {\mathcal{B}}_0$ so that there is a relation, $\mathit{X}={\mathrm{\chi }}^{-1}(\mathit{x},t;t_0)$ , whereby ${\mathcal{U}}_0(\mathit{X},t_0)$ is the reference geometry domain of $\mathcal{U}(\mathit{x},t)$ and ${\mathrm{\chi }}^{-1}:{\mathcal{B}}_t\to {\mathcal{B}}_0$ is the mapping to trace back the material point $\mathit{X}$ .

Figure 3.

Kinematics of modelling: the REV of a spatial point $\mathit{x}$ as a geometric domain $\mathcal{U}$ with a closed surface $\mathcal{S}$ .

In general, the surface $\mathcal{S}$ of the REV is not a material surface because of the existence of mass fluxes passing through the surface $\mathcal{S}$ to supply the new material for growth [49]. In this study, it is assumed that there is no mass flux passing through the surface $\mathcal{S}$ so that $\mathcal{S}$ is a material surface. This is consistent with the definition of the first-order kinematical growth [49–52] and the interpretation of the G & R model proposed by Humphrey and Rajagopal [53] that an actual growth process including mass transport may be mathematically modelled as a local mass generation and absorption process. Under such an assumption, the new material deposition and transport in cell wall can be modelled as a local volumetric growth inside the domain $\mathcal{U}(\mathit{x},t)$ matching the local irreversible isochoric expansion for maintaining a certain geometry or structural stability condition, which leads to a non-associated flow rule with the contributions from both expansion and volumetric growth as shown later in section 5.5. Moreover, non-symmetry of the Cauchy stress tensor in an open system [49,50] may be avoided tactically. Hence, growth and remodelling of the cell wall may be expressed in a framework analogous to that of a closed system.

Consequently, the kinematics may be described in a way similar to the conventional continuum theory of a closed system. The mapping $\mathrm{\chi }:{\mathcal{B}}_0\to {\mathcal{B}}_t$ is defined as a deformation which maps the reference configuration, ${\mathcal{B}}_0$ , onto the current configuration, ${\mathcal{B}}_t$ , so that the current position of the material point $\mathit{X}$ is written as $\mathit{x}=\mathrm{\chi }(\mathit{X},t)\in {\mathcal{B}}_t$ , where $t_0$ is dropped from the expression for conciseness. The deformation gradient is defined as $\mathit{F}(\mathit{X},t)=(\partial \mathrm{\chi }/\partial \mathit{X})$ with its Jacobian denoted as $J=det\mathit{F}>0$ .

Let a unit vector, ${\mathit{a}}_{i0}(\mathit{X})$ , denotes the initial fibre direction of the $i\mathrm{th}$ family $(i\in \{1,\dots ,n_f\})$ of fibres attached to a solid particle at the material point $\mathit{X}$ on ${\mathcal{B}}_0$ . In theory, $n_f$ , the total number of the families of fibres with different directions in the REV, could be any finite number depending on the feature of material, e.g., symmetry. In this study, the fibre may be referred to either the physical CMF or a pure geometrical direction, e.g., the normal direction of the cell wall surface. For clarity, we shall explicit use “CMF” where the fibre refers to physical CMF. Noting that, from the point of view of laminated composite, this study models all CMFs in the REV as a discrete set of finite number families of CMFs with distinguishing directions, which is capable to capture the major features of the cell wall of interest. In reality, either from perspective of single lamina or multiple plies of cell wall, CMFs may have a continuous dispersion of fibre angles, which can be described using a probability density function of a continuous random orientational variable. The reader may refer to the plastic spin theory using the ODF proposed by Dafalias [38], the soft matter model with distributed collagen fibre orientations presented by Gasser et al. [33], and the fibres orientation density distribution for fibrous tissues growth and remodelling by Lanir [27].

The current spatial fibre direction of ${\mathit{a}}_{i0}(\mathit{X})$ at the same material point is as follows:

$${\tilde{\mathit{a}}}_i=\frac{{\mathit{a}}_i}{\sqrt{{\lambda }_{{\mathit{a}}_i}}}\forall i\in \{1,\dots ,n_f\},$$

(1)

where the deformed fibre and its stretch are described by:

$${\mathit{a}}_i=\mathit{F}(\mathit{X},t){\mathbf{\Omega }}_i^{-1}(\mathit{X},t){\mathit{a}}_{i0}(\mathit{X})\mathrm{and}{\lambda }_{{\mathit{a}}_i}=\mathit{C}:({\mathbf{\Omega }}_i^{-1}{\mathit{A}}_{i0}{\mathbf{\Omega }}_i)=\mathrm{tr}{\mathit{A}}_i,$$

(2)

in which “:” denotes the double contraction, “ $\mathrm{tr}$ ” is the trace operator on ${\mathcal{B}}_t$ , the right Cauchy–Green tensor is defined as:

$$\mathit{C}={\mathit{F}}^{\mathrm{T}}\mathit{F}.$$

(3)

${\mathit{A}}_{i0}={\mathit{a}}_{i0}\otimes {\mathit{a}}_{i0}$ and ${\mathit{A}}_i={\mathit{a}}_i\otimes {\mathit{a}}_i$ are the two structural tensors, and ${\mathbf{\Omega }}_i^{-1}(\mathit{X},t)$ represents a rotation of the fibre due to remodelling. The tensor ${\mathbf{\Omega }}_i$ itself is an orthogonal tensor as the representation of rotation. The superscript $\mathrm{T}$ stands for the operator of matrix transpose while ⊗ is the tensor product. To differentiate from the fixed fibre direction ${\mathit{a}}_{i0}(\mathit{X})$ , a time-dependent fibre direction on ${\mathcal{B}}_0$ is introduced as:

$${\mathit{a}}_i^r(\mathit{X},t)={\mathbf{\Omega }}_i^{-1}(\mathit{X},t){\mathit{a}}_{i0}(\mathit{X})\forall i\in \{1,\dots ,n_f\},$$

(4)

with the corresponding structural tensor, ${\mathit{A}}_i^r={\mathit{a}}_i^r\otimes {\mathit{a}}_i^r$ , to represent the evolution of fibre direction with remodelling.

Strictly speaking, we may consider equation (4) as a definition of the remodelling configuration of the $i\mathrm{th}$ family of CMFs, denoted as ${\mathcal{B}}_i^r$ . ${\mathcal{B}}_i^r$ is identical to ${\mathcal{B}}_0$ except the fibre direction so that it has the same mapping $\mathrm{\chi }:{\mathcal{B}}_i^r\to {\mathcal{B}}_t$ and the same deformation gradient $\mathit{F}$ from ${\mathcal{B}}_i^r$ to ${\mathcal{B}}_t$ . In other words, if the unit tensor $\mathit{I}$ (i.e., $[\mathit{I}{]}_{\mathrm{ij}}={\delta }_{\mathrm{ij}}$ , where ${\delta }_{\mathrm{ij}}$ is the Kronecker delta) is the metric tensor on ${\mathcal{B}}_t$ , the metric tensor on ${\mathcal{B}}_i^r$ is $\mathit{C}$ , which is helpful for clarification of some concepts of the remodelling spin as we may see in the later discussion.

Noting that, to avoid confusion in formulation, the summation rule is not applied to the indexes indicating the families of fibres (e.g., the index “ $i$ ” in equations (1) and (4)) throughout this study.

As the base of kinematical growth modelling, the assumption of the kinematical decomposition of the total deformation gradient,

$$\mathit{F}={\mathit{F}}^e{\mathit{F}}^g,$$

(5)

which was proposed initially for plasticity [54] and then found its place in biomechanics [55], is used in this study, whereby the superscripts “ $e$ ” and “ $g$ ” refer to elastic and growth parts, respectively. Bear in mind that ${\mathit{F}}^g$ includes the contributions from both the stress-driven isochoric expansion and the matching volumetric growth resulting from material deposition and transport. The decomposition (5) naturally introduces an intermediate configuration, $\overline{\mathcal{B}}$ , i.e., the growth configuration generated by ${\mathit{F}}^g$ , which is a collection of local domains since in general ${\mathit{F}}^g$ is not a compatible tensorial field. For this reason, $\overline{\mathcal{B}}$ as a whole may be considered as a fictitious configuration.

The counterparts of the vector ${\mathit{a}}_i$ and tensor ${\mathit{A}}_i$ on $\overline{\mathcal{B}}$ , which may be obtained by a pull-back from ${\mathcal{B}}_t$ to $\overline{\mathcal{B}}$ , are denoted as:

$${\overline{\mathit{a}}}_i(\mathit{X},t)={\mathit{F}}^{e^{-1}}(\mathit{X},t){\mathit{a}}_i\equiv {\mathit{F}}^g(\mathit{X},t){\mathbf{\Omega }}_i^{-1}(\mathit{X},t){\mathit{a}}_{i0}(\mathit{X})\equiv {\mathit{F}}^g(\mathit{X},t){\mathit{a}}_i^r(\mathit{X},t)$$

(6)

and ${\overline{\mathit{A}}}_i={\overline{\mathit{a}}}_i\otimes {\overline{\mathit{a}}}_i$ , respectively. Noting that the hypothesis of no growth in a CMF direction due to the CMF’s inertia in growth is adopted in this study. As shown in later discussion, Spencer’s deviatoric stress tensor serves as a proper driving force of growth which imposes a certain constraint on ${\mathit{F}}^g$ so that $\mathrm{tr}{\overline{\mathit{A}}}_i\equiv 1$ and the vector ${\overline{\mathit{a}}}_i$ is automatically a unit vector on $\overline{\mathcal{B}}$ . Hence, ${\mathit{a}}_{i0}$ and ${\mathit{a}}_i^r$ on ${\mathcal{B}}_0$ , ${\overline{\mathit{a}}}_i$ on $\overline{\mathcal{B}}$ , and ${\tilde{\mathit{a}}}_i$ on ${\mathcal{B}}_t$ are all treated as the unit vectors.

Equation (6) indicates that the tensor ${\mathbf{\Omega }}_i^{-1}(\equiv {\mathbf{\Omega }}_i^{\mathrm{T}})$ may be interpreted as the deviation of the fibre direction, ${\overline{\mathit{a}}}_i$ , from a virtual direction, ${\mathit{F}}^g{\mathit{a}}_{i0}$ , complying with the perfectly bonding assumption. By contrast, the elastic deformation of a fibre still obeys the perfectly bonding condition since ${\mathit{a}}_i={\mathit{F}}^e{\overline{\mathit{a}}}_i$ .

For convenience, the kinematical relations of the different configurations, particles, and fibre directions during a growth process are shown in Figure 3.

2.2. Kinematical relations in rate formats and the growth spin

Aiming at an incremental description of growth, kinematical relations in rate formats are discussed here, which naturally introduces the concept of the growth spin of the wall matrix.

Let $\mathit{L}=(\partial \stackrel{·}{\mathrm{\chi }}/\partial \mathit{x})=\stackrel{·}{\mathit{F}}{\mathit{F}}^{-1}$ be the spatial velocity gradient on ${\mathcal{B}}_t$ , where the symbol of superimposed dot denotes the material time derivative. An additive decomposition of the velocity gradient, which is the rate format of equation (5), may be expressed as:

$$\mathit{L}={\mathit{L}}^e+{\mathit{L}}^g,$$

(7)

where elastic and growth parts of $\mathit{L}$ are defined as:

$${\mathit{L}}^e=\stackrel{·}{{\mathit{F}}^e}{\mathit{F}}^{e^{-1}}\mathrm{and}{\mathit{L}}^g={\mathit{F}}^e{\overline{\mathit{L}}}^g{\mathit{F}}^{e^{-1}},$$

(8)

respectively, in which the counterpart of ${\mathit{L}}^g$ on $\overline{\mathcal{B}}$ is:

$${\overline{\mathit{L}}}^g={\stackrel{·}{\mathit{F}}}^g{\mathit{F}}^{g^{-1}}.$$

(9)

In consistence with equation (9), a pull-back of ${\mathit{L}}^e$ from ${\mathcal{B}}_t$ to $\overline{\mathcal{B}}$ yields:

$${\overline{\mathit{L}}}^e={\mathit{F}}^{e^{-1}}{\mathit{L}}^e{\mathit{F}}^e={\mathit{F}}^{e^{-1}}{\stackrel{·}{\mathit{F}}}^e.$$

(10)

Furthermore, $\mathit{L}$ , ${\mathit{L}}^e$ , and ${\mathit{L}}^g$ may be decomposed into symmetric and skew-symmetric parts as:

$$\mathit{L}=\mathit{D}+\mathit{W}$$

(11)

and

$${\mathit{L}}^e={\mathit{D}}^e+{\mathit{W}}^e\mathrm{and}{\mathit{L}}^g={\mathit{D}}^g+{\mathit{W}}^g,$$

(12)

whereby the rate of deformation tensor, $\mathit{D}$ , and the material spin tensor, $\mathit{W}$ , of the whole continuum of the REV are defined as:

$$\mathit{D}=\mathrm{sym}\mathit{L}{=}^{\mathrm{def}}\frac{1}{2}(\mathit{L}+{\mathit{L}}^{\mathrm{T}})\mathrm{and}\mathit{W}=\mathrm{skew}\mathit{L}{=}^{\mathrm{def}}\frac{1}{2}(\mathit{L}-{\mathit{L}}^{\mathrm{T}}),$$

(13)

respectively. In a similar way, the spatial elastic and growth rates of deformation may be defined as:

$${\mathit{D}}^e=\mathrm{sym}{\mathit{L}}^e=\frac{1}{2}\left({\mathit{L}}^e+{\mathit{L}}^{e^{\mathrm{T}}}\right)\mathrm{and}{\mathit{D}}^g=\mathrm{sym}{\mathit{L}}^g=\frac{1}{2}\left({\mathit{L}}^g+{\mathit{L}}^{g^{\mathrm{T}}}\right),$$

(14)

respectively. Moreover, the growth spin, ${\mathit{W}}^g$ , and rigid body spin, $\mathbf{\omega }$ , of the matrix, i.e., crosslink network including the hemicellulose and pectins in the REV, are defined as:

$${\mathit{W}}^g=\mathrm{skew}{\mathit{L}}^g=\frac{1}{2}\left({\mathit{L}}^g-{\mathit{L}}^{g^{\mathrm{T}}}\right)\mathrm{and}\mathbf{\omega }=\mathit{W}-{\mathit{W}}^g,$$

(15)

respectively (see Figure 4).

Figure 4.

Growth spin ${\mathit{W}}^g$ and remodelling spin ${\mathit{W}}_i^g$ : (i) The reference configuration $(a)$ is mapped in sequence to the configurations $(b)$ , $(c)$ , $(d)$ , and $(e)$ by ${\mathit{D}}^g$ , ${\mathit{L}}^g$ , ${\mathit{L}}^g+\mathbf{\omega }$ , and $\mathit{L}$ , respectively and (ii) The reference fibre direction $(1)$ is mapped to the direction $(2)$ by ${\mathit{L}}^g$ complying with the perfectly bonding assumption and then re-orientates to the direction $(3)$ by the remodelling ${\mathit{W}}_i^g$ .

Correspondingly, pulling ${\mathit{D}}^e$ and ${\mathit{D}}^g$ from ${\mathcal{B}}_t$ back to $\overline{\mathcal{B}}$ yields:

$${\overline{\mathit{D}}}^e=\mathrm{sym}{\overline{\mathit{L}}}^e=\frac{1}{2}\left(\overline{\mathit{C}}{\overline{\mathit{L}}}^e+{\overline{\mathit{L}}}^{e^{\mathrm{T}}}\overline{\mathit{C}}\right)\mathrm{and}{\overline{\mathit{D}}}^g=\mathrm{sym}{\overline{\mathit{L}}}^g=\frac{1}{2}\left(\overline{\mathit{C}}{\overline{\mathit{L}}}^g+{\overline{\mathit{L}}}^{g^{\mathrm{T}}}\overline{\mathit{C}}\right),$$

(16)

where we have the relations, ${\overline{\mathit{D}}}^e={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{D}}^e{\mathit{F}}^e$ and $\overline{\mathit{C}}={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{F}}^e$ . And the pull-back of ${\mathit{W}}^g$ from ${\mathcal{B}}_t$ to $\overline{\mathcal{B}}$ reads:

$${\overline{\mathit{W}}}^g=\mathrm{skew}{\overline{\mathit{L}}}^g=\frac{1}{2}\left(\overline{\mathit{C}}{\overline{\mathit{L}}}^g-{\overline{\mathit{L}}}^{g^{\mathrm{T}}}\overline{\mathit{C}}\right).$$

(17)

Noting that ${\overline{\mathit{D}}}^g+\overline{{\mathit{W}}^g}(=\overline{\mathit{C}}{\overline{\mathit{L}}}^g)\ne {\overline{\mathit{L}}}^g.$

2.3. Remodelling spin on ${\mathcal{B}}_t$

Bear in mind that ${\mathit{W}}^g$ is defined as the growth spin of the wall matrix, with which CMFs convect. We shall separately introduce the spin of CMFs due to the presence of remodelling. To differentiate from ${\mathit{W}}^g$ , such a spin of CMFs is called the remodelling spin.

For the $i\mathrm{th}$ family of CMFs $(i\in \{1,\dots ,n_f\})$ , an additive decomposition on ${\mathcal{B}}_t$ is introduced as (see Figure 4):

$$\underset{\mathrm{growth}\mathrm{spin}\mathrm{of}\mathrm{matrix}}{\underbrace{{\mathit{W}}^g(=\mathrm{skew}{\mathit{L}}^g)}}=\underset{\mathrm{constitutive}\mathrm{spin}\mathrm{of}i\mathrm{th}\mathrm{f}.\mathrm{CMFs}}{\underbrace{{\mathbf{\omega }}_i}}+\underset{\mathrm{remodelling}\mathrm{spin}\mathrm{of}i\mathrm{th}\mathrm{f}.\mathrm{CMFs}}{\underbrace{{\mathit{W}}_i^g}},$$

(18)

which is in a similar form as equation (12) in Dafalias [56] so that it clearly indicates the connection/analogy between cell wall and crystal plasticity, whereby the remodelling spin, ${\mathit{W}}_i^g$ , of the $i\mathrm{th}$ family of CMFs on ${\mathcal{B}}_t$ is defined by a transformation,

$${\mathit{W}}_i^g={\mathit{F}}^{-\mathrm{T}}\mathrm{skew}{\mathit{W}}_{0i}^g{\mathit{F}}^{-1},$$

(19)

in which:

$${\mathit{W}}_{0i}^g=-{\stackrel{·}{\mathbf{\Omega }}}_i^{\mathrm{T}}{\mathbf{\Omega }}_i\equiv {\mathbf{\Omega }}_i^{\mathrm{T}}{\stackrel{·}{\mathbf{\Omega }}}_i,$$

(20)

is the corresponding remodelling spin on ${\mathcal{B}}_0$ and its skew-symmetric part is defined as:

$$\mathrm{skew}{\mathit{W}}_{0i}^g=\frac{1}{2}({\mathit{CW}}_{0i}^g+{\mathit{W}}_{0i}^g\mathit{C}).$$

(21)

In order to interpret the form of equation (21), we may recall the remodelling configuration ${\mathcal{B}}_i^r$ defined in section 2.1. The remodelling spin ${\mathit{W}}_{0i}^g$ , as shown in equation (20), actually is a “velocity gradient” on ${\mathcal{B}}_i^r$ corresponding to the “deformation gradient,” ${\mathbf{\Omega }}_i^{-1}$ . As aforementioned, the metric tensor on ${\mathcal{B}}_i^r$ is $\mathit{C}$ . Hence, the $\mathrm{skew}$ operator on ${\mathcal{B}}_i^r$ is in a form of $\mathrm{skew}(·)=\frac{1}{2}(\mathit{C}(·)-(·{)}^{\mathrm{T}}\mathit{C})$ , which explains equation (21).

Given ${\mathit{W}}^g$ and ${\mathit{W}}_i^g$ , the corresponding constitutive spin, ${\mathbf{\omega }}_i$ , of the $i\mathrm{th}$ family of CMFs is defined by equation (18) as the difference between ${\mathit{W}}^g$ of the matrix and ${\mathit{W}}_i^g$ of the $i\mathrm{th}$ family of CMFs. ${\mathit{W}}_i^g$ may be interpreted as the remodelling of the $i\mathrm{th}$ family of CMFs with respect to ${\mathit{W}}^g$ of the matrix. Thus ${\mathbf{\omega }}_i$ is the rate of rotation of the $i\mathrm{th}$ family of CMFs taking into account the convected rotation with the growth spin of the matrix and the remodelling spin of CMFs.

In the specific case that there is no remodelling, i.e., ${\mathit{W}}_i^g=0,$ equation (18) indicates that ${\mathit{W}}^g={\mathbf{\omega }}_i,$ which retains the perfectly bonding assumption for the growth of the whole continuum. Based on the above discussion, $\mathbf{\omega }$ of the matrix is a rigid body spin which should be eliminated from the constitutive model of the continuum. Therefore, the objective rate of a tensor, say $\mathit{A}$ , associated with the continuum, with respect to $\mathbf{\omega }$ , may be defined as:

$$\stackrel{◦}{\mathit{A}}=\stackrel{·}{\mathit{A}}-\mathbf{\omega }\mathit{A}+\mathit{A}\mathbf{\omega },$$

(22)

where the superposed ° follows the notation proposed by Dafalias [48].

In order to define the objective rate of tensors associated with the $i\mathrm{th}$ family of CMFs $(i\in \{1,\dots ,n_f\})$ , we use the general definition of Lie derivative [57–60]:

$$L_{\hat{\mathbf{\upsilon }}}(·)={\hat{\mathrm{\chi }}}_{*}(\frac{\mathrm{D}}{\mathrm{D}t}{\hat{\mathrm{\chi }}}^{*}(·)),$$

(23)

in which $\hat{\mathrm{\chi }}$ stands for the motion composed of the mapping $\mathrm{\chi }$ of the continuum and the remodelling rotation ${\mathbf{\Omega }}_i^{-1}$ while the subscript “ $\hat{\mathbf{\upsilon }}$ ” stands for the velocity field of $\hat{\mathrm{\chi }}$ , ${\hat{\mathrm{\chi }}}_{*}$ and ${\hat{\mathrm{\chi }}}^{*}$ are the corresponding pull-back and push-forward operators, respectively.

For the specification of $L_{\hat{\mathbf{\upsilon }}}(·)$ , let us first consider the motion $\hat{\mathrm{\chi }}$ of the unconstrained structural tensor, ${\mathit{A}}_i={\mathit{a}}_i\otimes {\mathit{a}}_i$ , defined in equation (2)₁. So that we have the specified pull-back and push-forward as follows:

$${\hat{\mathrm{\chi }}}_{*}(·)=(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{-\mathrm{T}}(·)(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{-1}\mathrm{and}{\hat{\mathrm{\chi }}}^{*}(·)=(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{\mathrm{T}}(·)(\mathit{F}{\mathbf{\Omega }}_i^{-1}).$$

(24)

Substituting equation (24) into equation (23) yields the objective rate of a tensor associated with the $i\mathrm{th}$ family of CMFs, say ${\mathit{B}}_i$ which is defined on ${\mathcal{B}}_t$ and could be either a geometric or a stress tensor, as:

$${\stackrel{{\nabla }^c}{\mathit{B}}}_i{=}^{\mathrm{def}}{L}_{\hat{\mathbf{\upsilon }}}{\mathit{B}}_i={\stackrel{·}{\mathit{B}}}_i-(\mathit{L}-{\hat{\mathit{W}}}_i^g){\mathit{B}}_i-{\mathit{B}}_i(\mathit{L}-{\hat{\mathit{W}}}_i^g{)}^{\mathrm{T}},$$

(25)

where the superscript “ $c$ ” stands for the current configuration, and ${\hat{\mathit{W}}}_i^g$ is defined and may be decomposed as:

$${\hat{\mathit{W}}}_i^g{=}^{\mathrm{def}}{\mathit{FW}}_{0i}^g{\mathit{F}}^{-1}\equiv {\mathit{W}}_i^g+{\tilde{\mathit{W}}}_i^g,$$

(26)

in which we recognise the remodelling spin ${\mathit{W}}_i^g$ and define a symmetric tensor,

$${\tilde{\mathit{W}}}_i^g={\mathit{F}}^{-\mathrm{T}}\mathrm{sym}{\mathit{W}}_{0i}^g{\mathit{F}}^{-1}.$$

(27)

Whereby, in consistency with equation (21), the symmetry part of ${\mathit{W}}_{0i}^g$ is defined as:

$$\mathrm{sym}{\mathit{W}}_{0i}^g=\frac{1}{2}({\mathit{CW}}_{0i}^g-{\mathit{W}}_{0i}^g\mathit{C}).$$

(28)

It is understood that in general $\mathrm{sym}{\mathit{W}}_{0i}^g\ne 0$ if ${\mathit{W}}_{0i}^g\ne 0$ . Equation (28) indicates that, given a tensor $\mathit{C}$ , $\mathrm{sym}{\mathit{W}}_{0i}^g$ may be computed either directly from ${\mathit{W}}_{0i}^g$ or from ${\mathit{W}}_i^g$ via solving ${\mathit{W}}_{0i}^g$ with equations (19) and (21).

The relation ${\mathit{W}}_i^g={\mathit{W}}^g-{\mathbf{\omega }}_i$ in equation (18) and $\mathbf{\omega }=\mathit{W}-{\mathit{W}}^g$ in equation (15)₂ indicate that:

$${\mathit{W}}_i^g=\mathit{W}-\mathbf{\omega }-{\mathbf{\omega }}_i.$$

(29)

Furthermore, using the decomposition in equation (11), $\mathit{L}=\mathit{D}+\mathit{W}$ , equation (29) becomes:

$${\mathit{W}}_i^g=\mathit{L}-\mathit{D}-\mathbf{\omega }-{\mathbf{\omega }}_i.$$

(30)

Substituting equation (30) into equation (26) and then conducting some algebra rearrangment, we obtain:

$$\mathit{L}-{\hat{\mathit{W}}}_i^g=\mathbf{\omega }+\mathit{D}+{\mathbf{\omega }}_i-{\tilde{\mathit{W}}}_i^g.$$

(31)

With this relationship, the comparison between equations (25) and (22) indicates that:

$${\stackrel{{\nabla }^c}{\mathit{B}}}_i={\stackrel{◦}{\mathit{B}}}_i-(\mathit{D}+{\mathbf{\omega }}_i-{\tilde{\mathit{W}}}_i^g){\mathit{B}}_i-{\mathit{B}}_i(\mathit{D}-{\mathbf{\omega }}_i-{\tilde{\mathit{W}}}_i^g).$$

(32)

In the case of no remodelling, we have ${\mathbf{\omega }}_i={\mathit{W}}^g$ and ${\tilde{\mathit{W}}}_i^g=0$ . Consequently, equation (32) degenerates to an expression of the continuum.

For further understanding the remodelling spin ${\mathit{W}}_i^g$ and the objective rate $\stackrel{{\nabla }^c}{(·)}$ , let us consider a special case that the structural tensor ${\mathit{A}}_i$ is assigned to ${\mathit{B}}_i$ , i.e., ${\mathit{B}}_i:={\mathit{A}}_i$ . The push-forward ${\hat{\mathrm{\chi }}}^{*}$ in equation (24)₂ and the transformation in equation (2)₁ indicate that ${\hat{\mathrm{\chi }}}^{*}{\mathit{A}}_i={\mathit{A}}_{i0}$ . Hence, $\mathrm{D}/\mathrm{D}t({\hat{\mathrm{\chi }}}^{*}{\mathit{A}}_i)=0$ , which, using the definition (23), comes to the conclusion:

$${\stackrel{{\nabla }^c}{\mathit{A}}}_i=L_{\hat{\mathbf{\upsilon }}}{\mathit{A}}_i\equiv 0\forall i\in \{1,\dots ,n_f\}.$$

(33)

Equation (33) may be considered as a kinematical remodelling law on ${\mathcal{B}}_t$ based on the motion of ${\mathit{A}}_i$ .

Since $\mathrm{tr}{\mathit{A}}_i\ne 1$ under elastic deformations, the objective rate of the normalised structural tensor:

$${\tilde{\mathit{A}}}_i={\tilde{\mathit{a}}}_i\otimes {\tilde{\mathit{a}}}_i=(\mathrm{tr}{\mathit{A}}_i{)}^{-1}{\mathit{A}}_i,$$

(34)

may provide a more clear and direct description of remodelling.

In this case, the motion $\hat{\mathrm{\chi }}$ may be specified by a mapping derived from equations (1) and (2):

$${\tilde{\mathit{a}}}_i=(\frac{1}{\sqrt{{\lambda }_{{\mathit{a}}_i}}}\mathit{F}(\mathit{X},t){\mathbf{\Omega }}_i^{-1}(\mathit{X},t)){\mathit{a}}_{i0}(\mathit{X}).$$

(35)

Hence, the specification of the pull-back and push-forward is as follows:

$${\hat{\mathrm{\chi }}}_{*}(·)={\lambda }_{{\mathit{a}}_i}(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{-\mathrm{T}}(·)(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{-1}\mathrm{and}{\hat{\mathrm{\chi }}}^{*}(·)=\frac{1}{{\lambda }_{{\mathit{a}}_i}}(\mathit{F}{\mathbf{\Omega }}_i^{-1}{)}^{\mathrm{T}}(·)(\mathit{F}{\mathbf{\Omega }}_i^{-1}).$$

(36)

Substituting equation (36) into equation (23) yields the objective rate of a normalised tensor, ${\tilde{\mathit{B}}}_i={\lambda }_{{\mathit{a}}_i}^{-1}{\mathit{B}}_i$ , as:

$${\stackrel{{\nabla }^u}{\tilde{\mathit{B}}}}_i={\stackrel{.}{\stackrel{~}{\mathit{B}}}}_i-\left(\mathit{L}-\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right)\mathit{I}-{\widehat{\mathit{W}}}_i^g\right){\tilde{\mathit{B}}}_i-{\tilde{\mathit{B}}}_i{\left(\mathit{L}-\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right)\mathit{I}-{\widehat{\mathit{W}}}_i^g\right)}^{\mathrm{T}},$$

(37)

where the superscript “ $u$ ” stands for the unit vector at ${\mathcal{B}}_t$ , and we use the following result:

$${\stackrel{.}{\mathbf{\lambda }}}_{{\mathit{a}}_i}=2{\mathbf{\lambda }}_{{\mathit{a}}_i}\mathrm{tr}\left([\mathit{D}-{\tilde{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right)\equiv 2{\mathbf{\lambda }}_{{\mathit{a}}_i}\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right),$$

(38)

which is discussed in Appendix 1 (Calculation of ${\stackrel{·}{\lambda }}_{{\mathit{a}}_i}$ ).

The other way to derive equation (37) is to use the relation:

$${\stackrel{{\nabla }^u}{\tilde{\mathit{B}}}}_i=\stackrel{{\nabla }^c}{({\lambda }_{{\mathit{a}}_i}^{-1}{\mathit{B}}_i)}=-\stackrel{{\nabla }^c}{{\lambda }_{{\mathit{a}}_i}}{\lambda }_{{\mathit{a}}_i}^{-1}{\tilde{\mathit{B}}}_i+{\lambda }_{{\mathit{a}}_i}^{-1}\stackrel{{\nabla }^c}{{\mathit{B}}_i}.$$

(39)

Since ${\lambda }_{{\mathit{a}}_i}=\mathrm{tr}{\mathit{A}}_i$ , by calculating the trace of equation (25) with ${\mathit{B}}_i$ replaced by ${\mathit{A}}_i$ and then exchanging the sequence of trace operator and differentiation operator, we obtain:

$$\stackrel{{\nabla }^c}{{\lambda }_{{\mathit{a}}_i}}={\stackrel{{\nabla }^c}{\mathrm{tr}\mathit{A}}}_i=\mathrm{tr}{\stackrel{·}{\mathit{A}}}_i-2\mathrm{tr}([\mathit{L}-{\hat{\mathit{W}}}_i^g]{\mathit{A}}_i),$$

(40)

which is consistent with equation (38). Substituting equations (40) and (25) into the right-hand side (r.h.s) of equation (39) yields:

$${\stackrel{{\nabla }^u}{\tilde{\mathit{B}}}}_i=-{\lambda }_{{\mathit{a}}_i}^{-1}\left\{\stackrel{·}{\mathrm{tr}{\mathit{A}}_i}-2\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\mathit{A}}_i\right)\right\}{\tilde{\mathit{B}}}_i+{\lambda }_{{\mathit{a}}_i}^{-1}\left\{{\stackrel{.}{\mathit{B}}}_i-(\mathit{L}-{\widehat{\mathit{W}}}_i^g){\mathit{B}}_i-{\mathit{B}}_i{(\mathit{L}-{\widehat{\mathit{W}}}_i^g)}^{\mathrm{T}}\right\},$$

which, after some manipulations, becomes:

$${\stackrel{{\nabla }^u}{\tilde{\mathit{B}}}}_i={\stackrel{.}{\stackrel{~}{\mathit{B}}}}_i+2\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right){\tilde{\mathit{B}}}_i-(\mathit{L}-{\widehat{\mathit{W}}}_i^g){\tilde{\mathit{B}}}_i-{\tilde{\mathit{B}}}_i{(\mathit{L}-{\widehat{\mathit{W}}}_i^g)}^{\mathrm{T}}.$$

(41)

Equation (41) is identical to equation (37) so we prove this transformation-based method to obtain ${\stackrel{{\nabla }^u}{\tilde{\mathit{B}}}}_i$ .

It is of interest to observe that there is analogy between the relation, ${\tilde{\mathit{B}}}_i={\lambda }_{{\mathit{a}}_i}^{-1}{\mathit{B}}_i$ , and the relation, $\mathbf{\sigma }={\overline{J}}^{-1}\mathit{\tau }$ , of Cauchy stress and Kirchhoff stress, where $\overline{J}=\det {\mathit{F}}^e$ . That is why the objective rates (equations (25) and (41)) (and their relation) show analogies to the objective (e.g., Truesdell) rates of Cauchy stress and Kirchhoff stress (and their relation) [58,59]. This observation may be helpful for the clarification of the objective rates defined here.

Again, we consider a special case that the structural tensor ${\tilde{\mathit{A}}}_i$ is assigned to ${\tilde{\mathit{B}}}_i$ , i.e., ${\tilde{\mathit{B}}}_i:={\tilde{\mathit{A}}}_i$ . Using the same reasoning for deriving equation (33), we may obtain the other kinematical remodelling law on ${\mathcal{B}}_t$ :

$${\stackrel{{\nabla }^u}{\tilde{\mathit{A}}}}_i=0\forall i\in \{1,\dots ,n_f\}.$$

(42)

Equations (33) and (42) are consistent with the differential geometry concept of the invariance of a tensor field (e.g., ${\mathit{A}}_i$ or ${\tilde{\mathit{A}}}_i$ in here) under some transformation, which is expressed by the vanishing of Lie derivative of that tensor field (see Schutz [57], p. 86).

One of the application of equation (42) is to construct a structure-preserving increment of the unit vector ${\tilde{\mathit{a}}}_i$ :

$$\delta {\tilde{\mathit{a}}}_i=\left(\mathit{L}-\mathrm{tr}\left([\mathit{L}-{\widehat{\mathit{W}}}_i^g]{\tilde{\mathit{A}}}_i\right)\mathit{I}-{\widehat{\mathit{W}}}_i^g\right){\tilde{\mathit{a}}}_i.$$

(43)

It is straightforward to prove that:

$$\mathrm{tr}[\{{\tilde{\mathit{a}}}_i+\delta {\tilde{\mathit{a}}}_i\}\otimes \{{\tilde{\mathit{a}}}_i+\delta {\tilde{\mathit{a}}}_i\}]=0.$$

2.4. Remodelling spins on $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$

It is interesting to understand the counterparts of the kinematical remodelling law (42) on $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$ , respectively. For this purpose, we use the same pull-back as equation (10) to transform equation (18) from ${\mathcal{B}}_t$ to $\overline{\mathcal{B}}$ :

$$\overline{{\mathit{W}}^g}={\overline{\mathbf{\omega }}}_i+{\overline{\mathit{W}}}_i^g\forall i\in \{1,\dots ,n_f\},$$

(44)

where:

$${\overline{\mathbf{\omega }}}_i={\mathit{F}}^{e^{\mathrm{T}}}{\mathbf{\omega }}_i{\mathit{F}}^e\mathrm{and}{\overline{\mathit{W}}}_i^g={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{W}}_i^g{\mathit{F}}^e,$$

(45)

in which ${\overline{\mathit{W}}}_i^g$ is the remodelling spin on $\overline{\mathcal{B}}$ .

The objective rates on $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$ associated with the $i\mathrm{th}$ family of CMFs may be derived using the Lie derivative (23) in a process similar to what discussed in section 2.3. Specifically, for the objective rate on $\overline{\mathcal{B}}$ , the corresponding motion $\hat{\mathrm{\chi }}$ may be specified by a mapping $({\mathit{F}}^g(\mathit{X},t){\mathbf{\Omega }}_i^{-1}(\mathit{X},t)){\mathit{a}}_{i0}(\mathit{X})$ . Moreover, the objective rate on ${\mathcal{B}}_0$ may be obtained using a motion $\hat{\mathrm{\chi }}$ defined by the mapping in equation (4). However, rather than explicitly using the Lie derivative (23), here we show the other method, which is in a way similar to equation (39), to obtain the objective rates on $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$ through transformations from ${\mathcal{B}}_t$ to $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$ , respectively.

Let us consider a tensor associated with the $i\mathrm{th}$ family of CMFs, say ${\overline{\mathit{B}}}_i$ which is defined on $\overline{\mathcal{B}}$ and has a relation, ${\mathit{B}}_i={\mathit{F}}^e{\overline{\mathit{B}}}_i{\mathit{F}}^{e^{\mathrm{T}}}$ , with the aforementioned tensor ${\mathit{B}}_i$ on ${\mathcal{B}}_t$ . Using equations (8) and (10), this relation has the rate format as:

$${\stackrel{·}{\mathit{B}}}_i={\mathit{F}}^e({\stackrel{·}{\overline{\mathit{B}}}}_i+\overline{{\mathit{L}}^e}{\overline{\mathit{B}}}_i+{\overline{\mathit{B}}}_i{\overline{\mathit{L}}}^{e^{\mathrm{T}}}){\mathit{F}}^{e^{\mathrm{T}}}.$$

(46)

Substituting equation (46) into equation (25) yields:

$${\stackrel{{\nabla }^c}{\mathit{B}}}_i={\mathit{F}}^e({\stackrel{·}{\overline{\mathit{B}}}}_i+\overline{{\mathit{L}}^e}{\overline{\mathit{B}}}_i+{\overline{\mathit{B}}}_i{\overline{\mathit{L}}}^{e^{\mathrm{T}}}){\mathit{F}}^{e^{\mathrm{T}}}-(\mathit{L}-{\hat{\mathit{W}}}_i^g){\mathit{B}}_i-{\mathit{B}}_i(\mathit{L}-{\hat{\mathit{W}}}_i^g{)}^{\mathrm{T}},$$

which, using equations (7) and (8) and some algebra manipulation, may be transformed to:

$${\mathit{F}}^{e^{-1}}{\stackrel{{\nabla }^c}{\mathit{B}}}_i{\mathit{F}}^{e^{-\mathrm{T}}}={\stackrel{·}{\overline{\mathit{B}}}}_i-({\overline{\mathit{L}}}^g-{\mathit{F}}^{e^{-1}}{\hat{\mathit{W}}}_i^g{\mathit{F}}^e){\overline{\mathit{B}}}_i-{\overline{\mathit{B}}}_i({\overline{\mathit{L}}}^g-{\mathit{F}}^{e^{-1}}{\hat{\mathit{W}}}_i^g{\mathit{F}}^e{)}^{\mathrm{T}}.$$

(47)

We may defined an objective rate of ${\overline{\mathit{B}}}_i$ , denoted as ${\stackrel{{\nabla }^g}{\overline{\mathit{B}}}}_i$ , via the relation:

$${\stackrel{{\nabla }^g}{\overline{\mathit{B}}}}_i={\mathit{F}}^{e^{-1}}{\stackrel{{\nabla }^c}{\mathit{B}}}_i{\mathit{F}}^{e^{-\mathrm{T}}},$$

(48)

in which the superscript “ $g$ ” stands for the growth configuration, $\overline{\mathcal{B}}$ . Using this definition and substituting the first relation in equation (26), ${\hat{\mathit{W}}}_i^g={\mathit{FW}}_{0i}^g{\mathit{F}}^{-1}$ , into the r.h.s of equation (47), we obtain an expression,

$${\stackrel{{\nabla }^g}{\overline{\mathit{B}}}}_i={\stackrel{·}{\overline{\mathit{B}}}}_i-({\overline{\mathit{L}}}^g-{\mathit{F}}^g{\mathit{W}}_{0i}^g{\mathit{F}}^{g^{-1}}){\overline{\mathit{B}}}_i-{\overline{\mathit{B}}}_i({\overline{\mathit{L}}}^g-{\mathit{F}}^g{\mathit{W}}_{0i}^g{\mathit{F}}^{g^{-1}}{)}^{\mathrm{T}}.$$

(49)

Alternatively, we may substitute the second relation in equation (26), ${\hat{\mathit{W}}}_i^g={\mathit{W}}_i^g+{\tilde{\mathit{W}}}_i^g$ , into the r.h.s of equation (47) and then obtain the other expression of the objective rate,

$${\stackrel{{\nabla }^g}{\overline{\mathit{B}}}}_i={\stackrel{.}{\overline{\mathit{B}}}}_i-\left({\overline{\mathit{L}}}^g-{\overline{\mathit{C}}}^{-1}({\overline{\mathit{W}}}_i^g+{\widetilde{\overline{\mathit{W}}}}_i^g)\right){\overline{\mathit{B}}}_i-{\overline{\mathit{B}}}_i{\left({\overline{\mathit{L}}}^g-{\overline{\mathit{C}}}^{-1}({\overline{\mathit{W}}}_i^g+{\widetilde{\overline{\mathit{W}}}}_i^g)\right)}^{\mathrm{T}},$$

(50)

where, in consistency with equation (45)₂ and with the help of the definition (27), we define:

$${\widetilde{\overline{\mathit{W}}}}_i^g={\mathit{F}}^{e^{\mathrm{T}}}{\tilde{\mathit{W}}}_i^g{\mathit{F}}^e\equiv {\mathit{F}}^{g^{-\mathrm{T}}}\mathrm{sym}{\mathit{W}}_{0i}^g{\mathit{F}}^{g^{-1}}.$$

(51)

Furthermore, for studying the objective rate on ${\mathcal{B}}_0$ , we consider a tensor associated with the $i\mathrm{th}$ family of CMFs, say ${\mathit{B}}_i^r$ which is defined on ${\mathcal{B}}_0$ and has a relation, ${\mathit{B}}_i={\mathit{FB}}_i^r{\mathit{F}}^{\mathrm{T}}$ , with the tensor ${\mathit{B}}_i$ on ${\mathcal{B}}_t$ .

Substituting the relation, ${\mathit{B}}_i={\mathit{FB}}_i^r{\mathit{F}}^{\mathrm{T}}$ , into equation (25) and going through a process similar to the above discussion from equations (46)–(49) lead to a relation:

$$\stackrel{{\nabla }^r}{{\mathit{B}}_i^r}={\mathit{F}}^{-1}{\stackrel{{\nabla }^c}{\mathit{B}}}_i{\mathit{F}}^{-\mathrm{T}},$$

(52)

where the objective rate of ${\mathit{B}}_i^r$ is obtained as:

$$\stackrel{{\nabla }^r}{{\mathit{B}}_i^r}={\stackrel{·}{\mathit{B}}}_i^r+{\mathit{W}}_{0i}^g{\mathit{B}}_i^r-{\mathit{B}}_i^r{\mathit{W}}_{0i}^g,$$

(53)

in which the superscript “ $r$ ” stands for the remodelling configuration. This is consistent with the definition of ${\mathit{W}}_{0i}^g$ as the remodelling spin on ${\mathcal{B}}_0$ (see equation (20)).

Similar to ${\tilde{\mathit{a}}}_i$ on ${\mathcal{B}}_t$ , ${\overline{\mathit{a}}}_i$ and ${\mathit{a}}_i^r$ are the unit vectors representing the time-dependent fibre directions on $\overline{\mathcal{B}}$ and ${\mathcal{B}}_0$ , respectively. Thus, the kinematical remodelling laws in terms of the corresponding geometrical structural tensors, ${\overline{\mathit{A}}}_i$ and ${\mathit{A}}_i^r$ , are:

$${\stackrel{{\nabla }^g}{\overline{\mathit{A}}}}_i=0\mathrm{and}\stackrel{{\nabla }^r}{{\mathit{A}}_i^r}=0\forall i\in \{1,\dots ,n_f\},$$

(54)

respectively, which may be derived either directly from equation (33) with the transformations (48) and (52) or from the Lie derivative (23) with the definition of the push-forwards,

$${\hat{\mathrm{\chi }}}^{*}(·)=({\mathit{F}}^g{\mathbf{\Omega }}_i^{-1}{)}^{\mathrm{T}}(·)({\mathit{F}}^g{\mathbf{\Omega }}_i^{-1})\mathrm{for}{\overline{\mathit{A}}}_i\mathrm{and}{\hat{\mathrm{\chi }}}^{*}(·)={\mathbf{\Omega }}_i(·){\mathbf{\Omega }}_i^{-1}\mathrm{for}{\mathit{A}}_i^r,$$

(55)

respectively.

Similar to equation (43), equations in (54) give the incremental motion of CMFs. For example, equation (54)₂ gives $\delta {\mathit{A}}_i^r=-{\mathit{W}}_{0i}^g{\mathit{A}}_i^r+{\mathit{A}}_i^r{\mathit{W}}_{0i}^g,$ which is the increment of the structural tensor ${\mathit{A}}_i^r$ on the configuration ${\mathcal{B}}_0$ (or ${\mathcal{B}}_i^r$ ) solely due to the remodelling spin ${\mathit{W}}_{0i}^g$ .

3. The framework of constitutive model

3.1. Elastic constitutive law and the model of growth and remodelling in general forms

With the kinematics of growth and remodelling defined, we consider the constitutive modelling framework covering elastic responses and growth and remodelling behaviours of the cell wall.

The state variables of the cell wall material may be defined on the current configuration ${\mathcal{B}}_t$ in terms of Cauchy stress, $\mathbf{\sigma }$ , and a set of hardening variables consisting of (i) a second-order tensor, $\mathit{\alpha }$ , and a scalar, $k$ , for growth and (ii) a set of second-order tensors, ${\mathit{\alpha }}_i$ , for the remodelling of each and every $i\mathrm{th}$ family of CMFs $(i\in \{1,\dots ,n_f\})$ .

Accordingly, a set of constitutive equations, which are adopted and extended from the viscoplasticity model of metallic materials proposed by Dafalias [48], fully describe the growth and remodelling behaviours of interest:

• 1. The tangential elastic relation of the continuum,

$$\stackrel{◦}{\mathit{\sigma }}=\mathbb{C}:{\mathit{D}}^e,$$

(56)

• 2. The coupled irreversible isochoric expansion and volumetric growth of the continuum and the associated hardening behaviour,

$${\mathit{D}}^g=\mathrm{\Phi }{\mathit{N}}^g,{\mathit{W}}^g=\mathrm{\Phi }{\mathit{M}}^g,$$

(57)

and

$$\stackrel{◦}{\mathit{\alpha }}=\sum _{I=1}^{n_g}\left({\mathrm{\Phi }}_I{\hat{\mathit{\alpha }}}_I\right),\stackrel{·}{k}=\sum _{I=1}^{n_g}\left({\overline{\mathrm{\Phi }}}_I{\hat{k}}_I\right),$$

(58)

• 3. The remodelling of the $i\mathrm{th}$ family of CMFs and the associated hardening behaviour,

$${\mathit{W}}_i^g={\mathrm{\Phi }}^r{\mathit{M}}_i^g\mathrm{and}{\stackrel{{\nabla }^u}{\mathit{\alpha }}}_i=\sum _{I=1}^{n_r}\left({\tilde{\mathrm{\Phi }}}_I{\hat{\mathit{\alpha }}}_{\mathrm{iI}}\right)\forall i\in \{1,\dots ,n_f\}.$$

(59)

The four-order tensor $\mathbb{C}$ represents the tangent elastic moduli of the stress rate $\stackrel{◦}{\mathbf{\sigma }-5}$ which, via equation (65) below, may be derived from hyperelasticity as shown later in section 5.1.3. The functions $\mathrm{\Phi }$ , ${\mathrm{\Phi }}_I$ , ${\overline{\mathrm{\Phi }}}_I$ , ${\mathrm{\Phi }}^r$ , and ${\tilde{\mathrm{\Phi }}}_I$ are the non-negative scalar overstress functions while $n_g$ and $n_r$ are the numbers of the terms of the sequences of the overstress functions for growth and remodelling, respectively. Moreover, ${\mathit{N}}^g$ , ${\mathit{M}}^g$ , ${\hat{\mathit{\alpha }}}_I$ , ${\hat{k}}_I$ , ${\mathit{M}}_i^g$ , and ${\hat{\mathit{\alpha }}}_{\mathrm{iI}}$ define the admissible directions of ${\mathit{D}}^g$ , ${\mathit{W}}^g$ , $\stackrel{◦}{\mathit{\alpha }-5}$ , $\stackrel{·}{k}$ , ${\mathit{W}}_i^g$ , and $\stackrel{{\nabla }^c}{{\mathit{\alpha }}_i}$ , respectively. Noting that ${\mathit{M}}^g$ and ${\mathit{M}}_i^g$ are independent from each other. The invariance requirements under the superposed rigid body rotation render all of the tensor and scalar-valued constitutive functions of equations (56)–(59) the isotropic functions of their variables, $\mathbf{\sigma }$ , $\mathit{\alpha }$ , ${\mathit{\alpha }}_i$ , and $k$ .

The objective rate in equations (56) and (58)₁ is introduced here based on the postulate that the rates of changes of $\mathbf{\sigma }$ and $\mathit{\alpha }$ are in reference to $\mathbf{\omega }$ , instead of $\mathit{W}$ [48]. This may become more clear if we re-formulate equations (56) and (58)₁ using a rate of change with respect to $\mathit{W}$ , i.e., the well-known Jaumann rate. By substituting $\mathit{W}-{\mathit{W}}^g$ for $\mathbf{\omega }$ according to equation (15) and using equations (14), (15), (22), and (56)–(59), equation (56) may be re-expressed as:

$$\stackrel{{\nabla }^J}{\mathbf{\sigma }}=\mathbb{C}:{\mathit{D}}^e-\mathrm{\Phi }({\mathit{M}}^g\mathbf{\sigma }-\mathbf{\sigma }{\mathit{M}}^g),$$

(60)

where:

$$\stackrel{{\nabla }^J}{(·)}=\stackrel{·}{(·)}-\mathit{W}(·)+(·)\mathit{W}\equiv \stackrel{°}{(·)}-{\mathit{W}}^g(·)+(·){\mathit{W}}^g,$$

(61)

is the Jaumann rate. In the same way, hardening law (58)₁ may be re-expressed as:

$$\stackrel{{\nabla }^J}{\mathit{\alpha }}=\sum _{I=1}^{n_r}\left({\mathrm{\Phi }}_I{\hat{\mathit{\alpha }}}_I\right)-\mathrm{\Phi }({\mathit{M}}^g\mathit{\alpha }-\mathit{\alpha }{\mathit{M}}^g).$$

(62)

It is observed that, as the introduction of the growth spin ${\mathit{W}}^g$ , the growth spin direction, ${\mathit{M}}^g$ , plays the role of a structural tensor which introduces the additional anisotropic terms into the constitutive model (60) and (62), where the conventional Jaumann rate is used.

Similarly, using elastic Truesdell rate of Cauchy stress tensor [59],

$$\stackrel{{\nabla }^{eT}}{(•)}=\stackrel{·}{(•)}-\mathrm{tr}{\mathit{L}}^e(•)-{\mathit{L}}^e(•)-(•){\mathit{L}}^{e^{\mathrm{T}}}\equiv \stackrel{◦}{(•)}-\mathrm{tr}{\mathit{D}}^e(•)-{\mathit{D}}^e(•)-(•){\mathit{D}}^e,$$

(63)

equation (56) may also be re-expressed as:

$$\stackrel{{\nabla }^{\mathrm{eT}}}{\mathbf{\sigma }}={\mathbb{C}}^{\sigma T}:{\mathit{D}}^e,$$

(64)

where the elastic tangent moduli for the elastic Truesdell rate is:

$${\mathbb{C}}^{\sigma T}=\mathbb{C}-{\mathbb{C}}^{\prime }-\mathbf{\sigma }\otimes \mathit{I},$$

(65)

in which ${\mathbb{C}}^{\prime }$ is defined by ${\mathbb{C}}^{\prime }:{\mathit{D}}^e={\mathit{D}}^e\mathbf{\sigma }+\mathbf{\sigma }{\mathit{D}}^e$ [59]. Accordingly, the evolution of $\mathit{\alpha }$ in equation (58)₁ is re-written as:

$${\mathit{\alpha }}^{{\nabla }^{\mathbf{eT}}}=\stackrel{◦}{\mathit{\alpha }-4}-\mathrm{tr}({\mathit{D}}^{\mathbf{e}})\mathit{\alpha }-{\mathit{D}}^{\mathbf{e}}\mathit{\alpha }-\mathit{\alpha }{\mathit{D}}^{\mathbf{e}}.$$

(66)

The Truesdell rate is useful where hyperelastic formulation is applied since such a formulation naturally leads to a rate format elastic constitutive law using the Truesdell rate.

3.2. Specification of the hardening laws to the growing cell wall

In the continuum modelling of growth and remodelling, hardening laws should be paid significant attention as they may play a unique role to represent the underlying mechanical and non-mechanical mechanisms of living cell walls. This observation was demonstrated by a mechanism-based hardening model proposed for understanding the dynamic regulation of cell wall growth at micromechanical level [17], where growing and non-growing cell walls were clearly differentiated by one single hardening law. Moreover, the mechanism-based hardening model leaded to a scalar phenomenological hardening model (equation (85) in Huang et al. [17]) which may indicate the proper form of specifications of equations (58) and (59)₂ for the cell wall growth.

Based on such a consideration, we propose a method for the specification of hardening laws as follows. The hardening model proposed by Dafalias (equation (10) in Dafalias [48]) for viscoplasticity is adopted and modified mathematically for the specification of equations (58) and (59)₂, and then the comparison between this type of hardening model and the hardening law proposed in Huang et al. [17] is conducted to show the consistency between them, even though they were proposed originally based on their very own and different mechanisms of the totally different materials. This consistency not only guides and justifies the specification of hardening laws here but also provides the physical interpretation as shown in the following discussion.

In reference to equation (58)₁, it is proposed that $n_g=3$ so that we have $\stackrel{◦}{\mathit{\alpha }}={\sum }_{I=1}^3({\mathrm{\Phi }}_I{\widehat{\mathit{\alpha }}}_I)$ . The detail of each term is discussed as follows:

First, for $I=1$ , it is proposed that ${\hat{\mathit{\alpha }}}_1$ is proportional to (and depending only on) ${\mathit{N}}^g$ , i.e.,

$${\hat{\mathit{\alpha }}}_1=H_1{\mathit{N}}^g,$$

(67)

with the corresponding overstress function ${\mathrm{\Phi }}_1=\mathrm{\Phi }$ , where $\mathrm{\Phi }$ is the overstress function of growth as shown in equation (57).

Second, for $I=2$ , with ${\mathrm{\Phi }}_2=\mathrm{\Phi }$ as well, ${\hat{\mathit{\alpha }}}_2$ is specified in a coupling form as:

$${\hat{\mathit{\alpha }}}_2=H_2\parallel {\mathit{N}}^g\parallel {\widehat{\tilde{\mathit{\alpha }}}}_2(\mathit{\alpha }),$$

(68)

where ${\widehat{\tilde{\mathit{\alpha }}}}_2(\mathit{\alpha })$ is a tensorial function of $\mathit{\alpha }$ , and $\parallel (·)\parallel =\sqrt{\mathrm{tr}[(·{)}^{\mathrm{T}}(·)]}$ is the norm of a tensor (·).

• Third, for the remaining term ${\hat{\mathit{\alpha }}}_3$ , we set:

$${\hat{\mathit{\alpha }}}_3={\widehat{\tilde{\mathit{\alpha }}}}_3(\mathit{\alpha }),$$

(69)

as a tensorial function solely in terms of $\mathit{\alpha }$ . Accordingly, it is requested that ${\mathrm{\Phi }}_3$ is explicitly independent of strain/stress so that ${\mathrm{\Phi }}_3$ is decoupled from $\mathrm{\Phi }$ . For this purpose, ${\mathrm{\Phi }}_3$ may be expressed solely as a functional of $\parallel \mathit{\alpha }\parallel$ , i.e., ${\mathrm{\Phi }}_3={\hat{\mathrm{\Phi }}}_3(\parallel \mathit{\alpha }\parallel ).$ An example of ${\hat{\mathrm{\Phi }}}_3$ using power law is:

$${\hat{\mathrm{\Phi }}}_3(\parallel \mathit{\alpha }\parallel )=H_3\parallel \mathit{\alpha }{\parallel }^{m_3},$$

(70)

where $m_3$ is a constant. $H_3$ as well as $H_1$ in equation (67) and $H_2$ in equation (68) are three scalar parameters.

Hence, using above assumptions of ${\hat{\mathit{\alpha }}}_I$ in equations (67)–(69) and those of ${\mathrm{\Phi }}_I(I\in \{1,2,3\})$ , equation (58)₁ may be specified as:

$$\stackrel{◦}{\mathit{\alpha }}=\underset{(1)\mathrm{strain}\mathrm{hardening}}{\underbrace{H_1{\mathit{D}}^g}}+\underset{(2)\mathrm{dynamic}\mathrm{recovery}}{\underbrace{H_2\parallel {\mathit{D}}^g\parallel {\widehat{\tilde{\mathit{\alpha }}}}_2(\mathit{\alpha })}}+\underset{(3)\mathrm{static}\mathrm{recovery}}{\underbrace{{\hat{\mathrm{\Phi }}}_3(\parallel \mathit{\alpha }\parallel ){\widehat{\tilde{\mathit{\alpha }}}}_3(\mathit{\alpha })}}.$$

(71)

According to equation (71), $\stackrel{◦}{\mathit{\alpha }-5}$ is consisting of three distinct terms. The strain hardening term, $H_1{\mathit{D}}^g$ , and the dynamic recovery term, $H_2\parallel {\mathit{D}}^g\parallel {\widehat{\tilde{\mathit{\alpha }}}}_2(\mathit{\alpha })$ , occur only when ${\mathit{D}}^g\ne 0$ (i.e., $\mathrm{\Phi }>0$ ), while the static recovery term, ${\mathrm{\Phi }}_3{\widehat{\tilde{\mathit{\alpha }}}}_3(\mathit{\alpha })$ , operates always as long as ${\mathrm{\Phi }}_3>0$ .

If we further specify that, ${\widehat{\tilde{\mathit{\alpha }}}}_2(\mathit{\alpha })={\widehat{\tilde{\mathit{\alpha }}}}_3(\mathit{\alpha })=\mathit{\alpha },$ and use the power law (70), equation (71) becomes:

$$\stackrel{◦}{\mathit{\alpha }}=H_1{\mathit{D}}^g+H_2\parallel {\mathit{D}}^g\parallel \mathit{\alpha }+H_3\parallel \mathit{\alpha }{\parallel }^{m_3}\mathit{\alpha }.$$

(72)

In a form similar to equation (72), equation (58)₂ may be proposed to have a specific expression as:

$$\stackrel{·}{k}=H_{k1}\parallel {\mathit{D}}^g\parallel +H_{k2}\parallel {\mathit{D}}^g\parallel k+H_{k3}{k}^{m_{k3}+1},$$

(73)

where $H_{k1}$ , $H_{k2}$ , $H_{k3}$ , and $m_{k3}$ are the scalar parameters.

As aforementioned, in order to clarify the role of each term in the r.h.s of equation (73) (or equations (71) and (72)) and justify the specification of the hardening laws (71)–(73), it is helpful to compare equation (73) with the mechanism-induced hardening law of cell wall growth proposed in equation (85) of Huang et al. [17], which is cited here as:

$$\stackrel{·}{\alpha }=\underset{(1)}{\underbrace{h_1\stackrel{·}{\lambda }}}+\underset{(2)}{\underbrace{h_2\stackrel{·}{\lambda }\alpha }}+\underset{(3)}{\underbrace{0}},$$

(74)

where $\alpha$ and $\stackrel{·}{\lambda }$ are the hardening variable and scalar growth rate, respectively, and $h_1$ and $h_2$ are the two constants. According to the discussion in Huang et al. [17], a comparison between equations (74) and (73) indicates that the first term in the r.h.s of equation (73) (or equations (71) and (72)) represents the strain hardening of the cell wall matrix due to the increase of tension in the crosslink network of the matrix. And the second term represents the effect of softening of the crosslink network due to two contributions: deposition of new load-bearing crosslinks coming in from outside and reconnection of the cut tethers and broken hydrogen bonds inside the cell wall. Hence, the second term features the key regulation mechanisms of growing cell walls as a living material. The third term in the r.h.s of equation (73) (or equations (71) and (72)), which does not have the counterpart in equation (74), may be interpreted as a chemical regulation term which is not explicitly dependent of mechanical regulation as there is no explicit dependency on strain/stress in the expression. Gravitropism provides a typical example demonstrating the coupling between chemical regulation and mechanical regulation of growth [61,62].

Therefore, in the light of equation (74), the specification of hardening laws in equations (71)–(73) can be justified by the physical interpretation. Furthermore, in a form similar to equation (72), equation (59)₂ of the hardening of a remodelling process may be specified as:

$${\stackrel{{\nabla }^u}{\mathit{\alpha }}}_i=H_{i1}{\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i)\frac{{\mathit{n}}_{\alpha }^i}{\parallel {\mathit{n}}_{\alpha }^i\parallel }+H_{i2}{\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i){\mathit{\alpha }}_i+H_{i3}\parallel {\mathit{\alpha }}_i{\parallel }^{m_{i3}}{\mathit{\alpha }}_i,$$

(75)

where the second-order tensor, ${\mathit{n}}_{\alpha }^i$ , is the driving force of remodelling and ${\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i)$ is the corresponding overstress function. The hardening internal variable, ${\mathit{\alpha }}_i$ , may be interpreted as a local back stress created by ${\mathit{W}}_i^g$ . Noting that the objective rate in equation (75) is the same as the one in the kinematical remodelling law (42). This is based on the consideration that ${\mathit{\alpha }}_i$ is associated with ${\tilde{\mathit{A}}}_i$ . More detailed discussion about equation (75) is presented in the later section of a case study.

4. Invariants and spin representation for the growth and remodelling of cell wall

Hitherto, we consider general situations where the number of families of fibres with different direction, $n_f$ , are in theory any finite number. From this section, we shall confine our attention to the specific case that $n_f=2$ . The reason of this is that cell wall may be considered as a laminated composite with a finite number of layers. And each layer may be considered as a lamina composed of one or two families of fibres. According to composite theory, we could model the material property of a representative layer and then obtain the whole multi-layer composite property using the lay-up and orientation of the layers. Therefore, it is reasonable to consider typical materials of two families of fibres.

4.1. Spencer’s deviatoric stress tensors as the driving forces of growth and remodelling

Let us consider a material with two families of fibres $(n_f=2)$ : (i) the two families of fibres could be both of physical CMFs, or (ii) one family is of physical CMFs while another is simply a pure geometric direction (e.g., a director attached to the continuum).

Under the assumption that CMFs are only extensible elastically, it is required that the driving force of growth should vanish in the direction of each family of CMFs [17]. This requirement naturally leads to the application of the so-called Spencer’s deviatoric stress tensor proposed by Spencer [30] and adopted and extended to be applied to cell wall growth modelling by Huang et al. [17,18].

Spencer’s deviatoric (Kirchhoff) stress, ${\mathit{\tau }}^{\mathrm{dev}}$ , may be derived from an expression,

$${\mathit{\tau }}^{\mathrm{dev}}=\mathit{\tau }+\alpha \mathit{I}+\sum _{i=1}^2\left({\beta }_i{\tilde{\mathit{A}}}_i\right),$$

(76)

and its constraints,

$$\mathrm{tr}{\mathit{\tau }}^{\mathrm{dev}}=0\mathrm{and}{\mathit{\tau }}^{\mathrm{dev}}:{\tilde{\mathit{A}}}_i=0\forall i\in \{1,2\}.$$

(77)

Whereby $\mathit{\tau }=\overline{J}\mathbf{\sigma }$ is the Kirchhoff stress while $\overline{J}=\det {\mathit{F}}^e$ . The index “ $i$ ” refers to the $i\mathrm{th}$ family of fibres as aforementioned. The first constraint in equation (77) makes sure ${\mathit{\tau }}^{\mathrm{dev}}$ satisfying the isochoric condition while the second constraint eliminates the stress component in both directions of the two families of CMFs. Hence, ${\mathit{\tau }}^{\mathrm{dev}}$ serves as the proper driving force of the isochoric wall expansion prohibiting inelastic extensions of CMFs.

Solving equations (76) and (77) yields the expression of the coefficients, $\alpha$ and ${\beta }_i$ $(i\in \{1,2\})$ in equation (76) as:

$$\alpha =\frac{1}{3}\mathit{\tau }:[-\mathit{I}+\frac{3}{3c+1}\sum _{i=1}^2\mathrm{dev}{\tilde{\mathit{A}}}_i],$$

(78)

$${\beta }_i=\frac{2}{3c+1}\mathit{\tau }:\left[-\frac{3}{2}\mathrm{dev}{\tilde{\mathit{A}}}_{i^{*}}+\frac{1}{c-1}\left({\tilde{\mathit{A}}}_i-{\tilde{\mathit{A}}}_{i^{*}}\right)\right]\forall i\in \{1,2\},$$

(79)

where $c={\tilde{\mathit{A}}}_1:{\tilde{\mathit{A}}}_2,$ $\mathrm{dev}(·)=(·)-\frac{1}{3}\mathrm{tr}(·)\mathit{I}$ is the deviatoric operator on ${\mathcal{B}}_t$ , and the index $i^{*}$ is the complementary index of an index $i$ in the set $\{1,2\}$ , i.e., $1^{*}=2$ and $2^{*}=1$ (numerically $i^{*}=3-i$ ).

For the purpose of growth modelling using finite element method, it is useful to have the counterpart of ${\mathit{\tau }}^{\mathrm{dev}}$ on $\overline{\mathcal{B}}$ where finite element formulation is more convenient in contrast to that on ${\mathcal{B}}_t$ . Pulling ${\mathit{\tau }}^{\mathrm{dev}}$ from ${\mathcal{B}}_t$ back to $\overline{\mathcal{B}}$ by a transformation,

$${\overline{\mathit{S}}}^{\mathrm{DEV}}={\mathit{F}}^{e^{-1}}{\mathit{\tau }}^{\mathrm{dev}}{\mathit{F}}^{e^{-\mathrm{T}}},$$

(80)

yields the expression of Spencer’s deviatoric PK2 stress, ${\overline{\mathit{S}}}^{\mathrm{DEV}}$ , as:

$${\overline{\mathit{S}}}^{\mathrm{DEV}}=\overline{\mathit{S}}+\overline{\alpha }{\overline{\mathit{C}}}^{-1}+\sum _{i=1}^2\left({\overline{\beta }}_i{\overline{\mathit{A}}}_i\right),$$

(81)

where $\overline{\mathit{S}}={\mathit{F}}^{e^{-1}}\mathit{\tau }{\mathit{F}}^{e^{-\mathrm{T}}}$ is the second Piola–Kirchhoff (PK2) stress on $\overline{\mathcal{B}}$ while $\overline{\alpha }\equiv \alpha$ and ${\overline{\beta }}_i\equiv {\beta }_i$ $(i\in \{1,2\})$ are the coefficients to be expressed in terms of variables on $\overline{\mathcal{B}}$ .

Noting that we have such a relation between the deviatoric parts of ${\overline{\mathit{A}}}_i$ and ${\tilde{\mathit{A}}}_i$ ,

$$\mathrm{DEV}{\overline{\mathit{A}}}_i={\lambda }_{{\mathit{a}}_i}{\mathit{F}}^{e^{-1}}\mathrm{dev}{\tilde{\mathit{A}}}_i{\mathit{F}}^{e^{-\mathrm{T}}},$$

(82)

where $\mathrm{DEV}(·)=(·)-\frac{1}{3}\mathrm{Tr}(·){\overline{\mathit{C}}}^{-1}$ and $\mathrm{Tr}(·)=\overline{\mathit{C}}:(·)$ within its r.h.s are the deviatoric and trace operators on $\overline{\mathcal{B}}$ , respectively. It is observed that ${\lambda }_{{\mathit{a}}_i}=\mathrm{Tr}{\overline{\mathit{A}}}_i=\mathrm{tr}{\mathit{A}}_i$ .

Hence, using relation (82), the scalars, $\overline{\alpha }$ and ${\overline{\beta }}_i$ , in equations (78) and (79) may be expressed in terms of the variables on $\overline{\mathcal{B}}$ as:

$$\overline{\alpha }=\frac{1}{3}\overline{\mathit{S}}:\left\{\overline{\mathit{C}}[-{\overline{\mathit{C}}}^{-1}+\frac{3}{3\overline{c}+1}\sum _{i=1}^2\mathrm{DEV}{\widetilde{\overline{\mathit{A}}}}_i]\overline{\mathit{C}}\right\}\equiv \alpha ,$$

(83)

$${\overline{\beta }}_i=\frac{2}{3\overline{c}+1}\overline{\mathit{S}}:\left\{\overline{\mathit{C}}[-\frac{3}{2}\mathrm{DEV}{\widetilde{\overline{\mathit{A}}}}_{i^{*}}+\frac{1}{\overline{c}-1}({\widetilde{\overline{\mathit{A}}}}_i-{\widetilde{\overline{\mathit{A}}}}_{i^{*}})]\overline{\mathit{C}}\right\}\equiv {\beta }_i\forall i\in \{1,2\},$$

(84)

where $\overline{c}={\widetilde{\overline{\mathit{A}}}}_1:{\widetilde{\overline{\mathit{A}}}}_2\equiv c$ , in which the weighted structural tensor ${\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$ is defined as:

$${\widetilde{\overline{\mathit{A}}}}_i={\widetilde{\overline{\mathit{a}}}}_i\otimes {\widetilde{\overline{\mathit{a}}}}_i\equiv {\lambda }_{{\mathit{a}}_i}^{-1}{\overline{\mathit{A}}}_i,$$

(85)

while ${\widetilde{\overline{\mathit{a}}}}_i={\lambda }_{{\mathit{a}}_i}^{-1/2}{\overline{\mathit{a}}}_i$ is the corresponding weighted vector.

It can be proven that Spencer’s deviatoric PK2 stress ${\overline{\mathit{S}}}^{\mathrm{DEV}}$ satisfies the constraints consistent with equation (77),

$$\mathrm{Tr}{\overline{\mathit{S}}}^{\mathrm{DEV}}=0\mathrm{and}\parallel {\overline{\mathit{S}}}^{\mathrm{DEV}}{\parallel }_{{\widetilde{\overline{\mathit{a}}}}_i}=0\forall i\in \{1,2\},$$

(86)

where $\parallel (·){\parallel }_{{\widetilde{\overline{\mathit{a}}}}_i}=[\overline{\mathit{C}}·(·)·\overline{\mathit{C}}]:{\widetilde{\overline{\mathit{A}}}}_i$ . Furthermore, if the Mandel-stress-style tensors corresponding to $\overline{\mathit{S}}$ and ${\overline{\mathit{S}}}^{\mathrm{DEV}}$ are introduced, respectively, as:

$$\overline{\mathit{\sum }}=\overline{\mathit{S}}\overline{\mathit{C}}\mathrm{and}{\overline{\mathit{\sum }}}^{\mathrm{DEV}}={\overline{\mathit{S}}}^{\mathrm{DEV}}\overline{\mathit{C}},$$

(87)

the constraints in equation (86) may be re-expressed as:

$$\mathrm{tr}{\overline{\mathit{\sum }}}^{\mathrm{DEV}}=0\mathrm{and}\mathrm{Tr}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\widetilde{\overline{\mathit{A}}}}_i)=0\forall i\in \{1,2\}.$$

(88)

Mandel-stress-style tensors play the important role in the constitutive modelling of growth [49] and finite element formulation [59]. We shall see the application of this type of tensors in formulations on $\overline{\mathcal{B}}$ in the subsequent discussion.

4.2. Invariants for the model of growth and remodelling

The purpose of this section is to list the invariants which may be potentially used for constructing yield functions of the growth and remodelling of the cell wall. Yield functions should be in terms of the invariants for their objectivity, i.e., invariance under rigid body motion or observer transformation [63].

Since this study involves a large number of tensors, it is unlikely and unnecessary to list the invariants of all of tensors and their combinations. Therefore, we select a typical set of second-order tensors including (i) a stress tensor, (ii) a geometric structural tensor, and (iii) a hardening tensor to demonstrate the forms of invariants. The result presented in this section may be applied to other set of tensors accordingly. For more discussion about the invariants related to Spencer’s deviatoric stress or general theory of the representation of anisotropic invariants, the reader may refers to the publications [30,64–68].

4.2.1. Invariants on $\mathcal{B}$ _t

In consistency with the aforementioned set of state variables of the constitutive model, the invariants of the continuum on ${\mathcal{B}}_t$ with two families of fibres may be expressed in terms of Kirchhoff stress $\mathit{\tau }$ , the geometrical structural tensors ${\tilde{\mathit{A}}}_i(i\in \{1,2\})$ , and the Spencer-deviatoric-stress-like hardening tensor, $\mathit{\alpha }$ , which is introduced in equation (58)₁ and here is requested to satisfy the conditions:

$$\mathrm{tr}\mathit{\alpha }=0\mathrm{and}\mathit{\alpha }:{\tilde{\mathit{A}}}_i^h=0\forall i\in \{1,2\},$$

(89)

where ${\tilde{\mathit{A}}}_i^h(\mathrm{tr}{\tilde{\mathit{A}}}_i^h=1\forall i\in \{1,2\})$ are the structural tensors characterising the hardening tensor $\mathit{\alpha }$ so they are similar to, but need not to be identical to, ${\tilde{\mathit{A}}}_i$ . We take this typical set of second-order tensors, $\{\mathit{\tau },{\tilde{\mathit{A}}}_i,\mathit{\alpha }\},$ to demonstrate the forms of invariants, which shall be applicable to other set of second-order tensors.

The invariants in terms of $\mathit{\tau }$ , ${\tilde{\mathit{A}}}_i(i\in \{1,2\})$ , and $\mathit{\alpha }$ are listed as follows:

• 1. $\mathit{\tau }$

$$\mathrm{tr}\mathit{\tau },\mathrm{tr}{\mathit{\tau }}^2,\mathrm{tr}{\mathit{\tau }}^3;$$

• 2. $\mathit{\tau },{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$$\mathrm{tr}(\mathit{\tau }{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\tau }}^2{\tilde{\mathit{A}}}_i);\mathrm{tr}(\mathit{\tau }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\tau }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]);$$

• 3. $\mathit{\tau },\mathit{\alpha }$

$$\mathrm{tr}(\mathit{\tau }\mathit{\alpha }),\mathrm{tr}({\mathit{\tau }}^2\mathit{\alpha }),\mathrm{tr}(\mathit{\tau }{\mathit{\alpha }}^2),\mathrm{tr}({\mathit{\tau }}^2{\mathit{\alpha }}^2);$$

• 4. $\mathit{\tau },\mathit{\alpha },{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$$\begin{array}{c}\mathrm{tr}(\mathit{\tau }\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}(\mathit{\tau }{\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\tau }}^2\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\tau }}^2{\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i);\hfill \\ \mathrm{tr}(\mathit{\tau }\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}(\mathit{\tau }{\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\tau }}^2\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\tau }}^2{\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]);\hfill \end{array}$$

• 5. $\mathit{\alpha },{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$$\mathrm{tr}{\mathit{\alpha }}^2,\mathrm{tr}{\mathit{\alpha }}^3;\mathrm{tr}(\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i);\mathrm{tr}(\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]).$$

Alternatively, the invariants may be presented in terms of ${\mathit{\tau }}^{\mathrm{dev}}$ (and $\mathit{\tau }$ where it is needed), ${\tilde{\mathit{A}}}_i(i\in \{1,2\})$ , and $\mathit{\alpha }$ as follows:

• 1. ${\mathit{\tau }}^{\mathrm{dev}}(\mathrm{and}\mathit{\tau })$

$$\mathrm{tr}\mathit{\tau };\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^3);$$

• 2. ${\mathit{\tau }}^{\mathrm{dev}}(\mathrm{and}\mathit{\tau }),{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$$\mathrm{tr}(\mathit{\tau }{\tilde{\mathit{A}}}_i),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2{\tilde{\mathit{A}}}_i);\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]);$$

• 3. ${\mathit{\tau }}^{\mathrm{dev}},\mathit{\alpha }$

$$\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2\mathit{\alpha }),\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}{\mathit{\alpha }}^2),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2{\mathit{\alpha }}^2);$$

• 4. ${\mathit{\tau }}^{\mathrm{dev}},\mathit{\alpha },{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$$\begin{array}{c}\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}{\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2{\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i);\hfill \\ \mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}{\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]);\hfill \\ \mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}(({\mathit{\tau }}^{\mathrm{dev}}{)}^2{\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]);\hfill \end{array}$$

• 5. $\mathit{\alpha },{\tilde{\mathit{A}}}_i(i\in \{1,2\})$

$${\mathit{\alpha }}^2,\mathrm{tr}{\mathit{\alpha }}^3;\mathrm{tr}(\mathit{\alpha }{\tilde{\mathit{A}}}_i),\mathrm{tr}({\mathit{\alpha }}^2{\tilde{\mathit{A}}}_i);\mathrm{tr}(\mathit{\alpha }·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]),\mathrm{tr}({\mathit{\alpha }}^2·[{\tilde{\mathit{a}}}_1\otimes {\tilde{\mathit{a}}}_2]).$$

4.2.2. Invariants on $\overline{\mathcal{B}}$

In a similar way, the invariants may be expressed in terms of variables defined on $\overline{\mathcal{B}}$ , which is useful for finite element formulation. For this purpose, in addition to the geometrical structural tensor ${\widetilde{\overline{\mathit{A}}}}_i$ defined in equation (85) and the Mandel-style stress tensors $\overline{\mathit{\sum }}$ and ${\overline{\mathit{\sum }}}^{\mathrm{DEV}}$ defined in equation (87), a Mandel-stress-like hardening tensor, $\overline{\mathit{\varkappa }}$ , is introduced here as:

$$\overline{\mathit{\varkappa }}=\overline{\mathit{\alpha }}\overline{\mathit{C}},$$

(90)

where:

$$\overline{\mathit{\alpha }}={\mathit{F}}^{e^{-1}}\mathit{\alpha }{\mathit{F}}^{e^{-\mathrm{T}}},$$

(91)

is a PK2-like hardening tensor transformed from $\mathit{\alpha }$ by a pull-back.

Thus, the invariants on $\overline{\mathcal{B}}$ are expressed in terms of ${\overline{\mathit{\sum }}}^{\mathrm{DEV}}$ (and $\overline{\mathit{\sum }}$ where it is needed), ${\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$ , and $\overline{\mathit{\varkappa }}$ as follows:

• 1. ${\overline{\mathit{\sum }}}^{\mathrm{DEV}}(\mathrm{and}\overline{\mathit{\sum }})$

$$\mathrm{tr}\overline{\mathit{\sum }}\equiv \mathrm{Tr}\overline{\mathit{S}}\equiv \mathrm{tr}\mathit{\tau };\mathrm{tr}{({\overline{\mathit{\sum }}}^{\mathrm{DEV}})}^2\equiv \mathrm{tr}{({\mathit{\tau }}^{\mathrm{dev}})}^2,\mathrm{tr}{({\overline{\mathit{\sum }}}^{\mathrm{DEV}})}^3\equiv \mathrm{tr}{({\mathit{\tau }}^{\mathrm{dev}})}^3;$$

• 2. ${\overline{\mathit{\sum }}}^{\mathrm{DEV}}(\mathrm{and}\overline{\mathit{\sum }}),{\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$

$$\begin{array}{l}\mathrm{Tr}(\overline{\mathbf{\sum }}{\widetilde{\overline{\mathit{A}}}}_i)\equiv \mathrm{tr}(\overline{\mathbf{\sum }}{\widetilde{\overline{\mathit{A}}}}_i\overline{\mathit{C}})\equiv \parallel \overline{\mathit{S}}{\parallel }_{{\widetilde{\overline{\mathit{a}}}}_i}\equiv \mathrm{tr}(\mathit{\tau }{\widetilde{\overline{\mathit{A}}}}_i);\mathrm{Tr}\left({({\overline{\mathbf{\sum }}}^{\mathrm{DEV}})}^2{\widetilde{\overline{\mathit{A}}}}_i\right)\equiv \mathrm{tr}\left({({\mathit{\tau }}^{\mathrm{dev}})}^2{\widetilde{\overline{\mathit{A}}}}_i\right),\\ \mathrm{Tr}\left({\overline{\mathbf{\sum }}}^{\mathrm{DEV}}\cdot [{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]\right)\equiv \mathrm{tr}\left({\mathit{\tau }}^{\mathrm{dev}}\cdot [{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]\right),\mathrm{Tr}\left({({\overline{\mathbf{\sum }}}^{\mathrm{DEV}})}^2\cdot [{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]\right)\equiv \mathrm{tr}\left({({\mathit{\tau }}^{\mathrm{dev}})}^2\cdot [{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]\right);\end{array}$$

• 3. ${\overline{\mathit{\sum }}}^{\mathrm{DEV}},\overline{\mathit{\varkappa }}$

  $$\mathrm{tr}({\overline{\mathbf{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}),\mathrm{tr}({\overline{\mathbf{\sum }}}^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}^2),\mathrm{tr}\left({({\overline{\mathbf{\sum }}}^{\mathrm{DEV}})}^2\overline{\mathit{\varkappa }}\right),\mathrm{tr}\left({({\overline{\mathbf{\sum }}}^{\mathrm{DEV}})}^2{\overline{\mathit{\varkappa }}}^2\right);$$

• 4. ${\overline{\mathit{\sum }}}^{\mathrm{DEV}},\overline{\mathit{\varkappa }},{\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$

$$\begin{array}{c}\mathrm{Tr}\left({\overline{\mathit{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}{\widetilde{\overline{\mathit{A}}}}_i\right),\mathrm{Tr}\left({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}^2{\widetilde{\overline{\mathit{A}}}}_i\right),\mathrm{Tr}\left({\left({\overline{\mathit{\sum }}}^{\mathrm{DEV}}\right)}^2\overline{\mathit{\varkappa }}{\widetilde{\overline{\mathit{A}}}}_i\right),\mathrm{Tr}\left({\left({\overline{\mathit{\sum }}}^{\mathrm{DEV}}\right)}^2{\overline{\mathit{\varkappa }}}^2{\widetilde{\overline{\mathit{A}}}}_i\right);\hfill \\ \mathrm{Tr}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]),\mathrm{Tr}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}^2·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]);\hfill \\ \mathrm{Tr}({({\overline{\mathit{\sum }}}^{\mathrm{DEV}})}^2\overline{\mathit{\varkappa }}·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]),\mathrm{Tr}({\left({\overline{\mathit{\sum }}}^{\mathrm{DEV}}\right)}^2{\overline{\mathit{\varkappa }}}^2·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]);\hfill \end{array}$$

• 5. $\overline{\mathit{\varkappa }},{\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$

$$\mathrm{tr}{\overline{\mathit{\varkappa }}}^2,\mathrm{tr}{\overline{\mathit{\varkappa }}}^3;\mathrm{Tr}(\overline{\mathit{\varkappa }}{\widetilde{\overline{\mathit{A}}}}_i),\mathrm{Tr}({\overline{\mathit{\varkappa }}}^2{\widetilde{\overline{\mathit{A}}}}_i);\mathrm{Tr}(\overline{\mathit{\varkappa }}·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]),\mathrm{Tr}({\overline{\mathit{\varkappa }}}^2·[{\widetilde{\overline{\mathit{a}}}}_1\otimes {\widetilde{\overline{\mathit{a}}}}_2]).$$

4.3. Spin representation for the model of growth and remodelling

In essence, spins are the second-order skew-symmetric tensors, which is naturally related to the concept of differential forms (see, e.g., p. 104 in Marsden and Hughes [60] or p. 115 in Schutz [57]). According to the definition of differential forms, given two symmetric second-order tensors, say $\mathit{A}$ and $\mathit{B}$ , we may compute their wedge product,

$$\mathit{A}\wedge \mathit{B}=\mathit{A}\otimes \mathit{B}-\mathit{B}\otimes \mathit{A},$$

(92)

and conduct the tensor contraction as shown in an index format,

$$[\mathit{A}\wedge \mathit{B}{]}_{\mathrm{iklj}}{\delta }_{\mathrm{kl}}=[\mathit{A}\otimes \mathit{B}-\mathit{B}\otimes \mathit{A}{]}_{\mathrm{iklj}}{\delta }_{\mathrm{kl}}=[\mathit{AB}-\mathit{BA}{]}_{\mathrm{ij}}.$$

(93)

It is straightforward to prove that $[\mathit{AB}-\mathit{BA}]$ is a skew-symmetric tensor.

Based on the form in equation (93) and referring to the general format of spin proposed by Dafalias [69], the set of skew-symmetric tensors may be taken as the basis of the representation of the growth spin on ${\mathcal{B}}_t$ is constructed from the set of state variables of the constitutive model, $\{{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_i,\mathit{\alpha }\}$ , listed as follows:

$$\begin{array}{c}{\mathit{\tau }}^{\mathrm{dev}}{\tilde{\mathit{A}}}_i-{\tilde{\mathit{A}}}_i{\mathit{\tau }}^{\mathrm{dev}};({\mathit{\tau }}^{\mathrm{dev}}{)}^2{\tilde{\mathit{A}}}_i-{\tilde{\mathit{A}}}_i({\mathit{\tau }}^{\mathrm{dev}}{)}^2;\hfill \\ {\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }-\mathit{\alpha }{\mathit{\tau }}^{\mathrm{dev}};({\mathit{\tau }}^{\mathrm{dev}}{)}^2\mathit{\alpha }-\mathit{\alpha }({\mathit{\tau }}^{\mathrm{dev}}{)}^2;{\mathit{\tau }}^{\mathrm{dev}}{\mathit{\alpha }}^2-{\mathit{\alpha }}^2{\mathit{\tau }}^{\mathrm{dev}};\hfill \\ \mathit{\alpha }{\mathit{\tau }}^{\mathrm{dev}}{\mathit{\alpha }}^2-{\mathit{\alpha }}^2{\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha };{\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }({\mathit{\tau }}^{\mathrm{dev}}{)}^2-({\mathit{\tau }}^{\mathrm{dev}}{)}^2\mathit{\alpha }{\mathit{\tau }}^{\mathrm{dev}}.\hfill \end{array}$$

Moreover, the corresponding set of tensors on $\overline{\mathcal{B}}$ constructed from the set of state variables, $\{{\overline{\mathit{\sum }}}^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i,\overline{\mathit{\varkappa }}\},$ is as follows:

$$\begin{array}{c}{\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\widetilde{\overline{\mathit{A}}}}_i\overline{\mathit{C}}-\overline{\mathit{C}}{\widetilde{\overline{\mathit{A}}}}_i{\overline{\mathit{\sum }}}^{\mathrm{DEV}};({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2{\widetilde{\overline{\mathit{A}}}}_i\overline{\mathit{C}}-\overline{\mathit{C}}{\widetilde{\overline{\mathit{A}}}}_i({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2;\hfill \\ {\overline{\mathit{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}-\overline{\mathit{\varkappa }}{\overline{\mathit{\sum }}}^{\mathrm{DEV}};({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2\overline{\mathit{\varkappa }}-\overline{\mathit{\varkappa }}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2;{\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}^2-{\overline{\mathit{\varkappa }}}^2{\overline{\mathit{\sum }}}^{\mathrm{DEV}};\hfill \\ \overline{\mathit{\varkappa }}{\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}^2-{\overline{\mathit{\varkappa }}}^2{\overline{\mathit{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }};{\mathit{\tau }}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2-({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2\overline{\mathit{\varkappa }}{\mathit{\tau }}^{\mathrm{DEV}}.\hfill \end{array}$$

The relation between the tensors in the above two sets is the same relation between Kirchhoff stress on ${\mathcal{B}}_t$ and Mandel stress on $\overline{\mathcal{B}}$ . To demonstrate this observation, first we show that for the first tensors in both sets, we have:

$$\underset{\mathrm{Mandel}-\mathrm{style}\mathrm{tensor}}{\underbrace{{\overline{\mathit{\sum }}}^{\mathrm{DEV}}{\widetilde{\overline{\mathit{A}}}}_i\overline{\mathit{C}}-\overline{\mathit{C}}{\widetilde{\overline{\mathit{A}}}}_i{\overline{\mathit{\sum }}}^{\mathrm{DEV}}}}=\underset{\mathrm{PK}2-\mathrm{style}\mathrm{tensor}}{\underbrace{\left[{\mathit{F}}^{e^{-1}}\underset{\mathrm{Kirchhoff}-\mathrm{style}\mathrm{tensor}}{\underbrace{({\mathit{\tau }}^{\mathrm{dev}}{\tilde{\mathit{A}}}_i-{\tilde{\mathit{A}}}_i{\mathit{\tau }}^{\mathrm{dev}})}}{\mathit{F}}^{e^{-\mathrm{T}}}\right]}}\overline{\mathit{C}},$$

(94)

where we use the definition (87) and a relation, ${\widetilde{\overline{\mathit{A}}}}_i={\mathit{F}}^{e^{-1}}{\tilde{\mathit{A}}}_i{\mathit{F}}^{e^{-\mathrm{T}}}$ , derived from equations (34) and (85). Second, using equations (87) and (90), it can be proved that for the third tensors in both sets, we have:

$$\underset{\mathrm{Mandel}-\mathrm{style}\mathrm{tensor}}{\underbrace{{\overline{\mathit{\sum }}}^{\mathrm{DEV}}\overline{\mathit{\varkappa }}-\overline{\mathit{\varkappa }}{\overline{\mathit{\sum }}}^{\mathrm{DEV}}}}=\underset{\mathrm{PK}2-\mathrm{style}\mathrm{tensor}}{\underbrace{\left[{\mathit{F}}^{e^{-1}}\underset{\mathrm{Kirchhoff}-\mathrm{style}\mathrm{tensor}}{\underbrace{({\mathit{\tau }}^{\mathrm{dev}}\mathit{\alpha }-\mathit{\alpha }{\mathit{\tau }}^{\mathrm{dev}})}}{\mathit{F}}^{e^{-\mathrm{T}}}\right]}}\overline{\mathit{C}}.$$

(95)

The growth spin and remodelling spins may be constructed as linear combinations of skew-symmetric tensors in the forms as shown in this section. We shall discuss the application in section 5.6.

4.4. Yield surfaces

As aforementioned, wall yielding or loosening is the necessary condition of wall expansion under turgor pressure [3,13]. The Lockhart equation indicated the threshold of turgor pressure for cell wall growth [14], which was generalised to the concept of yield surface of growth [15–17]. In consistency with the Lockhart equation, which was cited in Huang et al. [17], the yield surface of wall expansion, $f_Y$ , may be defined on ${\mathcal{B}}_t$ in a general form in terms of $\mathit{\tau }$ , ${\mathit{\tau }}^{\mathrm{dev}}$ , $\mathit{\alpha }$ , ${\tilde{\mathit{A}}}_i(i\in \{1,2\})$ as:

$$f_Y\left(\mathit{\tau },{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_i,\mathit{\alpha }\right)=0,$$

(96)

or equivalently on $\overline{\mathcal{B}}$ in terms of $\overline{\mathit{\sum }}$ , ${\overline{\mathit{\sum }}}^{\mathrm{DEV}}$ , $\overline{\mathit{\varkappa }}$ , ${\widetilde{\overline{\mathit{A}}}}_i(i\in \{1,2\})$ as:

$${\overline{f}}_Y(\overline{\mathit{\sum }},{\overline{\mathit{\sum }}}^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i,\overline{\mathit{\varkappa }})=0,$$

(97)

where $f_Y(\mathit{\tau },{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_i,\mathit{\alpha })\equiv {\overline{f}}_Y(\overline{\mathit{\sum }},{\overline{\mathit{\sum }}}^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i,\overline{\mathit{\varkappa }})$ . For the objectivity of a yield surface, the yield surfaces $f_Y$ and ${\overline{f}}_Y$ must be expressed in terms of the invariants as discussed above. Therefore, mathematically, what $f_Y(\mathit{\tau },{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_i,\mathit{\alpha })={\overline{f}}_Y(\overline{\mathit{\sum }},{\overline{\mathit{\sum }}}^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i,\overline{\mathit{\varkappa }})$ exactly implies is that $f_Y$ and ${\overline{f}}_Y$ are the identical functions in terms of the invariants which may be expressed in terms of variables on ${\mathcal{B}}_t$ and $\overline{\mathcal{B}}$ , respectively. For example, $f_Y(\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}{)}^2)\equiv {\overline{f}}_Y(\mathrm{tr}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2)$ where $\mathrm{tr}({\mathit{\tau }}^{\mathrm{dev}}{)}^2\equiv \mathrm{tr}({\overline{\mathit{\sum }}}^{\mathrm{DEV}}{)}^2$ .

An anisotropic yield surface with isotropic hardening law was reported in the cell wall growth model in Huang et al. [17]. The case study discussed subsequently uses an anisotropic yield surface adopted from the previous studies [17,30] with an isotropic hardening law for the wall growth while using a kinematic hardening law for the remodelling of CMFs. This gives a thorough demonstration of the proposed growth spin theory.

Noting that such a kinematic hardening laws requires a method to compute the objective rate of Spencer’s deviatoric stress tensor ${\mathit{\tau }}^{\mathrm{dev}}$ (and of $\mathit{\alpha }$ ), which we discuss in details in Appendix 1 (The objective rate of Spencer’s deviatoric stress tensor).

5. Case study: cell wall as a cylindrical wall of elastically incompressible material

A cell wall may be simplified as a closed cylindrical wall of elastically incompressible material subjected to an inner (turgor) pressure, $P$ (see Figure 5). The assumption of incompressibility is due to the property of space-filling pectins [70]. The general theory about the finite deformation elasticity of cylindrical wall structures may be referred to Ogden [63] while Soldatos [71] provides the extensive discussion about the mass growth of axisymmetric fibre-reinforced tubes.

Figure 5.

The model of cell wall as a cylindrical wall with two mechanically equivalent families of CMFs. The fibre angle $\vartheta (r,t)$ represents the CMFs’ temporal–spatial distribution.

In such an idealised geometry model, cell wall growth may include the expansive growth on the cylindrical side surface and tip growth on two caps of the closed wall. As the interest of this study is on the cell wall elongation in the longitudinal direction by expansive growth, the effect of tip growth can be ignored when compared to the wall elongation according to the experimental observation [9,35] and the theoretical analysis [15]. Therefore, the model may be further simplified to be an open cylindrical wall subjected to (i) the inner (turgor) pressure $P$ on the wall surface and (ii) the equivalent loads, $N$ , on two end sections of the open wall obtained by lumping the turgor pressure $P$ applied on the cap of the closed wall and then transferring to the end section of the open wall (see Figure 5).

The wall is considered as a laminated composite with a finite number of layers. From the point of view of composite lay-up, it is assumed that the wall has a balanced lay-up of the layers (say, [−45/45]) which, however, is not necessary to be symmetrical with respect to the centroid of the composite shape. In other words, the lay-up on the two sides of the wall cross-section centroid need not to be symmetrical as indicated in the literature [2,72]. Thus, we may assume that (i) each layer consists of two families of CMFs which are mechanically equivalent locally (so that a double helical distribution globally) and (ii) the whole structure is axisymmetric (see Figure 5). It is also assumed that the lay-up throughout the wall, i.e., the distribution of the fibre angles of CMFs along the wall thickness direction, is known. As reported in the literature [2,72], the CMFs in the innermost layer of the wall, which are the latest deposited by wall growth, are aligned in the circumferential direction of the wall. From the innermost surface to the outermost surface along the wall thickness direction, the distribution of the fibre directions may sequentially change from the circumferential direction towards a direction closer to the longitudinal direction of the wall, which is an interesting non-symmetrical lay-up from the points of view of growth modelling and mechanical design.

5.1. Governing equations

5.1.1. Kinematical equations

In the reference configuration ${\mathcal{B}}_0$ , the cylindrical geometry of the wall is defined in terms of the cylindrical coordinates $(R,\mathrm{\Theta },Z)$ by:

$$A\le R\le B,0\le \mathrm{\Theta }\le 2\mathrm{\pi },0\le Z\le L,$$

where $A$ and $B$ are the inner and outer radii, respectively, of the wall on ${\mathcal{B}}_0$ and $L$ is its length.

In terms of the spatial coordinate system, $(r,\theta ,z)$ , the geometry of the wall on the current configuration ${\mathcal{B}}_t$ experiencing growth and elastic deformation is given by (see Figure 5):

$$r=\mathrm{Rf}(R,t),\theta =\mathrm{\Theta },z={\lambda }_3(t)Z,$$

(98)

where $f$ is an unknown function that depends on both the growth and the elastic deformation, and ${\lambda }_3$ is the stretch in the longitudinal direction $Z$ at the current time $t$ . Thus, we have:

$$a(A,t)\le r(R,t)\le b(B,t),0\le z(t)\le l(t)={\lambda }_3(t)L,$$

where the current inner and outer radii, $a$ and $b$ , of the wall are unknown unless specified as part of boundary conditions.

Let $\{{\mathit{E}}_I\}$ and $\{{\mathit{e}}_i\}$ $(I,i\in \{1,2,3\})$ be the orthonormal bases of the cylindrical coordinate systems on ${\mathcal{B}}_0$ and ${\mathcal{B}}_t$ , respectively, where ${\mathit{E}}_3={\mathit{e}}_3$ are the unit vectors in the longitudinal principal direction of the wall. Hence, equation (98) may be expressed in a vector format as:

$$\mathit{x}=f(R,t)\mathcal{P}°(\mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3)+{\lambda }_3(t)(\mathit{X}·{\mathit{E}}_3){\mathit{e}}_3,$$

(99)

where $\mathcal{P}°(*)$ is the projection of a vector from $\{{\mathit{E}}_i\}$ to $\{{\mathit{e}}_i\}$ . The deformation gradient $\mathit{F}=\partial \mathit{x}/\partial \mathit{X}$ is then obtained as:

$$\mathit{F}=\frac{\partial \mathit{x}}{\partial \mathit{X}}=\mathcal{P}°(\mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3)\otimes \frac{\partial f(R,t)}{\partial \mathit{X}}+f(R)(\hat{\mathit{I}}-{\mathit{e}}_3\otimes {\mathit{E}}_3)+{\lambda }_3(t){\mathit{e}}_3\otimes {\mathit{E}}_3,$$

(100)

in which $[\hat{\mathit{I}}{]}_{\mathrm{iJ}}={\delta }_{\mathrm{iJ}}$ is an unit two-point tensor. From the expression, $R=\parallel \mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3\parallel$ , the derivative in the r.h.s of equation (100) can be worked out as:

$$\frac{\partial f(R,t)}{\partial \mathit{X}}=\frac{f^{\prime }(R,t)}{R}(\mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3),$$

(101)

where a prime indicates the derivative with respect to $R$ . So that the expression of $\mathit{F}$ in equation (100) becomes:

$$\begin{array}{c}\mathit{F}=R{f}^{\prime }(R,t)\frac{\mathcal{P}°(\mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3)}{R}\otimes \frac{\mathit{X}-(\mathit{X}·{\mathit{E}}_3){\mathit{E}}_3}{R}\hfill \\ +f(R,t)(\hat{\mathit{I}}-{\mathit{e}}_3\otimes {\mathit{E}}_3)+{\lambda }_3(t){\mathit{e}}_3\otimes {\mathit{E}}_3.\hfill \end{array}$$

(102)

Thus, $\mathit{F}$ is symmetric and its principal axes coincide with the cylindrical axes. The corresponding principal stretches of $\mathit{F}$ , ${\lambda }_1$ , ${\lambda }_2$ , and ${\lambda }_3$ , are then simply read off as:

$${\lambda }_1=R{f}^{\prime }(R,t)+f(R,t),{\lambda }_2=f(R,t)=\frac{r}{R},{\lambda }_3={\lambda }_3(t).$$

(103)

5.1.2. Equilibrium equations

The equilibrium equation on ${\mathcal{B}}_t$ without body force is:

$$\nabla ·\mathbf{\sigma }=\frac{\partial \mathbf{\sigma }}{\partial \mathit{x}}=0,$$

(104)

where ∇ is the gradient operator on ${\mathcal{B}}_t$ . Similar to equation (102), the Cauchy stress $\mathbf{\sigma }$ in the $\{{\mathit{e}}_i\}$ coordinate system may be expressed as:

$$\mathbf{\sigma }=({\sigma }_1-{\sigma }_2)\frac{\mathit{x}-(\mathit{x}·{\mathit{e}}_3){\mathit{e}}_3}{r}\otimes \frac{\mathit{x}-(\mathit{x}·{\mathit{e}}_3){\mathit{e}}_3}{r}+{\sigma }_2(\mathit{I}-{\mathit{e}}_3\otimes {\mathit{e}}_3)+{\sigma }_3{\mathit{e}}_3\otimes {\mathit{e}}_3,$$

(105)

where ${\sigma }_i(i\in \{1,2,3\})$ are the principal stresses of $\mathbf{\sigma }$ . The equilibrium equation in this case reads:

$$\frac{\mathrm{d}{\sigma }_1}{\mathrm{d}r}+\frac{1}{r}({\sigma }_1-{\sigma }_2)=0.$$

(106)

5.1.3. Elastic constitutive law of the growing cell wall

Let $\mathrm{\Psi }$ and ${\mathrm{\Psi }}_0$ be the free energy per unit mass on the configuration ${\mathcal{B}}_t$ and the free energy per unit fixed volume on ${\mathcal{B}}_0$ , respectively. We have the relation:

$${\mathrm{\Psi }}_0={\rho }_0\mathrm{\Psi },$$

(107)

where ${\rho }_0=J\rho$ is the mass density on ${\mathcal{B}}_0$ while $\rho$ is the mass density on ${\mathcal{B}}_t$ . The mass balance is simply expressed as:

$${\stackrel{·}{\rho }}_0(\mathit{X},t)=M_V\mathrm{or}\mathrm{equivalently}\stackrel{·}{\rho }(\mathit{x},t)=\frac{\mathrm{D}\rho (\mathit{x},t)}{\mathrm{D}t}=m_V-\rho \mathrm{tr}\mathit{D},$$

(108)

where $M_V$ and $m_V=J^{-1}{M}_V$ are the mass sources per unit volume on ${\mathcal{B}}_0$ and ${\mathcal{B}}_t$ , respectively, representing the effects of material deposition, transport, and integration, whereby mass fluxes are not considered explicitly for the reasons mentioned in section 2.1. Equation (108)₁ indicates that the density ${\rho }_0(\mathit{X},t)$ within ${\mathcal{B}}_0$ changes with growth, just as the fibre direction ${\mathit{a}}_i^r(\mathit{X},t)$ on ${\mathcal{B}}_0$ changes with remodelling.

The assumption of constant density on $\overline{\mathcal{B}}$ , which is one of the key features of kinematical growth model [49], may be expressed as:

$${\stackrel{·}{\rho }}_g\equiv \overline{J}(\frac{\mathrm{D}\rho }{\mathrm{D}t}+\rho \mathrm{tr}{\mathit{D}}^e)=0,$$

(109)

where ${\rho }_g=\overline{J}\rho$ is the mass density on $\overline{\mathcal{B}}$ . Substituting equation (109) into equation (108)₂ yields:

$$\mathrm{tr}{\mathit{D}}^g=\frac{m_V}{\rho }.$$

(110)

Equations (109) and (110) indicate that (i) the mass source purely contributes to the volumetric part of growth, $\mathrm{tr}{\mathit{D}}^g$ and (ii) the change of $\rho$ is purely elastic.

Based on this observation, a specification of ${\mathrm{\Psi }}_0$ may be introduced as:

$${\mathrm{\Psi }}_0({\mathit{F}}^e,{\mathit{F}}^g;{\mathit{a}}_i^r)=J_g\overline{\mathrm{\Psi }}({\mathit{F}}^e;{\overline{\mathit{a}}}_i),$$

(111)

where $J_g=det{\mathit{F}}^g$ , and $\overline{\mathrm{\Psi }}$ is the free energy density per unit volume on the configuration $\overline{\mathcal{B}}$ so that $\overline{\mathrm{\Psi }}={\rho }_g\mathrm{\Psi }$ . Here, we are concerned only with the volumetric part of growth (110) and the remodelling due to the additional fibre rotation ${\mathbf{\Omega }}_i^{-1}$ , so that the material response characterised by $\overline{\mathrm{\Psi }}$ is independent of ${\mathit{F}}^g$ . Moreover, the implications of objectivity require that $\overline{\mathrm{\Psi }}$ depends on ${\mathit{F}}^e$ through the right Cauchy–Green tensor $\overline{\mathit{C}}={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{F}}^e$ and the free energy function should be in terms of the invariants of a set of tensors, $\{\overline{\mathit{C}},{\overline{\mathit{A}}}_i(i\in 1,2)\}$ [17,33].

In the case of interest that the elastic response of the material is incompressible, i.e., $\overline{J}\equiv 1$ , a Lagrange multiplier, $p$ , is introduced to replace $\overline{\mathrm{\Psi }}$ by $\overline{\mathrm{\Psi }}-p\ln \overline{J}$ . Then, it can be deduced in a standard way (see, e.g., Ogden [63]) that the elastic constitutive equation of such a material is in a form as:

$$\mathbf{\sigma }={\mathit{F}}^e\frac{\partial \overline{\mathrm{\Psi }}}{\partial {\mathit{F}}^e}-p{\overline{J}}^{-1}\mathit{I}={\mathit{F}}^e\frac{\partial \overline{\mathrm{\Psi }}}{\partial {\mathit{F}}^e}-p\mathit{I}.$$

(112)

The corresponding tangent moduli ${\mathbb{C}}^{\sigma T}$ in equation (64) may be computed from the transformation,

$$[{\mathbb{C}}^{\sigma T}{]}_{\mathrm{ijkl}}=4{\overline{J}}^{-1}{F}_{\mathrm{iI}}^e{F}_{\mathrm{jJ}}^e{F}_{\mathrm{kK}}^e{F}_{\mathrm{lL}}^e{\left[\frac{{\partial }^2(\overline{\mathrm{\Psi }}-p\ln \overline{J})}{\partial \overline{\mathit{C}}\partial \overline{\mathit{C}}}\right]}_{\mathrm{IJKL}},$$

(113)

where $4{\partial }^2(\overline{\mathrm{\Psi }}-p\ln \overline{J})/\partial \overline{\mathit{C}}\partial \overline{\mathit{C}}$ is the elastic tangent moduli on $\overline{\mathcal{B}}$ . More details of the expressions of Cauchy stress and tangent moduli of the cell wall model may be referred to Huang et al. [17], which was adopted from Gasser et al. [33].

5.2. Elastic solution of Cauchy stress

Let us consider part of the wall between radii $a$ and $r$ in ${\mathcal{B}}_t$ (see Figure 5). The volume of this part is:

$$\mathrm{\pi }(r^2-a^2)l.$$

(114)

Due to the elastic incompressibility, the elastic response from $\overline{\mathcal{B}}$ to ${\mathcal{B}}_t$ is volume-preserving, and the corresponding material volume on $\overline{\mathcal{B}}$ should be equal to the one in equation (114). However, since the deformation gradient ${\mathit{F}}_g$ is in general incompatible and the configuration $\overline{\mathcal{B}}$ is only a collection of infinitesimal domains, the corresponding volume on $\overline{\mathcal{B}}$ cannot be calculated with a similar formula as equation (114). We denote this volume by $V(R,L,t)$ , and calculate it by a push-forward from ${\mathcal{B}}_0$ to $\overline{\mathcal{B}}$ with the Jacobian $J_g$ :

$$V(R,L,t)=2\mathrm{\pi }L{\int }_A^R\xi J_g(\xi ,t)\mathrm{d}\xi \equiv \mathrm{LS}(R,t),$$

(115)

where:

$$S(R,t)=2\mathrm{\pi }{\int }_A^R\xi J_g(\xi ,t)\mathrm{d}\xi .$$

(116)

The elastic incompressibility indicates that the two volumes in equations (114) and (115) are equal,

$$\mathrm{\pi }(r^2-a^2)l=V(R,L,t),$$

and hence,

$$r^2=a^2+\frac{V(R,L,t)}{\pi l}=a^2+\frac{S(R,t)}{\pi {\lambda }_3}.$$

(117)

The stretch ${\lambda }_2=f(R,t)=r/R$ is then obtained as:

$${\lambda }_2=\frac{1}{R}{(a^2+\frac{S(R,t)}{\mathrm{\pi }{\lambda }_3})}^{1/2}.$$

(118)

Consequently, the relations in equation (103) indicate that:

$${\lambda }_1=\frac{S^{\prime }(R,t)}{2\mathrm{\pi }R{\lambda }_2{\lambda }_3},\mathrm{or}\mathrm{equivalently}{\lambda }_1{\lambda }_2{\lambda }_3=J_g(R,t).$$

(119)

It is observed that, ${\lambda }_3$ (or equivalently the length $l$ ) and the inner radius of the wall, $a$ , in equations (118) and (119), are yet unknown when ${\mathit{F}}^g$ is known once the growth law is specified. To determine $a$ and ${\lambda }_3$ , we need to obtain an expression of Cauchy stress $\mathbf{\sigma }$ in terms of the elastic principle stretches and then to solve the problem using the equilibrium equation and boundary conditions.

In this case, ${\mathit{F}}^e$ is diagonal with the principal stretches denoted by ${\lambda }_{e1}(={\lambda }_{e2}^{-1}{\lambda }_{e3}^{-1})$ , ${\lambda }_{e2}$ , and ${\lambda }_{e3}$ . The free energy $\overline{\mathrm{\Psi }}$ depends only on the stretch ${\lambda }_{\mathrm{ej}}(j\in \{1,2,3\})$ , and we represent this dependence as:

$$\overline{\mathrm{\Psi }}={\overline{\mathrm{\Psi }}}_{*}({\lambda }_{e1}(={\lambda }_{e2}^{-1}{\lambda }_{e3}^{-1}),{\lambda }_{e2},{\lambda }_{e3})=\hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3}).$$

(120)

The reason of choosing ${\lambda }_{e2}$ and ${\lambda }_{e3}$ as the major variables is that the major driving force of growth is considered in an in-plane direction.

Hence, substituting equation (120) into equation (112) yields:

$${\sigma }_j={\lambda }_{\mathrm{ej}}\frac{\partial {\overline{\mathrm{\Psi }}}_{*}}{\partial {\lambda }_{\mathrm{ej}}}-p\forall j\in \{1,2,3\}.$$

(121)

Using equation (121) and the incompressibility condition, it can be proven that the derivatives of $\hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})$ with respect to the major variables are:

$${\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}={\sigma }_2-{\sigma }_1$$

(122)

and

$${\lambda }_{e3}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e3}}={\sigma }_3-{\sigma }_1.$$

(123)

Substituting equation (122) into the equilibrium equation (106) yields:

$$\frac{\mathrm{d}{\sigma }_1}{\mathrm{d}r}=\frac{{\lambda }_{e2}}{r}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}.$$

(124)

${\sigma }_1$ must satisfy the traction boundary conditions: (1) ${\sigma }_1=P$ on the inner surface, $r=a$ and (2) ${\sigma }_1=0$ on the outer surface, $r=b$ . By integrating equation (124) with respect to $r$ from $a$ to $b$ and using the two-point traction boundary conditions of ${\sigma }_1$ , we obtain:

$$P={\int }_a^b{\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}\frac{\mathrm{d}r}{r}.$$

(125)

Equation (117) indicates that:

$$r={(a^2+\frac{S(R)}{\mathrm{\pi }{\lambda }_3})}^{1/2},b={(a^2+\frac{S(B)}{\mathrm{\pi }{\lambda }_3})}^{1/2}.$$

(126)

Hence, $b$ in the integral (125) is determined by $a$ and ${\lambda }_3$ if the growth is specified.

In addition to equation (125), we need the equilibrium condition in $Z$ (longitudinal) direction to completely determine $a$ and ${\lambda }_3$ . For this purpose, integrating equation (124) from $a$ to $r$ yields:

$${\sigma }_1(r)={\int }_a^r{\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}\frac{\mathrm{d}\xi }{\xi }.$$

(127)

Substituting equation (127) into equation (123) leads to the expression of ${\sigma }_3$ ,

$${\sigma }_3={\lambda }_{e3}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e3}}+{\int }_a^r{\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}\frac{\mathrm{d}\xi }{\xi }.$$

(128)

Hence, the resultant force, $N^{*}$ , on the end section of the open wall in the longitudinal direction, is obtained by integrating ${\sigma }_3$ over the end section area as:

$$N^{*}{=}^{\mathrm{def}}2\mathrm{\pi }{\int }_a^b{\sigma }_{3}r\mathrm{d}r=2\mathrm{\pi }{\int }_a^b\{{\lambda }_{e3}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e3}}+{\int }_a^r{\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}\frac{\mathrm{d}\xi }{\xi }\}r\mathrm{d}r.$$

(129)

In addition, given the turgor pressure $P$ , the total force applying on the end cap of the closed wall is $N=\mathrm{\pi }{a}^{2}P$ , which is equivalent to the resultant force $N^{*}$ in equation (129) (refer to Figure 5):

$$N^{*}=N.$$

(130)

Given the growth (i.e., ${\mathit{F}}^g$ ), the equilibrium conditions (125) and (129) together with equations (126) and (130) are sufficient to determine $a$ and ${\lambda }_3$ , which completely solves the elastic deformation of the wall. Consequently, equations (127), (128), and (122), which is re-written in a form as:

$${\sigma }_2={\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}+{\int }_a^r{\lambda }_{e2}\frac{\partial \hat{\mathrm{\Psi }}({\lambda }_{e2},{\lambda }_{e3})}{\partial {\lambda }_{e2}}\frac{\mathrm{d}\xi }{\xi },$$

(131)

provide the principal stresses ${\sigma }_1,{\sigma }_2,{\sigma }_3$ of $\mathbf{\sigma }$ in terms of ${\lambda }_{e2}$ and ${\lambda }_{e3}$ .

5.3. Specification of the expressions of τ^dev

With $\mathbf{\sigma }$ (and then $\mathit{\tau }=\overline{J}\mathbf{\sigma }=\mathrm{Diag}[{\tau }_1,{\tau }_2,{\tau }_3]$ ) obtained, Spencer’s deviatoric stress tensor can be calculated according to equations (76), (78), and (79).

As aforementioned, it is assumed that each layer of the wall has two families of CMFs which are mechanically equivalent. This means the driving force of growth of such a layer as a continuum should vanish in both directions of the two families of CMFs due to the inertia of CMFs in growth. However, each family has its own remodelling behaviour depending on its fibre direction and the axis of re-orientation. The driving force of the remodelling should vanish in these two directions as well since any force in these two directions does not contribute to the re-orientation of this family of CMFs and should be subtracted from the total stress.

Therefore, we consider two specifications of the expression of ${\mathit{\tau }}^{\mathrm{dev}}$ : one with two physically equivalent families of CMFs, the other with one family of CMFs and one pure geometry direction (i.e., the axis of re-orientation). The former is for the modelling of the growth of the wall as a continuum, and the latter is for the modelling of the remodelling of each and every family of CMFs. We discuss the two specifications in details separately.

5.3.1. Specification 1: two mechanically equivalent families of CMFs for the growth of the cell wall

Let us write the radial axis of the wall as $\mathit{r}=(\mathit{x}-(\mathit{x}·{\mathit{e}}_3){\mathit{e}}_3)/\parallel \mathit{x}-(\mathit{x}·{\mathit{e}}_3){\mathit{e}}_3\parallel$ and the circumferential axis as $\mathit{t}={\mathit{e}}_3\times \mathit{r}$ , where × is the cross product of vectors. Since the basis of a local coordinate system, $\{\mathit{r},\mathit{t},{\mathit{e}}_3\}$ , is identical to the principal directions of stress, ${\mathit{\tau }}^{\mathrm{dev}}$ may be expressed in a diagonal form in this local coordinate system.

Using the local coordinate system $\{\mathit{r},\mathit{t},{\mathit{e}}_3\}$ , the two families of CMFs of the representative layer are represented by two unit vectors:

$${\tilde{\mathit{a}}}_i^f=\{0,\cos {\vartheta }_i,\sin {\vartheta }_i{\}}^{\mathrm{T}}i\in \{1,2\},$$

(132)

where ${\vartheta }_i$ is the angle of $i\mathrm{th}$ family of CMFs on the $\mathit{t}\wedge {\mathit{e}}_3$ plane (i.e., the plane expanded by $\mathit{t}$ and ${\mathit{e}}_3$ ) measured from the circumferential axis $\mathit{t}$ towards the longitudinal axis ${\mathit{e}}_3$ around the axis $\mathit{r}$ .

The two mechanically equivalent families of CMFs may be specified using a single specific angle $\vartheta (r,t)\in [0,\mathrm{\pi }/2]$ ,

$${\vartheta }_1=\vartheta \mathrm{and}{\vartheta }_2=\mathrm{\pi }-\vartheta .$$

(133)

Thus, equation (132) becomes:

$${\tilde{\mathit{a}}}_1^f=\{0,\cos \vartheta ,\sin \vartheta {\}}^{\mathrm{T}}\mathrm{and}{\tilde{\mathit{a}}}_2^f=\{0,-\cos \vartheta ,\sin \vartheta {\}}^{\mathrm{T}}.$$

(134)

Let us assign the values of the two specific vectors ${\tilde{\mathit{a}}}_i^f$ to the general unit vectors ${\tilde{\mathit{a}}}_i$ in equation (1),

$${\tilde{\mathit{a}}}_i:={\tilde{\mathit{a}}}_i^f\forall i\in \{1,2\},$$

(135)

which characterise Spencer’s deviatoric stress in equation (76) via the structural tensor ${\tilde{\mathit{A}}}_i={\tilde{\mathit{a}}}_i\otimes {\tilde{\mathit{a}}}_i$ in equation (34). Consequently, the components of the structural tensor ${\tilde{\mathit{A}}}_i(i\in \{1,2\})$ are obtained as:

$${\tilde{\mathit{A}}}_1=\frac{1}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+\cos 2\vartheta \hfill & \sin 2\vartheta \hfill \\ \multicolumn{1}{c}{0}& \sin 2\vartheta \hfill & 1-\cos 2\vartheta \hfill \end{array}\right],{\tilde{\mathit{A}}}_2=\frac{1}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+\cos 2\vartheta \hfill & -\sin 2\vartheta \hfill \\ \multicolumn{1}{c}{0}& -\sin 2\vartheta \hfill & 1-\cos 2\vartheta \hfill \end{array}\right].$$

(136)

The corresponding deviatoric tensors are:

$$\mathrm{dev}{\tilde{\mathit{A}}}_1=\frac{1}{6}\left[\begin{array}{ccc}-2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+3\cos 2\vartheta \hfill & 3\sin 2\vartheta \hfill \\ \multicolumn{1}{c}{0}& 3\sin 2\vartheta \hfill & 1-3\cos 2\vartheta \hfill \end{array}\right]$$

(137)

and

$$\mathrm{dev}{\tilde{\mathit{A}}}_2=\frac{1}{6}\left[\begin{array}{ccc}-2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+3\cos 2\vartheta \hfill & -3\sin 2\vartheta \hfill \\ \multicolumn{1}{c}{0}& -3\sin 2\vartheta \hfill & 1-3\cos 2\vartheta \hfill \end{array}\right].$$

(138)

Substituting equation (136) into equation (76) yields a specification of ${\mathit{\tau }}^{\mathrm{dev}}$ for the continuum (i.e., matrix and CMFs) denoted as:

$${\mathit{\tau }}_c^{\mathrm{dev}}=\mathrm{Diag}[{\tau }_{c11}^{\mathrm{dev}},{\tau }_{c22}^{\mathrm{dev}},{\tau }_{c33}^{\mathrm{dev}}],$$

(139)

where the subscript “ $c$ ” stands for the continuum, and the diagonal components are:

$${\tau }_{c11}^{\mathrm{dev}}={\tau }_1+\alpha ,{\tau }_{c22}^{\mathrm{dev}}={\tau }_2+\alpha +{\beta }_1(1+\cos 2\vartheta ),{\tau }_{c33}^{\mathrm{dev}}={\tau }_3+\alpha +{\beta }_1(1-\cos 2\vartheta ).$$

(140)

Using equations (137) and (138) and $c={\tilde{\mathit{A}}}_1:{\tilde{\mathit{A}}}_2=({\tilde{\mathit{a}}}_1·{\tilde{\mathit{a}}}_2{)}^2=\stackrel{2}{\cos }2\vartheta ,$ the three coefficients in equations (78) and (79) can be worked out as:

$$\alpha =-\frac{1}{3}\mathrm{tr}\mathit{\tau }+\frac{2}{3c+1}\frac{1}{3}(-{\tau }_1+\frac{{\tau }_2+{\tau }_3}{2}+\frac{{\tau }_2-{\tau }_3}{2}3\cos 2\vartheta ),$$

(141)

$${\beta }_1={\beta }_2=-\frac{1}{3c+1}(-{\tau }_1+\frac{{\tau }_2+{\tau }_3}{2}+\frac{{\tau }_2-{\tau }_3}{2}3\cos 2\vartheta ).$$

(142)

Substituting equations (141) and (142) into equation (140) leads to the expression of the diagonal components (no summation of the index $i$ below):

$${\tau }_{\mathrm{cii}}^{\mathrm{dev}}=c_{i1}({\tau }_2-{\tau }_3)+c_{i2}({\tau }_3-{\tau }_1)+c_{i3}({\tau }_1-{\tau }_2)\forall i\in \{1,2,3\},$$

(143)

where the coefficients, $c_{\mathrm{ij}}(i,j\in \{1,2,3\})$ , are:

$$\begin{array}{c}{c}_{11}=\frac{\cos 2\vartheta }{3c+1},c_{12}=-\frac{c}{3c+1},c_{13}=\frac{c}{3c+1};\hfill \\ c_{21}=\frac{1}{6}(2-3B_1\cos 2\vartheta ),c_{22}=-\frac{1}{6}{B}_1,c_{23}=-\frac{1}{6}(2-B_1);\hfill \\ c_{31}=-\frac{1}{6}(2+3B_2\cos 2\vartheta ),c_{32}=\frac{1}{6}(2-B_2),c_{33}=\frac{1}{6}{B}_2;\hfill \end{array}$$

in which:

$$B_1=\frac{1+3\cos 2\vartheta }{3c+1},B_2=\frac{1-3\cos 2\vartheta }{3c+1}.$$

Equation (143) shows that, using the constitutive equations (122) and (123), the components of ${\mathit{\tau }}_c^{\mathrm{dev}}$ may be expressed as the functions of ${\lambda }_{e2}$ and ${\lambda }_{e3}$ only since the Lagrangian multiplier $p$ is eliminated automatically from equation (143).

${\mathit{\tau }}_c^{\mathrm{dev}}$ in equation (139) is in a form as the direct driving force of wall expansion. As examples, the expressions of ${\mathit{\tau }}_c^{\mathrm{dev}}$ in two cases are shown as follows:

• Case 1: $\vartheta =0$

$${\tau }_{c11}^{\mathrm{dev}}=\frac{{\tau }_1-{\tau }_3}{2},{\tau }_{c22}^{\mathrm{dev}}=0,{\tau }_{c33}^{\mathrm{dev}}=\frac{{\tau }_3-{\tau }_1}{2}.$$

• Case 2: $\vartheta =\mathrm{\pi }/4$

$${\tau }_{c11}^{\mathrm{dev}}=0,{\tau }_{c22}^{\mathrm{dev}}=\frac{{\tau }_2-{\tau }_3}{2},{\tau }_{c33}^{\mathrm{dev}}=\frac{{\tau }_3-{\tau }_2}{2}.$$

Bear in mind that wall expansion must be associated with volumetric growth to maintain wall thickness and strength. However, ${\mathit{\tau }}_c^{\mathrm{dev}}$ itself is a traceless deviatoric tensor. Thus, although ${\mathit{\tau }}_c^{\mathrm{dev}}$ is assumed as the driving force of growth, a role similar to the deviatoric stress in $J_2$ plasticity, the matching volumetric growth need to be represented by a non-traceless tensor, which is shown in the subsequent discussion about the flow rule of growth.

5.3.2. Specification 2: one family of CMFs and one pure geometrical direction for the remodelling of individual family of CMFs

Here, we consider the remodelling of the $i\mathrm{th}$ family, i.e., one of the two families, of CMFs $(i\in \{1,2\})$ of a representative layer of the cell wall. Thus, different from equation (135), the two unit vectors, ${\tilde{\mathit{a}}}_1$ and ${\tilde{\mathit{a}}}_2$ , are assigned as:

$${\tilde{\mathit{a}}}_1:={\tilde{\mathit{a}}}_i^f=\{0,\cos {\vartheta }_i,\sin {\vartheta }_i{\}}^{\mathrm{T}}\mathrm{and}{\tilde{\mathit{a}}}_2:=\{1,0,0{\}}^{\mathrm{T}}=\mathit{r},$$

(144)

where ${\tilde{\mathit{a}}}_1$ represents the $i\mathrm{th}$ family of CMFs while ${\tilde{\mathit{a}}}_2$ is a pure geometrical direction aligning with the radial direction $\mathit{r}$ to represent the re-orientation direction. Consequently, the structural tensors in equation (76) for defining Spencer’s deviatoric stress are obtained as:

$${\tilde{\mathit{A}}}_1=\frac{1}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+\cos 2{\vartheta }_i\hfill & \sin 2{\vartheta }_i\hfill \\ \multicolumn{1}{c}{0}& \sin 2{\vartheta }_i\hfill & 1-\cos 2{\vartheta }_i\hfill \end{array}\right]\mathrm{and}{\tilde{\mathit{A}}}_2=\left[\begin{array}{ccc}1\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \end{array}\right].$$

(145)

The correspondingly deviatoric tensors are:

$$\mathrm{dev}{\tilde{\mathit{A}}}_1=\frac{1}{6}\left[\begin{array}{ccc}-2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+3\cos 2{\vartheta }_i\hfill & 3\sin 2{\vartheta }_i\hfill \\ \multicolumn{1}{c}{0}& 3\sin 2{\vartheta }_i\hfill & 1-3\cos 2{\vartheta }_i\hfill \end{array}\right]$$

(146)

and

$$\mathrm{dev}{\tilde{\mathit{A}}}_2=\frac{1}{3}\left[\begin{array}{ccc}2\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -1\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -1\hfill \end{array}\right].$$

(147)

Using equations (146) and (147), the three coefficients in equations (78) and (79) are obtained as:

$$\alpha =-\frac{{\tau }_2+{\tau }_3}{2}+\frac{({\tau }_2-{\tau }_3)\cos 2{\vartheta }_i}{2},$$

(148)

$${\beta }_1=-({\tau }_2-{\tau }_3)\cos 2{\vartheta }_i,$$

(149)

$${\beta }_2=-{\tau }_1+\frac{{\tau }_2+{\tau }_3}{2}-\frac{1}{2}({\tau }_2-{\tau }_3)\cos 2{\vartheta }_i.$$

(150)

Substituting equations (145) and (148)–(150) into equation (76) yields the specification of Spencer’s deviatoric stress tensor for the $i\mathrm{th}$ family of CMFs, denoted as:

$${\mathit{\tau }}_i^{\mathrm{dev}}=\frac{{\tau }_3-{\tau }_2}{4}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -(1-\cos 4{\vartheta }_i)\hfill & \sin 4{\vartheta }_i\hfill \\ \multicolumn{1}{c}{0}& \sin 4{\vartheta }_i\hfill & (1-\cos 4{\vartheta }_i)\hfill \end{array}\right].$$

(151)

Equation (151) may be re-expressed in the other form,

$${\mathit{\tau }}_i^{\mathrm{dev}}=\frac{{\tau }_2-{\tau }_3}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& \stackrel{2}{\sin }2{\vartheta }_i\hfill & -\cos 2{\vartheta }_i\sin 2{\vartheta }_i\hfill \\ \multicolumn{1}{c}{0}& -\cos 2{\vartheta }_i\sin 2{\vartheta }_i\hfill & -\stackrel{2}{\sin }2{\vartheta }_i\hfill \end{array}\right],$$

(152)

which indicates that ${\mathit{\tau }}_i^{\mathrm{dev}}$ may be derived from two vectors, $\{0,1,0{\}}^{\mathrm{T}}$ and $\{0,\cos 2{\vartheta }_i,\sin 2{\vartheta }_i{\}}^{\mathrm{T}}$ , as follows:

$${\mathit{\tau }}_i^{\mathrm{dev}}=\frac{{\tau }_2-{\tau }_3}{2}[\left\{\begin{array}{c}0\hfill \\ \multicolumn{1}{c}{1}\\ \multicolumn{1}{c}{0}\end{array}\right\}\left\{\begin{array}{ccc}0\hfill & 1\hfill & 0\hfill \end{array}\right\}-\left\{\begin{array}{c}0\hfill \\ \multicolumn{1}{c}{\cos 2{\vartheta }_i}\\ \multicolumn{1}{c}{\sin 2{\vartheta }_i}\end{array}\right\}\left\{\begin{array}{ccc}0\hfill & \cos 2{\vartheta }_i\hfill & \sin 2{\vartheta }_i\hfill \end{array}\right\}].$$

(153)

${\mathit{\tau }}_i^{\mathrm{dev}}$ is assumed as the driving force of the remodelling of the $i\mathrm{th}$ family of CMFs of the representative layer.

5.4. Yield surface of growth

As aforementioned, cell wall growth is composed of two parts (i) a directly stress-driven irreversible isochoric in-plane expansion and (ii) the out-of-plane volumetric growth by material deposition and transport (see Figure 2(a)). Based on the biological evidences and postulates about wall yielding/loosening [3,13] and mechanics understandings, it is considered that the expansion plays the leading role in wall growth while material deposition and transport are matching the expansion accordingly for maintaining the stability of the wall. Therefore, it is reasonable to assume that a yield surface of wall growth only needs to explicitly account for the yielding condition of the wall expansion. The volumetric growth by material deposition and transport does not affect the yielding condition explicitly but may be taken into account in the flow rule. Under such an assumption, the yield surface of cell wall growth is merely the yielding condition of the stress-driven irreversible isochoric expansion of the wall, which obviously links to von Mises condition in plasticity using a deviatoric stress.

Therefore, it is assumed that the wall growth is driven by the tensor ${\mathit{\tau }}_c^{\mathrm{dev}}$ specified in equation (139) and the growth model of the continuum may be formulated using a set of state variables, $\{{\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1,{\tilde{\mathit{A}}}_2,k\},$ in which ${\tilde{\mathit{A}}}_i:={\tilde{\mathit{A}}}_i^f={\tilde{\mathit{a}}}_i^f\otimes {\tilde{\mathit{a}}}_i^f(i=1,2)$ are specified according to equation (135), and $k$ is a scalar hardening variable. The yield surface, $f_Y$ introduced in equation (96), with isotropic hardening, is proposed as:

$$f_Y\left({\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1^f,{\tilde{\mathit{A}}}_2^f\right)=J_Y\left({\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1^f,{\tilde{\mathit{A}}}_2^f\right)-(1+k)=0,$$

(154)

where the yield function, $J_Y$ on ${\mathcal{B}}_t$ , and its counterpart, ${\overline{J}}_Y$ on $\overline{\mathcal{B}}$ , may be specified, respectively, by adopting the anisotropic yield functions proposed in the previous studies [17,30] as:

$$\begin{array}{c}\multicolumn{1}{c}{J_Y({\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1^f,{\tilde{\mathit{A}}}_2^f)=J_Y^{*}({\tilde{J}}_2,{\tilde{J}}_3,{\tilde{J}}_4,{\tilde{J}}_5)=\sum _{I=2}^5\frac{{\tilde{J}}_I}{Y_I^2}}\\ \multicolumn{1}{c}{\equiv \sum _{I=2}^5\frac{{\overline{J}}_I}{Y_I^2}={\overline{J}}_Y^{*}({\overline{J}}_2,{\overline{J}}_3,{\overline{J}}_4,{\overline{J}}_5)={\overline{J}}_Y({\overline{\mathit{\sum }}}_c^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_1^f,{\widetilde{\overline{\mathit{A}}}}_2^f),}\end{array}$$

(155)

in which, similar to equation (87)₂, a Mandel-stress-style tensor is defined as:

$${\overline{\mathit{\sum }}}_c^{\mathrm{DEV}}{=}^{\mathrm{def}}({\mathit{F}}^{e^{-1}}{\mathit{\tau }}_c^{\mathrm{dev}}{\mathit{F}}^{e^{-\mathrm{T}}})\overline{\mathit{C}}\equiv {\mathit{F}}^{e^{-1}}{\mathit{\tau }}_c^{\mathrm{dev}}{\mathit{F}}^e,$$

(156)

the weighted structural tensor ${\widetilde{\overline{\mathit{A}}}}_i^f(i\in \{1,2\})$ is defined in a way similar to equation (85),

$${\widetilde{\overline{\mathit{A}}}}_i^f={\widetilde{\overline{\mathit{a}}}}_i^f\otimes {\widetilde{\overline{\mathit{a}}}}_i^f,$$

(157)

while ${\widetilde{\overline{\mathit{a}}}}_i^f={\mathit{F}}^{e^{-1}}{\tilde{\mathit{a}}}_i$ , the invariants are referring to the lists in section 4.2 with $\{{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_i,{\overline{\mathit{\sum }}}^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i\}$ replaced by $\{{\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_i^f,{\overline{\mathit{\sum }}}_c^{\mathrm{DEV}},{\widetilde{\overline{\mathit{A}}}}_i^f\},$ respectively, as:

$$\begin{array}{c}{\tilde{J}}_2=\mathrm{tr}\left[{({\mathit{\tau }}_c^{\mathrm{dev}})}^2\right]\equiv \mathrm{tr}\left[{({\overline{\mathit{\sum }}}_c^{\mathrm{DEV}})}^2\right]={\overline{J}}_2;\hfill \\ {\tilde{J}}_3=\mathrm{tr}\left({({\mathit{\tau }}_c^{\mathrm{dev}})}^2{\tilde{\mathit{A}}}_1^f\right)\equiv \mathrm{Tr}\left({({\overline{\mathit{\sum }}}_c^{\mathrm{DEV}})}^2{\widetilde{\overline{\mathit{A}}}}_2^f\right)={\overline{J}}_3;\hfill \\ {\tilde{J}}_4=\mathrm{tr}\left({({\mathit{\tau }}_c^{\mathrm{dev}})}^2{\tilde{\mathit{A}}}_2^f\right)\equiv \mathrm{Tr}\left({({\overline{\mathit{\sum }}}_c^{\mathrm{DEV}})}^2{\widetilde{\overline{\mathit{A}}}}_2^f\right)={\overline{J}}_4;\hfill \\ {\tilde{J}}_5=\mathrm{tr}({\mathit{\tau }}_c^{\mathrm{dev}}·[{\tilde{\mathit{a}}}_1^f\otimes {\tilde{\mathit{a}}}_2^f])\equiv \mathrm{Tr}({\overline{\mathit{\sum }}}_c^{\mathrm{DEV}}·[{\widetilde{\overline{\mathit{a}}}}_1^f\otimes {\widetilde{\overline{\mathit{a}}}}_2^f])={\overline{J}}_5;\hfill \end{array}$$

the scalars, $Y_2$ , $Y_3=Y_4$ , and $Y_5$ , are the yield parameters with the dimension of stress, and the hardening variable $k$ is governed by the evolution law (58)₂ in a general form and may be specified in a form as equation (73).

Noting that yielding and softening of the growing cell wall have their own biological mechanisms. The yielding is due to the proactive loosening of cell wall crosslinks to enable growth while the softening is part of the dynamic regulation of the cell wall to balance the growth and strength [2,3].

5.5. Flow rules of growth on $\overline{\mathcal{B}}$ _t and $\overline{\mathcal{B}}$

Let us consider the growth law in equation (57), ${\mathit{D}}^g=\mathrm{\Phi }{\mathit{N}}^g$ and ${\mathit{W}}^g=\mathrm{\Phi }{\mathit{M}}^g$ , on ${\mathcal{B}}_t$ . In consistency with the yield surface in equation (154), the overstress function $\mathrm{\Phi }$ may be defined as:

$$\mathrm{\Phi }=〈J_Y({\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1^f,{\tilde{\mathit{A}}}_2^f)-(1+k)〉,$$

(158)

where <*> are the Macaulay brackets.

As aforementioned, Spencer’s deviatoric stress tensor ${\mathit{\tau }}_c^{\mathrm{dev}}$ in equation (139) is assumed as the driving force of the wall growth, but its traceless feature indicates that the associated volumetric growth must be represented by other non-traceless tensors. Thus, taking into account the axisymmetric feature, the flow directions of the growth are supposed as:

$${\mathit{N}}^g={\mathit{N}}_{\mathrm{s}}^g+{\mathit{N}}_{\mathrm{d}}^g,{\mathit{M}}^g=0,$$

(159)

where the two parts of ${\mathit{N}}^g$ are:

$${\mathit{N}}_{\mathrm{s}}^g=\frac{{\mathit{\tau }}_c^{\mathrm{dev}}}{\parallel {\mathit{\tau }}_c^{\mathrm{dev}}\parallel }\mathrm{and}{\mathit{N}}_{\mathrm{d}}^g=-(\mathit{r}·{\mathit{N}}_{\mathrm{s}}^g·\mathit{r})[\mathit{r}\otimes \mathit{r}].$$

(160)

Equation (159) ₁ may be considered as a non-associated flow rule on ${\mathcal{B}}_t$ . It may be interpreted that ${\mathit{N}}_{\mathrm{s}}^g$ represents the direction of the stress-driven wall expansion while ${\mathit{N}}_{\mathrm{d}}^g$ represents the direction of the matching volumetric growth resulting from the new material deposition and transport (Figure 2(a)). Noting that $\mathrm{tr}{\mathit{N}}_{\mathrm{s}}^g=0$ and $\mathit{r}·{\mathit{N}}^g·\mathit{r}=0$ , which mean (i) the directly stress-driven wall expansion direction ${\mathit{N}}_{\mathrm{s}}^g$ is isochoric so it tends to reduce the wall thickness (and strength) during expansion and (ii) the matching volumetric growth is supplied by the new material deposition and transport in the wall thickness direction $\mathit{r}$ for keeping a constant wall thickness on $\overline{\mathcal{B}}$ (see Figure 2). This is how the cell wall growth differs from the current engineering AM: deposition and transport of new material to maintain the in-plane wall expansion rather than to create an out-of-plane extrusion.

Moreover, we have the pull-back formulae of (i) rate of deformation, $\overline{{\mathit{D}}^g}={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{D}}^g{\mathit{F}}^e$ , and (ii) spins,

$$\overline{\mathbf{\omega }}={\mathit{F}}^{e^{\mathrm{T}}}\mathbf{\omega }{\mathit{F}}^e\mathrm{and}\overline{{\mathit{W}}^g}={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{W}}^g{\mathit{F}}^e,$$

(161)

which are consistent with equation (45). Thus, it is straightforward to obtain the counterpart of the growth law (57) on $\overline{\mathcal{B}}$ as:

$$\overline{{\mathit{D}}^g}=\mathrm{\Phi }{\overline{\mathit{N}}}^g\mathrm{and}\overline{\mathit{W}}=\mathrm{\Phi }{\overline{\mathit{M}}}^g,$$

(162)

where the growth directions on $\overline{\mathcal{B}}$ are:

$${\overline{\mathit{N}}}^g={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{N}}^g{\mathit{F}}^e\mathrm{and}{\overline{\mathit{M}}}^g={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{M}}^g{\mathit{F}}^e.$$

(163)

As the counterpart of equation (159), the specification of the growth directions on $\overline{\mathcal{B}}$ is calculated from equation (163) as:

$${\overline{\mathit{N}}}^g=\frac{{\mathit{F}}^{e^{\mathrm{T}}}{\mathit{\tau }}_c^{\mathrm{dev}}{\mathit{F}}^e}{\parallel {\mathit{\tau }}_c^{\mathrm{dev}}\parallel }-(\mathit{r}·{\mathit{N}}_{\mathrm{s}}^g·\mathit{r})[({\mathit{F}}^{e^{\mathrm{T}}}\mathit{r})\otimes ({\mathit{F}}^{e^{\mathrm{T}}}\mathit{r})],{\overline{\mathit{M}}}^g=0.$$

(164)

Noting that the r.h.s of equation (164)₁ is an expression in terms of the variables on ${\mathcal{B}}_t$ . To obtain the alternative expression in terms of the variables on $\overline{\mathcal{B}}$ , the two equations in the growth law (57) are summed up as ${\mathit{L}}^g=\mathrm{\Phi }({\mathit{N}}^g+{\mathit{M}}^g)$ , which is then pulled back to $\overline{\mathcal{B}}$ as:

$${\overline{\mathit{L}}}^g=\mathrm{\Phi }{\mathit{F}}^{e^{-1}}({\mathit{N}}^g+{\mathit{M}}^g){\mathit{F}}^e,$$

(165)

using the relation, ${\overline{\mathit{L}}}^g={\mathit{F}}^{e^{-1}}{\mathit{L}}^g{\mathit{F}}^e$ . Substituting the specification (159) into equation (165) yields:

$${\overline{\mathit{L}}}^g=\mathrm{\Phi }\{\frac{{\overline{\mathit{\sum }}}_c^{\mathrm{DEV}}}{\parallel {\overline{\mathit{\sum }}}_c^{\mathrm{DEV}}\parallel }-(\mathit{r}·{\mathit{N}}_{\mathrm{s}}^g·\mathit{r})[({\mathit{F}}^{e^{-1}}\mathit{r})\otimes ({\mathit{F}}^{e^{\mathrm{T}}}\mathit{r})]\},$$

(166)

where we use the definition of the Mandel-stress-style tensor ${\overline{\mathit{\sum }}}_c^{\mathrm{DEV}}$ in equation (156). From the point of view of the conventional continuum, equation (166) is a non-associated flow rule on $\overline{\mathcal{B}}$ .

5.6. Flow rules of remodelling on $\overline{\mathcal{B}}$ _t and $\overline{\mathcal{B}}$

Hitherto, we complete the growth model of the continuum with isotropic hardening using the set of specified state variables, $\{{\mathit{\tau }}_c^{\mathrm{dev}},{\tilde{\mathit{A}}}_1:={\tilde{\mathit{A}}}_1^f,{\tilde{\mathit{A}}}_2:={\tilde{\mathit{A}}}_2^f,k\}$ accounting for the two families of CMFs in one layer. For modelling the remodelling of individual family of CMFs in the same layer with kinematic hardening, we use a set of specified state variables of the $i\mathrm{th}$ family of CMFs, $\{{\mathit{\tau }}_i^{\mathrm{dev}},{\tilde{\mathit{A}}}_1:={\tilde{\mathit{A}}}_i^f,{\tilde{\mathit{A}}}_2:=\mathit{r}\otimes \mathit{r},{\mathit{\alpha }}_i\},$ where ${\mathit{\tau }}_i^{\mathrm{dev}}$ is defined in equation (151), ${\tilde{\mathit{A}}}_i={\tilde{\mathit{a}}}_i\otimes {\tilde{\mathit{a}}}_i(i=1,2)$ are specified by equation (144), and ${\mathit{\alpha }}_i$ is a second-order hardening tensor associated with the $i\mathrm{th}$ family of CMFs which is defined later. Thus, we shall be able to formulate the remodelling of each and every family of CMFs using the invariants and spin representation presented in section 4.2 and 4.3 with $\{{\mathit{\tau }}^{\mathrm{dev}},{\tilde{\mathit{A}}}_1,{\tilde{\mathit{A}}}_2,\mathit{\alpha }\}$ replaced by $\{{\mathit{\tau }}_i^{\mathrm{dev}},{\tilde{\mathit{A}}}_i^f,\mathit{r}\otimes \mathit{r},{\mathit{\alpha }}_i\}.$

We consider the remodelling law in equation (59)₁, ${\mathit{W}}_i^g={\mathrm{\Phi }}^r{\mathit{M}}_i^g$ , of the $i\mathrm{th}$ family of CMFs. The overstress function ${\mathrm{\Phi }}^r$ is discussed later in section 5.7. Based on the discussion of the representation of spins in section 4.3, here we show two different formulations to construct the flow direction ${\mathit{M}}_i^g$ of remodelling.

5.6.1. Flow direction of remodelling on $\overline{\mathcal{B}}$ _t : formulation 1

In this formulation, Spencer’s deviatoric stress tensor, ${\mathit{\tau }}_i^{\mathrm{dev}}$ , and the associated geometrical structural tensor, ${\tilde{\mathit{A}}}_i^f={\tilde{\mathit{a}}}_i^f\otimes {\tilde{\mathit{a}}}_i^f$ , of the $i\mathrm{th}$ family of CMFs, which is assigned to ${\tilde{\mathit{A}}}_1$ in equation (145)₁ as:

$${\tilde{\mathit{A}}}_1:={\tilde{\mathit{A}}}_i^f=\frac{1}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+\cos 2{\vartheta }_i\hfill & \sin 2{\vartheta }_i\hfill \\ \multicolumn{1}{c}{0}& \sin 2{\vartheta }_i\hfill & 1-\cos 2{\vartheta }_i\hfill \end{array}\right],$$

(167)

are used to construct ${\mathit{M}}_i^g$ in a form referring to equation (94) as:

$${\mathit{M}}_i^g={\mathit{\tau }}_i^{\mathrm{dev}}{\tilde{\mathit{A}}}_i^f-{\tilde{\mathit{A}}}_i^f{\mathit{\tau }}_i^{\mathrm{dev}}.$$

(168)

Therefore, substituting equations (151) and (167) into equation (168) yields:

$${\mathit{M}}_i^g=\frac{{\tau }_3-{\tau }_2}{2}\sin 2{\vartheta }_i\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -1\hfill \\ \multicolumn{1}{c}{0}& 1\hfill & 0\hfill \end{array}\right]\forall i\in \{1,2\}.$$

(169)

5.6.2. Flow direction of remodelling on $\overline{\mathcal{B}}$ _t : formulation 2

Different from the formulation (168), in this formulation, Spencer’s deviatoric stress tensor, ${\mathit{\tau }}_i^{\mathrm{dev}}$ , and the stress-like hardening tensor, ${\mathit{\alpha }}_i$ , introduced in equation (59)₂ for the $i\mathrm{th}$ family of CMFs, are used to construct ${\mathit{M}}_i^g$ referring to equation (95) as:

$${\mathit{M}}_i^g=\frac{{\mathit{\tau }}_i^{\mathrm{dev}}{\mathit{\alpha }}_i-{\mathit{\alpha }}_i{\mathit{\tau }}_i^{\mathrm{dev}}}{\parallel {\mathit{\tau }}_i^{\mathrm{dev}}-{\mathit{\alpha }}_i\parallel },$$

(170)

where, being consistent with equation (89) and similar to ${\mathit{\tau }}_i^{\mathrm{dev}}$ , the hardening tensor ${\mathit{\alpha }}_i$ satisfies the constraints,

$$\mathrm{tr}{\mathit{\alpha }}_i=0\mathrm{and}{\mathit{\alpha }}_i:{\hat{\mathit{A}}}_i^h=0\forall i\in \{1,2\},$$

(171)

in which ${\hat{\mathit{A}}}_i^h(\mathrm{tr}{\hat{\mathit{A}}}_i^h=1\forall i\in \{1,2\})$ are the structural tensors characterising the hardening tensor ${\mathit{\alpha }}_i$ . Noting that ${\mathit{\tau }}_i^{\mathrm{dev}}$ and ${\mathit{\alpha }}_i$ are characterised by the structural tensors ${\tilde{\mathit{A}}}_i^f$ and ${\hat{\mathit{A}}}_i^h$ , respectively. It may be interpreted that ${\hat{\mathit{A}}}_i^h$ is the background geometrical structural tensor associated with the hardening state of the matrix resulting from the remodelling of the $i\mathrm{th}$ family of CMFs characterised by the state of ${\tilde{\mathit{A}}}_i^f$ . Thus, ${\hat{\mathit{A}}}_i^h$ need not to be identical to ${\tilde{\mathit{A}}}_i^f$ .

For clarity, let us construct ${\mathit{\alpha }}_i$ explicitly in a way analogous to that of ${\mathit{\tau }}_i^{\mathrm{dev}}$ . Similar to equation (145)₁, ${\hat{\mathit{A}}}_i^h$ may be expressed as:

$${\hat{\mathit{A}}}_i^h=\frac{1}{2}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 1+\cos 2{\hat{\vartheta }}_i\hfill & \sin 2{\hat{\vartheta }}_i\hfill \\ \multicolumn{1}{c}{0}& \sin 2{\hat{\vartheta }}_i\hfill & 1-\cos 2{\hat{\vartheta }}_i\hfill \end{array}\right],$$

(172)

where ${\hat{\vartheta }}_i$ , analogous to ${\vartheta }_i$ , is an in-plane angle on the $\mathit{t}\wedge {\mathit{e}}_3$ plane. Then, analogous to equation (145), two structural tensors, ${\tilde{\mathit{A}}}_j^{\alpha }$ $(j\in \{1,2\})$ in which the superscript ‘ $\alpha$ ’ stands for the hardening variable, are proposed as:

$${\tilde{\mathit{A}}}_1^{\alpha }={\hat{\mathit{A}}}_i^h\mathrm{and}{\tilde{\mathit{A}}}_2^{\alpha }=\left[\begin{array}{ccc}1\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill \end{array}\right].$$

(173)

By introduction of a so-called back stress, ${\mathit{\tau }}^b=\mathrm{Diag}[{\tau }_1^b,{\tau }_2^b,{\tau }_3^b]$ , the tensor ${\mathit{\alpha }}_i$ may be expressed in a form analogous to ${\mathit{\tau }}^{\mathrm{dev}}$ in equation (76):

$${\mathit{\alpha }}_i={\mathit{\tau }}^b+{\alpha }^b\mathit{I}+\sum _{j=1}^2\left({\beta }_j^b{\tilde{\mathit{A}}}_j^{\alpha }\right),$$

(174)

where ${\alpha }^b$ , ${\beta }_1^b$ , and ${\beta }_2^b$ may be calculated using equations (78) and (79) with $\mathit{\tau }$ and ${\tilde{\mathit{A}}}_j$ replaced by ${\mathit{\tau }}^b$ and ${\tilde{\mathit{A}}}_j^{\alpha }$ , respectively. Accordingly, by substituting equation (173) into equation (174), ${\mathit{\alpha }}_i$ may be obtained in a form analogous to equation (152):

$${\mathit{\alpha }}_i={\hat{\alpha }}_{23}^i\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -\sin 2{\hat{\vartheta }}_i\hfill & \cos 2{\hat{\vartheta }}_i\hfill \\ \multicolumn{1}{c}{0}& \cos 2{\hat{\vartheta }}_i\hfill & \sin 2{\hat{\vartheta }}_i\hfill \end{array}\right],$$

(175)

in which:

$${\hat{\alpha }}_{23}^i=\frac{{\tau }_3^b-{\tau }_2^b}{2}\sin 2{\hat{\vartheta }}_i.$$

(176)

Therefore, substituting equations (152) and (175) into equation (170) yields the flow direction of remodelling,

$${\mathit{M}}_i^g={\hat{\alpha }}_{23}^i\sin 2{\vartheta }_i\sin (2{\vartheta }_i-2{\hat{\vartheta }}_i)\frac{{\tau }_3-{\tau }_2}{\parallel {\mathit{\tau }}_i^{\mathrm{dev}}-{\mathit{\alpha }}_i\parallel }\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -1\hfill \\ \multicolumn{1}{c}{0}& 1\hfill & 0\hfill \end{array}\right]\forall i\in \{1,2\}.$$

(177)

The comparison between equations (169) and (177) indicates that the flow directions of remodelling may be expressed in a united form as:

$${\mathit{M}}_i^g=f_i^r\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & -1\hfill \\ \multicolumn{1}{c}{0}& 1\hfill & 0\hfill \end{array}\right]\forall i\in \{1,2\},$$

(178)

where the function $f_i^r$ has the expression,

$$f_i^r=\frac{{\tau }_3-{\tau }_2}{2}\sin 2{\vartheta }_i,$$

(179)

for formulation 1, or,

$$f_i^r={\hat{\alpha }}_{23}^i\sin 2{\vartheta }_i\sin (2{\vartheta }_i-2{\hat{\vartheta }}_i)\frac{{\tau }_3-{\tau }_2}{\parallel {\mathit{\tau }}_i^{\mathrm{dev}}-{\mathit{\alpha }}_i\parallel },$$

(180)

for formulation 2. This completes the formulation of the flow direction of remodelling on ${\mathcal{B}}_t$ .

5.6.3. Flow direction of remodelling on $\overline{\mathcal{B}}$ : formulations 1 and 2

Regarding the yield function ${\overline{J}}_Y$ in equation (155) and the growth flow direction ${\overline{\mathit{L}}}^g$ in equation (166), it is interesting to formulate the counterparts of equations (168) and (170) on $\overline{\mathcal{B}}$ using Mandel-stress-style tensors. This may be obtained using a pull-back transformation similar to equation (163)₂,

$${\overline{\mathit{M}}}_i^g={\mathit{F}}^{e^{\mathrm{T}}}{\mathit{M}}_i^g{\mathit{F}}^e,$$

(181)

where ${\overline{\mathit{M}}}_i^g$ is the flow direction of remodelling on $\overline{\mathcal{B}}$ .

First, a Mandel-stress-style tensor corresponding to ${\mathit{\tau }}_i^{\mathrm{dev}}$ is defined as:

$${\overline{\mathit{\sum }}}_i^{\mathrm{DEV}}{=}^{\mathrm{def}}({\mathit{F}}^{e^{-1}}{\mathit{\tau }}_i^{\mathrm{dev}}{\mathit{F}}^{e^{-\mathrm{T}}})\overline{\mathit{C}}\equiv {\mathit{F}}^{e^{-1}}{\mathit{\tau }}_i^{\mathrm{dev}}{\mathit{F}}^e.$$

(182)

Second, similar to equations (90) and (91), two hardening tensors of remodelling are introduced as:

$${\overline{\mathit{\varkappa }}}_i{=}^{\mathrm{def}}{\overline{\mathit{\alpha }}}_i\overline{\mathit{C}}\equiv {\mathit{F}}^{e^{-1}}{\mathit{\alpha }}_i{\mathit{F}}^e$$

(183)

and

$${\overline{\mathit{\alpha }}}_i={\mathit{F}}^{e^{-1}}{\mathit{\alpha }}_i{\mathit{F}}^{e^{-\mathrm{T}}}.$$

(184)

Again, it may be interpreted that ${\overline{\mathit{\alpha }}}_i$ is a hardening tensor analogous to the PK2 stress tensor while ${\overline{\mathit{\varkappa }}}_i$ is a Mandel-stress-like hardening tensor on $\overline{\mathcal{B}}$ .

Consequently, the formula (170) of the formulation 2 on ${\mathcal{B}}_t$ may be pulled back to $\overline{\mathcal{B}}$ using equations (181)–(183) as:

$${\overline{\mathit{M}}}_i^g=\frac{-2\mathrm{skew}\left({\overline{\mathit{\sum }}}_i^{\mathrm{DEV}}{\overline{\mathit{\varkappa }}}_i\right)}{\parallel {\overline{\mathit{\sum }}}_i^{\mathrm{DEV}}-{\overline{\mathit{\varkappa }}}_i\parallel }\forall i\in \{1,2\},$$

(185)

where the operator, $\mathrm{skew}(·)=\frac{1}{2}(\overline{\mathit{C}}(·)-(·{)}^{\mathrm{T}}\overline{\mathit{C}})$ , is defined in the same format as in equation (17). As expected, the flow direction ${\overline{\mathit{M}}}_i^g$ is expressed in terms of the Mandel-stress-style tensors, ${\overline{\mathit{\sum }}}_i^{\mathrm{DEV}}$ and ${\overline{\mathit{\varkappa }}}_i$ .

In a similar way, the pull-back of equation (168) of the formulation 1 yields:

$${\overline{\mathit{M}}}_i^g=-2\mathrm{skew}({\overline{\mathit{\sum }}}_i^{\mathrm{DEV}}({\mathit{F}}^{e^{-1}}{\tilde{\mathit{A}}}_i^f{\mathit{F}}^{e^{-\mathrm{T}}})\overline{\mathit{C}})\forall i\in \{1,2\}.$$

(186)

5.7. Hardening law of remodelling

Different from the isotropic hardening law of growth in equation (158), a kinematic hardening law is assumed for the remodelling. Hence, the evolution of the hardening tensor ${\mathit{\alpha }}_i$ , which is proposed in a general form of equation (59)₂, may be specified in a form consistent with equation (75) as,

$${\stackrel{{\nabla }^s}{\mathit{\alpha }}}_i=H_{i1}{\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i)\frac{{\mathit{n}}_{\alpha }^i}{\parallel {\mathit{n}}_{\alpha }^i\parallel }+H_{i2}{\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i){\mathit{\alpha }}_i+H_{i3}\parallel {\mathit{\alpha }}_i{\parallel }^{m_{i3}}{\mathit{\alpha }}_i,$$

(187)

where $\stackrel{{\nabla }^u}{(·)}$ in equation (75) is replaced by $\stackrel{{\nabla }^s}{(·)}$ according to the features of ${\mathit{\alpha }}_i$ shown in equation (171), so that the objective rate is structure-preserving.

The driving force ${\mathit{n}}_{\alpha }^i$ in equation (187) is specified in a form similar to that of the $J_2$ theory [44,58],

$${\mathit{n}}_{\alpha }^i={\mathit{\tau }}_i^{\mathrm{dev}}-{\mathit{\alpha }}_i.$$

(188)

By substituting equations (152) and (175) into equation (188), the components of ${\mathit{n}}_{\alpha }^i$ are obtained as:

$${\mathit{n}}_{\alpha }^i=\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& n_{\alpha 22}^i\hfill & n_{\alpha 23}^i\hfill \\ \multicolumn{1}{c}{0}& n_{\alpha 23}^i\hfill & -n_{\alpha 22}^i\hfill \end{array}\right],$$

(189)

where:

$$n_{\alpha 22}^i=-\frac{{\tau }_3-{\tau }_2}{2}\stackrel{2}{\sin }2{\vartheta }_i+{\hat{\alpha }}_{23}^i\sin 2{\hat{\vartheta }}_i,n_{\alpha 23}^i=\frac{{\tau }_3-{\tau }_2}{2}\sin 2{\vartheta }_i\cos 2{\vartheta }_i-{\hat{\alpha }}_{23}^i\cos 2{\hat{\vartheta }}_i.$$

(190)

Accordingly, the norm of ${\mathit{n}}_{\alpha }^i$ is obtained as:

$$\parallel {\mathit{n}}_{\alpha }^i\parallel =\sqrt{2}\sqrt{{\left(\frac{{\tau }_3-{\tau }_2}{2}\sin 2\vartheta \right)}^2+{\hat{\alpha }}_{23}^2-({\tau }_3-{\tau }_2){\hat{\alpha }}_{23}\sin 2\vartheta \cos 2(\vartheta -\hat{\vartheta })}.$$

(191)

With ${\mathit{n}}_{\alpha }^i$ obtained, the overstress function ${\mathrm{\Phi }}^r({\mathit{n}}_{\alpha }^i)$ in the hardening law (187) (and in the flow rule (59)₁) may be defined as:

$${\mathrm{\Phi }}^r=〈\parallel {\mathit{n}}_{\alpha }^i\parallel -Y_r〉,$$

(192)

where $Y_r$ is a remodelling parameter. Equation (192) implies a remodelling of CMFs not explicitly depending on the wall expansion. Alternatively, equation (192) may be modified as ${\mathrm{\Phi }}^r=〈\parallel {\mathit{n}}_{\alpha }^i\parallel -Y_r〉\mathrm{\Phi }$ by multiplying the overstress function of the wall expansion, $\mathrm{\Phi }$ , so that the remodelling is explicitly associated with the wall expansion.

Hitherto, the objective rate of ${\mathit{\alpha }}_i$ is determined from the hardening law (187). To calculate the material time derivative, ${\stackrel{·}{\mathit{\alpha }}}_i$ , for a complete update of ${\mathit{\alpha }}_i$ at time $t$ , we refer to equation (214) in Appendix 1 (The objective rate of Spencer’s deviatoric stress tensor) and propose that the objective rate is related with the material time derivative in such a form,

$${\stackrel{{\nabla }^s}{\mathit{\alpha }}}_i={\stackrel{·}{\mathit{\alpha }}}_i-\delta {\mathit{\alpha }}_i,$$

(193)

where $\delta {\mathit{\alpha }}_i$ may be expressed in a format analogous to either equation (223) or equation (215). As mentioned in Appendix 1 (The objective rate of Spencer’s deviatoric stress tensor), for analytical analysis, it is more convenient to use equation (215). Therefore, referring to the definition of ${\mathit{\alpha }}_i$ in equation (174), we write the general form of $\delta {\mathit{\alpha }}_i$ as follows:

$$\delta {\mathit{\alpha }}_i=[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}+\sum _{j=1}^2(\frac{\partial {\alpha }^b}{\partial {\tilde{\mathit{A}}}_j^{\alpha }}:\delta {\tilde{\mathit{A}}}_j^{\alpha })\mathit{I}+\sum _{l,j=1}^2\left((\frac{\partial {\beta }_l^b}{\partial {\tilde{\mathit{A}}}_j^{\alpha }}:\delta {\tilde{\mathit{A}}}_j^{\alpha }){\tilde{\mathit{A}}}_l^{\alpha }\right)+\sum _{j=1}^2({\beta }_j^b\delta {\tilde{\mathit{A}}}_j^{\alpha }),$$

(194)

where, analogous to $[\delta \mathit{\tau }{]}^{\mathrm{dev}}$ in equation (217), $[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}$ may be calculated by:

$$[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}=\delta {\mathit{\tau }}^b+(\frac{\partial {\alpha }^b}{\partial {\mathit{\tau }}^b}:\delta {\mathit{\tau }}^b)\mathit{I}+\sum _{j=1}^2\left((\frac{\partial {\beta }_j^b}{\partial {\mathit{\tau }}^b}:\delta {\mathit{\tau }}^b){\tilde{\mathit{A}}}_j^{\alpha }\right),$$

(195)

in which, analogous to equation (216), the increment of ${\mathit{\tau }}^b$ is:

$$\delta {\mathit{\tau }}^b={\mathit{L}}^e{\mathit{\tau }}^b+{\mathit{\tau }}^b{\mathit{L}}^{e^{\mathrm{T}}}.$$

(196)

According to the analogies between $\{{\mathit{\tau }}^b,{\hat{\vartheta }}_i\}$ and $\{\mathit{\tau },{\vartheta }_i\},$ the specification of the components of $[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}$ in the r.h.s of equation (194) is obtained by a direct analogy with equation (151) as:

$$[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}=\frac{\delta {\tau }_3^b-\delta {\tau }_2^b}{4}\left[\begin{array}{ccc}0\hfill & 0\hfill & 0\hfill \\ \multicolumn{1}{c}{0}& -(1-\cos 4{\hat{\vartheta }}_i)\hfill & \sin 4{\hat{\vartheta }}_i\hfill \\ \multicolumn{1}{c}{0}& \sin 4{\hat{\vartheta }}_i\hfill & (1-\cos 4{\hat{\vartheta }}_i)\hfill \end{array}\right].$$

(197)

For the remaining part in the r.h.s of equation (194), we first calculate the increment, $\delta {\tilde{\mathit{A}}}_j^{\alpha }$ , in the last term of r.h.s in a way analogous to equation (218) as:

$$\delta {\tilde{\mathit{A}}}_j^{\alpha }=\left[L^e-\mathrm{tr}\left({\mathit{L}}^e{\tilde{\mathit{A}}}_j^{\alpha }\right)\mathit{I}\right]{\tilde{\mathit{A}}}_j^{\alpha }+{\tilde{\mathit{A}}}_j^{\alpha }{\left[{\mathit{L}}^e-\mathrm{tr}\left({\mathit{L}}^e{\tilde{\mathit{A}}}_j^{\alpha }\right)\mathit{I}\right]}^{\mathrm{T}}.$$

(198)

Second, the coefficients, ${\beta }_j^b$ $(j\in \{1,2\})$ , also in the last term of r.h.s are calculated as the analogies to equations (149) and (150),

$${\beta }_1^b=-({\tau }_2^b-{\tau }_3^b)\cos 2{\hat{\vartheta }}_i,$$

(199)

$${\beta }_2^b=-{\tau }_1^b+\frac{{\tau }_2^b+{\tau }_3^b}{2}-\frac{1}{2}({\tau }_2^b-{\tau }_3^b)\cos 2{\hat{\vartheta }}_i.$$

(200)

Third, we calculate the second term in the r.h.s of equation (194),

$$\sum _{j=1}^2(\frac{\partial {\alpha }^b}{\partial {\tilde{\mathit{A}}}_j^{\alpha }}:\delta {\tilde{\mathit{A}}}_j^{\alpha })={\mathit{\tau }}^b:\delta {\tilde{\mathit{A}}}_1^{\alpha },$$

(201)

and part of the expression inside the third term,

$$\sum _{j=1}^2(\frac{\partial {\beta }_i^b}{\partial {\tilde{\mathit{A}}}_j^{\alpha }}:\delta {\tilde{\mathit{A}}}_j^{\alpha })=2{\mathit{\tau }}^b:[-\frac{3}{2}\mathrm{dev}\delta {\tilde{\mathit{A}}}_{i^{*}}^{\alpha }-(\delta {\tilde{\mathit{A}}}_i^{\alpha }-\delta {\tilde{\mathit{A}}}_{i^{*}}^{\alpha })],$$

(202)

which gives the whole third term in the r.h.s of equation (194),

$$\sum _{l,j=1}^2\left((\frac{\partial {\beta }_i^b}{\partial {\tilde{\mathit{A}}}_j^{\alpha }}:\delta {\tilde{\mathit{A}}}_j^{\alpha }){\tilde{\mathit{A}}}_l^{\alpha }\right)=-(2{\mathit{\tau }}^b:\delta {\tilde{\mathit{A}}}_1^{\alpha }){\tilde{\mathit{A}}}_1^{\alpha }-({\mathit{\tau }}^b:\delta {\tilde{\mathit{A}}}_1^{\alpha }){\tilde{\mathit{A}}}_2^{\alpha }.$$

(203)

Substituting equations (201) and (203) into equation (194) yields:

$$\delta {\mathit{\alpha }}_i=[\delta {\mathit{\tau }}^b{]}^{\mathrm{dev}}+({\mathit{\tau }}^b:\delta {\tilde{\mathit{A}}}_1^{\alpha })\{\mathit{I}-2{\tilde{\mathit{A}}}_1^{\alpha }-{\tilde{\mathit{A}}}_2^{\alpha }\}+{\beta }_1^b\delta {\tilde{\mathit{A}}}_1^{\alpha }.$$

(204)

This completes the calculation of $\delta {\mathit{\alpha }}_i$ in equation (194). Similar to equation (222), it can be proven that $\mathrm{tr}(\delta {\mathit{\alpha }}_i)=0$ .

With ${\stackrel{{\nabla }^s}{\mathit{\alpha }}}_i$ and $\delta {\mathit{\alpha }}_i$ obtained in equations (187) and (204), respectively, equation (193) may be re-arranged as:

$${\stackrel{·}{\mathit{\alpha }}}_i={\stackrel{{\nabla }^s}{\mathit{\alpha }}}_i+\delta {\mathit{\alpha }}_i.$$

(205)

Substituting equations (187) and (204) into the r.h.s of equation (205) yields the temporal evolution of ${\mathit{\alpha }}_i$ $(i\in \{1,2\})$ .

Noting that equation (205) includes a set of two independent scalar equations sufficient to compute the evolution of the two independent scalar variables, ${\hat{\alpha }}_{23}^i$ (or $({\tau }_3^b-{\tau }_2^b)$ ) and ${\hat{\vartheta }}_i$ , which gives the evolution of the tensor ${\mathit{\alpha }}_i$ as shown in equation (175).

6. Conclusion

The recent advances in cell wall bring the researchers closer to meeting two grand challenges in this field: (i) to connect cell wall structure to mechanics at the molecular and cellular scales and (ii) to connect cell wall structure to cell growth by identifying the molecular sites of wall loosening and the molecular movements of wall polymers as the cell wall surface enlarges [5]. In order to reflect those advances on continuum theory and to address the challenges from the point of view of continuum mechanics, this study proposes the concept of growth spin and the corresponding model of growth and remodelling of cell wall.

First, the model relaxes the perfectly bonding assumption, which has been widely used in both elastic [33] and inelastic [30] fibre-reinforced material models including cell wall growth models [17,73] with inspiration from Spencer’s work [30,64], so that the continuum model may accommodate the new observation and concepts about cell wall structure at nano-/micro-scale, especially the emerging concept that CMFs make a load-bearing network via close physical contacts rather than a tethered network of CMFs.

Second, it is demonstrated that in the proposed model, the three major components—i.e., CMFs, hemicellulose, and pectins—may have their own growth/remodelling spins, just similar to the plastic spins [38,48,56,69,74] of the slip/twining systems in crystal plasticity, which renders a platform of formulation not only for the antisymmetric part of growth deformation gradient (a term missing from most of the existing kinematical growth models) but also for the remodelling of the individual family of CMFs. This feature provides the more sophisticated and flexible kinematical description of cell wall aiming at quantitatively modelling the postulates of (i) the selective control of CMFs and (ii) flexible means to modulate the directionality of cell growth [5,75].

Third, the corresponding isotropic and kinematic hardening laws are proposed for growth and remodelling, respectively, which further shows the potential of the present model to capture the selective control and flexible regulation of the directionality of cell growth. In addition to generalising the strain-induced scalar hardening law reported in Huang et al. [17] into a tensorial format, a hardening term not explicitly depending on strain is introduced to account for chemical regulation. Together, the proposed hardening model provides a new modelling framework for addressing one of the key questions of cell wall growth: the balance between hardening and softening of the wall material for a dynamical control of growth while maintaining the material strength.

Fourth, this study thoroughly demonstrates the application of Spencer’s deviatoric stress tensors for the modelling of cell wall growth and remodelling. In addition to showing the role of Spencer’s deviatoric stress tensors in the yield surfaces and flow rules, the kinematic hardening law of remodelling shows the calculation and application of the objective rate of Spencer’s deviatoric stress tensor.

Fifth, to the best of the author’s knowledge, unlike the macroscopic plasticity theory which has the micro-structural interpretation using crystal plasticity theory, the kinematical growth of soft matter does not have such a micro-structural counterpart. This caused the difficulty in the interpretation of underlining mechanisms and may be one of the reasons that the growth configuration $\overline{\mathcal{B}}$ was considered as a “fictitious” configuration [41]. Based on the recent advances on the nano-/micro-scale mechanics of cell wall, this study may shed a new light into the micro-structural interpretation of the kinematical growth.

In summary, the growth spin model suggests a new framework of constitutive modelling for reflecting the recent advances in the research of cell wall growth and addressing some long-term challenges, especially the integration of the components and the mechanisms to maintain the integrity during growth, from the point of view of continuum theory. The formulation in this study is conducted on both the current configuration and intermediate (growth) configuration, which could be applied to finite element formation in future work.