%%   As
%% discussed in Section \ref{sect1}, the solid tissue phase has only
%% a mass source/sink associated with it, and no flux. The fluid
%% tissue phase has only a mass flux, and no
%% source/sink\footnote{Considering the case of the lymphatic fluid,
%% this implies that lymph glands are assumed not to be present.}.

%% Since the solid phase, $\mathrm{s}$, does not undergo transport,
%% its motion is specified entirely by $\Bvarphi(\bX,t)$. We describe
%% the remaining species $\mathrm{f},\alpha,\dots,\omega$ as
%% convecting with the solid phase and diffusing with respect to it.
%% They therefore have a velocity relative to $\mathrm{s}$. Since the
%% remaining species convect with $\mathrm{s}$, it implies a local
%% homogenization of deformation. The modelling assumption is made
%% that at each point, $\bX$, the individual phases undergo the same
%% deformation.

%% The concentrations,
%% $\rho_0^\iota$, change as a result of mass transport and
%% inter-conversion of species, implying that the total density in
%% the reference configuration, $\rho_0$, changes with time. They are
%% parameterized as $\rho_0^\iota(\bX,t)$.

%% We define concentrations of the species
%% $\rho_0^\iota=\bar{\rho}_0^\iota f^\iota$ as masses per unit
%% volume in $\Omega_0$. The intrinsic species density is
%% $\bar{\rho}_0^\iota$, and $f^\iota$ is the volume fraction of
%% $\iota$, for $\iota = \mathrm{s,f},\alpha,\dots,\omega$. In an
%% experiment it is far easier to measure the concentration,
%% $\rho_0^\iota$, rather than the intrinsic species density,
%% $\bar{\rho}_0^\iota$\footnote{In our experiments we have measured
%% the mass concentration $\rho_0^\iota$ of collagen in engineered
%% tendons grown \emph{in vitro}. These results will be presented
%% elsewhere \citep{Calve:04}.}. The concentrations also have
%% the property $\sum\limits_{\iota}\rho_0^\iota = \rho_0$, the total
%% material density of the tissue, with the sum being over all
%% species $\mathrm{s,f},\alpha,\dots,\omega$.

%% \noindent Since the solid phase of the tissue does not undergo
%% mass transport, there is no associated flux. Localizing the result
%% gives

%% \begin{equation}
%% \frac{\partial\rho_0^\mathrm{s}}{\partial t} = \Pi^\mathrm{s},
%% \label{massballocA}
%% \end{equation}
%% \noindent where the explicit dependence upon position and time has
%% been suppressed.

%% The fluid phase of the tissue, $\mathrm{f}$, may be thought of as
%% the interstitial or lymphatic fluid that perfuses the tissue. As
%% explained above, we do not consider sources of fluid in the region
%% of interest. The fluid therefore enters and leaves $\Omega_0$ as a
%% flux, $\bM^\mathrm{f}$. The balance of mass in integral form is
%% \begin{equation}
%% \frac{\mathrm{d}}{\mathrm{d}t} \int\limits_{\Omega_0}
%% \rho_0^\mathrm{f} (\bX,t)\mathrm{d}V =
%% -\int_{\partial\Omega_0}\bM^\mathrm{f}(\bX,t)\cdot\bN \mathrm{d}A,
%% \label{massbalintB}
%% \end{equation}

%% Applying the Divergence Theorem to the surface
%% integral and localizing the result gives
%% \begin{equation}
%% \frac{\partial\rho_0^\mathrm{f}}{\partial t} = -
%% \Bnabla\cdot\bM^\mathrm{f}, \label{massballocB}
%% \end{equation}

%% \noindent where $\Bnabla(\bullet)$ is the gradient operator
%% defined on $\Omega_0$, and $\Bnabla\cdot(\bullet)$ denotes the
%% divergence of a vector or tensor argument on $\Omega_0$.

%% For the precursor and byproduct species,
%% $\iota=\alpha,\dots,\omega$, the balance of mass in integral form
%% is
%% \begin{equation}
%% \frac{\mathrm{d}}{\mathrm{d}t} \int\limits_{\Omega_0} \rho_0^\iota
%% (\bX,t)\mathrm{d}V = \int\limits_{\Omega_0} \Pi^\iota
%% (\bX,t)\mathrm{d}V
%% -\int_{\partial\Omega_0}\bM^\iota(\bX,t)\cdot\bN \mathrm{d}A.
%% \label{massbalintI}
%% \end{equation}

%% \noindent In local form it is
%% \begin{equation}
%% \frac{\partial\rho_0^\iota}{\partial t} = \Pi^\iota -
%% \Bnabla\cdot\bM^\iota,\;\forall\,\iota=\alpha,\dots,\omega.
%% \label{massballocI}
%% \end{equation}

%% \noindent Of course, this last equation is the general form of
%% mass balance for any species $\iota$, recalling that in
%% particular, $\bM^\mathrm{s} = \bzero$ and $\Pi^\mathrm{f} = 0$.
%% This form will be used in the development that follows.

%% The fluxes,
%% $\bM^\iota,\;\forall\,\iota=\mathrm{f},\alpha,\dots,\omega$,
%% represent mass transport of the fluid, of precursors to the
%% reaction site, and of byproducts from sites of tissue breakdown.
%% The sources,
%% $\Pi^\iota,\;\forall\iota=\mathrm{s},\alpha,\dots,\omega$, arise
%% from inter-conversion of species. The sources/sinks in
%% (\ref{massballocA}) and (\ref{massballocI}) are therefore related,
%% as tissue and byproducts are formed by consuming precursors (amino
%% acids and nutrients, for instance). To maintain a degree of
%% simplicity in this initial exposition, we will restrict our
%% description of tissue breakdown to the reverse of this reaction.

%% %% The behaviour of the entire system can be determined by summing
%% %% \mbox{Equation~(\ref{massbalance1})} over all species $\iota$.
%% %% Additionally, sources and sinks satisfy the relation

%% %% \begin{equation}
%% %% \sum\limits_\iota\Pi^\iota = 0, \label{sourcebalance}
%% %% \end{equation}

%% %% \noindent which is consistent \citep{growthpaper} with the Law of Mass
%% %% Action for reaction rates and with Mixture Theory
%% %% \citep{TruesdellNoll:65}.


%% The body undergoes deformation, $\Bvarphi(\bX,t)$, and has a
%% material velocity field $\bV(\bX,t) =
%% \partial\Bvarphi(\bX,t)/\partial t$. In discussing momentum and energy, it proves convenient to define a
%% material velocity of species $\iota$ relative to the solid phase
%% as $\bV^\iota = (1/\rho_0^\iota)\bF\bM^\iota$. Recall (from
%% Section \ref{sect2}) that the remaining species are described as
%% deforming with the solid phase and diffusing relative to it.
%% Therefore $\bF$ is common to all species. The spatial velocity
%% corresponding to $\bV^\iota$ is $\bv^\iota =
%% (1/\rho^\iota)\bm^\iota = \bV^\iota$, by the Piola transform.
%% Since fluxes are defined relative to the solid tissue phase, which
%% does not diffuse, the total material velocity of the solid phase
%% is $\bV$, and for each of the remaining species it is $\bV +
%% \bV^\iota$, $\iota = \mathrm{f},\alpha,\dots,\omega$. Formally, we
%% can write the material velocity as  $\bV + \bV^\iota$, $\iota =
%% \mathrm{s},\mathrm{f},\alpha,\dots,\omega$ with the understanding
%% that $\bV^\mathrm{s} = \bzero$. Likewise, $\Pi^\mathrm{f} = 0$.
%% This convention has been adopted in the remainder of the paper.

%% We first write the balance of linear momentum in
%% the reference configuration, $\Omega_0$. 


%% %% The total first Piola-Kirchhoff stress tensor, $\bP$, is the sum of
%% %% the partial stresses $\bP^{\,\iota}$ (borne by a species $\iota$) over all
%% %% the species present.

%% %% With the introduction 
%% %%   of these quantities, the 
%% %% balance of linear momentum in local form over $\Omega_0$ for solid
%% %% collagen and fluid is,

%% The mass fluxes, $\bM^\iota$, and mass sources, $\Pi^\iota$ make
%% important contributions to the balance of linear momentum, as shown
%% below.

%% %% \noindent where $\bg$ is the body force per unit mass, and $\bq^\iota$
%% %% is an interaction term denoting the force per unit mass exerted upon
%% %% $\iota$ by all other species present. 

%% Summing over all species, the balance
%% of linear momentum for the system is obtained:

%% %%Additionally,
%% %%recognising that the rate of change of momentum of the entire tissue
%% %%is affected only by external agents and is independent of internal
%% %%interactions, the following relation arises.

%% %% \begin{equation}
%% %% \sum\limits_{\iota =
%% %%   \mathrm{c}}^\mathrm{f}\left(\rho^\iota_0\bq^\iota+\Pi^\iota
%% %% \bV^\iota 
%% %% \right)= 0. \label{qrelation}
%% %% \end{equation}

%% %% \noindent This is also consistent with Classical Mixture Theory
%% %% \citep{TruesdellNoll:65}. See \citet{growthpaper} for further
%% %% details on balance of linear momentum, and the formulation of
%% %% balance of angular momentum. We only note here that the latter
%% %% principle leads to a symmetric partial Cauchy stress,
%% %% $\Bsigma^\iota$ for each species in contrast with the unsymmetric
%% %% Cauchy stress of \cite{EpsteinMaugin:2000}.

%% Having established (\ref{interforcebalance}) we return to the
%% balance of linear momentum for a single species
%% (\ref{linmombalI1}) in order to simplify it. 

%% The balance of linear momentum for a single species in the current
%% configuration, $\Omega_t$, is obtained via similar arguments and
%% the Reynolds Transport Theorem:
%% \begin{eqnarray}
%% \rho^\iota\frac{\partial}{\partial t}\left(\bv+\bv^\iota\right)
%% &=& \rho^\iota\left(\bg+\bq^\iota\right) +
%% \Bnabla_x\cdot\Bsigma^\iota\nonumber\\
%% & & - \left(\Bnabla_x\left(\bv+\bv^\iota\right)\right)\bm^\iota -
%% \rho^\iota\left(\Bnabla_x\left(\bv+\bv^\iota\right)\right)\bv,
%% \label{ballinmomcurrI}
%% \end{eqnarray}

%% \noindent where $\Bsigma^\iota =
%% (\mathrm{det}\bF)^{-1}\bP^\iota\bF^\mathrm{T}$ is the partial
%% Cauchy stress of species $\iota$.

%% \begin{eqnarray}
%% & &\sum\limits_{\iota}\rho^\iota_0\frac{\partial e^\iota}{\partial
%% t} =\nonumber\\
%% & &\qquad\sum\limits_{\iota}\left(
%% \bP^\iota\colon\dot{\bF}+\bP^\iota\colon\Bnabla\bV^\iota -
%% \Bnabla\cdot\bQ^\iota +\rho_0^\iota r^\iota +
%% \rho^\iota_0\tilde{e}^\iota - \Bnabla e^\iota\cdot\bM^\iota\right)
%% \label{energysummation}
%% \end{eqnarray}

%% \begin{eqnarray}
%% \sum\limits_{\iota}\rho^\iota_0\frac{\partial e^\iota}{\partial t}
%% &=& \sum\limits_{\iota}\left(
%% \bP^\iota\colon\dot{\bF}+\bP^\iota\colon\Bnabla\bV^\iota -
%% \Bnabla\cdot\bQ^\iota +\rho_0^\iota r^\iota- \Bnabla
%% e^\iota\cdot\bM^\iota\right)\nonumber\\
%% &&-
%% \sum\limits_{\iota}\left(\rho_0^\iota\bq^\iota\cdot(\bV+\bV^\iota)
%% - \Pi^\iota\left(e^\iota +
%% \frac{1}{2}\Vert\bV+\bV^\iota\Vert^2\right)\right).
%% \label{energybaltot1}
%% \end{eqnarray}


%% \begin{eqnarray}
%% \sum\limits_{\iota}\frac{\mathrm{d}}{\mathrm{d}t}
%% \int\limits_{\Omega_0}\rho_0^\iota \left ( e^\iota + \frac{1}{2}
%% \Vert\bV+\bV^\iota\Vert^2 \right ) \mathrm{d}V =
%% \sum\limits_{\iota}\int\limits_{\Omega_0} \left(\rho_0^\iota\bg
%% \cdot\left(\bV+\bV^\iota\right) + \rho_0^\iota r^\iota
%% \right)\mathrm{d}V &
%% &\nonumber\\
%% +\sum\limits_{\iota}\int\limits_{\partial
%% \Omega_0}\left(\left(\bV+\bV^\iota\right)\cdot\bP^\iota -
%% \bM^\iota\left(e^\iota +\frac{1}{2}
%% \Vert\bV+\bV^\iota\Vert^2\right) -
%% \bQ^\iota\right)\cdot\bN\mathrm{d}A.\quad& & \label{energysys}
%% \end{eqnarray}

%% \noindent Substituting for $\sum\limits_\iota
%% \rho^\iota_0\tilde{e}^\iota$ from (\ref{energycond1}), 

%% The formulation up to this point has introduced some elements of
%% coupling between mass transport, mechanics and thermodynamics.  Mass
%% transport and mechanics are further coupled due to the kinematics of
%% growth. Local volumetric changes take place as species concentrations
%% evolve. As concentration increases, the material of a species swells,
%% and conversely, shrinks as concentration decreases.

%% Finite strain kinematics treats the total deformation gradient as
%% arising from a geometrically-necessary elastic deformation
%% accompanying growth, as well as a separate elastic deformation due
%% to an external stress. The deformation gradient is subject to a
%% split reminiscent of the classical decomposition of multiplicative
%% plasticity

%% This
%% sequence of maps is pictured in Figure \ref{growthkinematicsfig}.

%%  and rationale
%% underlying it is elucidated below.
%% At a continuum point the reference concentration of each species
%% admits the notion of an ``original'' state in which the
%% concentration of a species is $\rho_\mathrm{org}^\iota(\bX)$. This
%% is a state that may never be attained in a physical system.
%% However, if attained, the corresponding species would be
%% stress-free in the absence of deformation. Neglecting other
%% possible kinematics (such as plasticity) and microstructural
%% details, the set of quantities
%% $\{\rho_0^\mathrm{s},\dots,\rho_0^\omega\}$, and the temperature,
%% $\theta$, fully specify the reference state of the material at a
%% point. As mass transport alters the reference density to its value
%% $\rho_0^\iota(\bX,t)$, the species swells if $\rho_0^\iota >
%% \rho_\mathrm{org}^\iota$, and shrinks if $\rho_0^\iota <
%% \rho_\mathrm{org}^\iota$. Assuming that these volume changes are
%% isotropic leads to the following growth kinematics: For each
%% species, one can define a ``growth deformation gradient tensor'',
%% $\bF^{\mathrm{g}^\iota} :=
%% \frac{\rho_0^\iota}{\rho_\mathrm{org}^\iota}{\bf 1}$, where ${\bf
%% 1}$ is the second-order isotropic tensor. The tensor
%% $\bF^{\mathrm{g}^\iota}$ is analogous to the plastic deformation
%% gradient of multiplicative plasticity. 

%% With an external stress, there is further elastic
%% deformation, $\bar{\bF}^\mathrm{e}$, common to all species.

%% The kinematic relations are:
%% \begin{equation}
%% \bF =
%% \bar{\bF}^\mathrm{e}\tilde{\bF}^{\mathrm{e}^\iota}\bF^{\mathrm{g}^\iota},\quad
%% \bF^{\mathrm{g}^\iota} =
%% \frac{\rho_0^\iota}{\rho_\mathrm{org}^\iota}{\bf 1}.
%% \label{growthkinematicseq}
%% \end{equation}

%% \noindent Clearly, the elastic deformation gradients can be
%% combined to write $\bF^{\mathrm{e}^\iota} =
%% \bar{\bF}^\mathrm{e}\tilde{\bF}^{\mathrm{e}^\iota}$, the ``total''
%% elastic deformation gradient of species $\iota$.
%% \begin{figure}[ht]
%% \psfrag{A}{\small $\Omega_0$} \psfrag{B}{\small $\Omega^\ast$}
%% \psfrag{C}{\small $\Omega_t$} \psfrag{D}{\small $\Bvarphi$}
%% \psfrag{G}{\small $\bu^\ast$} \psfrag{E}{\small $\Bkappa$}
%% \psfrag{M}{\small $\bX$} \psfrag{I}{\small
%% $\bF^{\mathrm{g}^\iota}$} \psfrag{H}{\small
%% $\tilde{\bF}^{\mathrm{e}^\iota}$} \psfrag{J}{\small $\tilde{\bF}$}
%% \psfrag{Y}{\small $\bX^\ast$} \psfrag{K}{\small
%% $\bar{\bF}^\mathrm{e}$} \psfrag{X}{\small $\bx$} \psfrag{L}{\small
%% $\bF$} \centering
%% {\includegraphics[width=0.8\textwidth]
%%   {images/elucidation/kinematics}}
%% \caption{The kinematics of growth.} \label{growthkinematicsfig}
%% \end{figure}

%% %% As is customary in field theories of continuum physics, the
%% %% Clausius-Duhem inequality is obtained by multiplying the Entropy
%% %% Inequality (the Second Law of Thermodynamics) by the temperature
%% %% field, $\theta$, and subtracting it from the Balance of Energy (the
%% %% First Law of Thermodynamics). 

%% %%  Only the valid
%% %% constitutive laws relevant to the examples that follow are listed
%% %% here. For details, see \cite{growthpaper}.


%% \noindent Inequality (\ref{redentropyineq1}) represents a
%% fundamental restriction upon the physical processes during
%% biological growth. Any constitutive relations that are prescribed
%% must satisfy this restriction,

%% {\bf The Eshelby stress as a thermodynamic driving force}
%% \label{eshelby-force}

%% Combining the stress divergence and chemical potential gradient
%% contributions to the driving force for any species, and using the
%% mass-specific Helmholtz free energy, $\psi^\iota$, we write,

%% \begin{equation}
%% \bF^\mathrm{T}\Bnabla\cdot\bP^\iota - \rho^\iota_0\left(\Bnabla
%% e^\iota - \theta\Bnabla\eta^\iota\right) =
%% \Bnabla\cdot\left(\bF^\mathrm{T}\bP^\iota\right) -
%% \Bnabla\bF^\mathrm{T}\colon\bP^\iota
%%  - \rho^\iota_0\Bnabla\psi^\iota\vert_\theta.
%% \end{equation}

%% \noindent Regrouping terms this expression is
%% \begin{equation}
%% -\Bnabla\cdot\underbrace{\left(\rho^\iota_0\psi^\iota\vert_\theta\bone
%% - \bF^\mathrm{T}\bP^\iota\right)}_{\mbox{Eshelby
%% stress},\;\BXi^\iota} +
%% \left(\Bnabla\rho^\iota_0\right)\psi^\iota_0\vert_\theta -
%% \Bnabla\bF^\mathrm{T}\colon\bP^\iota.
%% \end{equation}

%% Thus, the divergence of the well-known Eshelby stress tensor is
%% also among the driving forces for mass transport. Also observe the
%% presence of a strain gradient-dependent driving force,
%% $-\Bnabla\bF^\mathrm{T}\colon\bP^\iota$ in the developments of
%% Sections \ref{sect5.2} and \ref{sect5.3}, independent of the
%% pressure gradient term for the fluid species.

%%  Recall that in (\ref{VIconstrel1}) we 
%% have taken $\tilde{\bD}^\iota$ to be positive
%% semi-definite\footnote{If
%% $\Vert\tilde{\bD}^\iota\Bnabla\bV\bF^{-1}/\rho^\iota_0\Vert << 1$,
%% we have $\bD^\iota \approx \tilde{\bD}^\iota$.}. 

%%  We proceed now to examine the four separate terms in the
%% thermodynamic driving force:
%% \begin{equation}
%% \boldmath{\sF}^\iota =
%% -\rho^\iota_0\bF^\mathrm{T}\frac{\partial\bV}{\partial t} +
%% \rho_0^\iota\bF^\mathrm{T}\bg +
%% \bF^\mathrm{T}\Bnabla\cdot\bP^\iota - \rho^\iota_0\left(\Bnabla
%% e^\iota - \theta\Bnabla\eta^\iota\right) \label{drivingforceI}
%% \end{equation}

%% \noindent The first two terms respectively represent the
%% influences of inertia and body force. Thus, the inertial effect is
%% to drive species~$\iota$ in the opposite direction to the body's
%% acceleration. 

%% The body force's influence is directed along itself.

%% The third term represents the stress divergence effect. . We
%% demonstrate this effect for the case of the fluid species in
%% Section \ref{sect5.2} below, for which it translates to the more
%% intuitive notion of transport along a fluid pressure gradient.


%% The fourth term in $\boldmath{\sF}^\iota$ admits the following
%% interpretation: The . An assumption inherent in the development
%% that began in Section \ref{sect2} is that any mass entering or
%% leaving $\Omega_0$ at a point $\bX$ on the boundary,
%% $\partial\Omega_0$, has the field values
%% $\rho_0^\iota,e^\iota,\eta^\iota,\theta$, and $\psi^\iota$
%% corresponding to $\bX$. Likewise, the incremental mass of species
%% $\iota$ created or absorbed via the source/sink $\Pi^\iota$ at
%% $\bX$ has the field values of that point. Consider a sufficiently
%% small neighborhood of a point, say
%% $\mathsf{N}(\bX)\subset\Omega_0$. Changing the mass of species
%% $\iota$ in $\mathsf{N}(\bX)$ by $\delta\mathsf{m}^\iota$ units
%% causes a change in the Helmholtz free energy of $\iota$ in
%% $\mathsf{N}(\bX)$ by $\delta\Psi^\iota =
%% \psi^\iota\delta\mathsf{m}^\iota$. By definition therefore,
%% $\psi^\iota =
%% \partial\Psi^\iota/\partial\mathsf{m}^\iota$. This derivative gives the \emph{chemical
%% potential}, $\mu^\iota$, of the transported species, $\iota$.
%% Thus, we have $\mu^\iota = e^\iota - \theta\eta^\iota$, and
%% $\Bnabla e^\iota -\theta\Bnabla\eta^\iota =
%% \Bnabla\mu^\iota\vert_\theta$. This last term in
%% ${\boldmath{\sF}}^\iota$ thus represents the thermodynamic driving
%% force due to a chemical potential gradient.


%% The gradient of internal energy in (\ref{drivingforceI}) leads to
%% a strain gradient-dependent term. A concentration gradient-driven
%% term arises from the gradient of mixing entropy. Together with the
%% other terms that were remarked upon above, they represent a
%% complete thermodynamic formulation of coupled mass transport and
%% mechanics. This is the central result of our paper.


%% It has recently come to our attention that the constitutive
%% relation for flux (\ref{MIconstrel}) is precisely the result
%% arrived at by \citet{DeGrootMazur:1984}, including the
%% identification of the chemical potential gradient term. However,
%% their approach involves a slightly different application of the
%% Second Law, and a less detailed treatment of the mechanics.

%% \noindent{\bf Remark 1}: The final version of the dissipation
%% inequality (\ref{dissipation1}), and the mass balance equation can
%% be manipulated to restrict the mathematical form of the mass
%% source. It is common to make the mass source depend upon the
%% strain energy density \citep{HarriganHamilton:1993} while
%% respecting the restriction imposed by the dissipation inequality.
%% This form is often used while modelling hard tissue. Such an
%% approach leads to strain-mediated mass transport. However, with a
%% strain-independent source, strain-mediated (or stress-mediated)
%% mass transport would not be obtained with such a formulation.

%% \noindent{\bf Remark 2}: We expect that evaluation of the
%% dissipation, $\sD$, using (\ref{dissipation1}) from field
%% quantities in a boundary value problem will provide a test of
%% soundness, and if necessary indications for improvement, of our
%% constitutive models.

%% In our earlier exposition \citep{growthpaper}, we
%% used simple first order chemical kinetics to define
%% $\Pi^{\mathrm{c}}$. Other forms, which have been studied in the
%% literature, can be used:
