~hnarayanan/phd-dissertation/trunk : revision 128

9

derived from general principles governing the behaviour of multiphase

10

mixtures. Specifically, Section~\ref{balance-laws} helps define the

11

system and formally introduces fundamental quantities characterising

12

it, before deriving the balance laws from fundamental

13

axioms. Section~\ref{kinematics-of-growth} presents the kinematics

14

associated with finite deformation growth. A fundamental axiom of

15

Thermodynamics, the entropy inequality, and the restrictions it places

16

on functional forms of constitutive relationships is the subject of

17

Section~\ref{entropy-inequality}. The chapter concludes with key

18

algorithmic considerations (Section~\ref{algorithmic-considerations})

19

which play an important role in the computational formulation

20

underlying the numerical experiments presented in

21

Chapter~\ref{numerical-simulations-1}.

12

it, before deriving the balance laws from fundamental axioms. Section%

13

~\ref{kinematics-of-growth} presents the kinematics associated with

14

finite deformation growth. A fundamental axiom of Thermodynamics, the

15

entropy inequality, and the restrictions it places on functional forms

16

of constitutive relationships is the subject of Section~\ref{entropy%

17

-inequality}. The chapter concludes with key algorithmic

18

considerations (Section~\ref{algorithmic-considerations}) which play

19

an important role in the computational formulation underlying the

20

numerical experiments presented in Chapter~\ref{numerical-simulations%

21

-1}.

22

23

\section{Balance laws for an open mixture}

24

\label{balance-laws}

57

58

The tissue consists of numerous species, of which the following

59

groupings are of importance for the models: A solid species,

60

consisting of solid \emph{collagen fibrils} and

61

\emph{cells},\footnote{At this point, the solid species is not

62

differentiated any further. This is a good approximation to the

63

physiological setting for tendons, which are relatively acellular

64

and whose dry mass consists of up to 75\% collagen

65

\citep{Nordinetal:2001}. When modelling tumour growth in a later

66

chapter (Section~\ref{tumour-growth}), where cell mechanics and

67

migration are significant \citep{namyetal:04}, the solid phase is

68

further distinguished.} denoted by $\mathrm{c}$, an extra-cellular

69

\emph{fluid} species, denoted by $\mathrm{f}$, consisting primarily of

70

water, and \emph{solute} species, consisting of precursors to

71

reactions, byproducts, nutrients, and other regulatory chemicals. A

72

generic solute will be denoted by $\mathrm{s}$. In the treatment that

73

follows, an arbitrary species will be denoted by $\iota$, where $\iota

74

= \mathrm{c,f,s}$.

60

consisting of solid \emph{collagen fibrils} and \emph{cells},%

61

\footnote{At this point, the solid species is not differentiated any

62

further. This is a good approximation to the physiological setting

63

for tendons, which are relatively acellular and whose dry mass

64

consists of up to 75\% collagen \citep{Nordinetal:2001}. When

65

modelling tumour growth in a later chapter (Section~\ref{tumour%

66

-growth}), where cell mechanics and migration are significant

67

\citep{namyetal:04}, the solid phase is further distinguished.}

68

denoted by $\mathrm{c}$, an extra-cellular \emph{fluid} species,

69

denoted by $\mathrm{f}$, consisting primarily of water, and

70

\emph{solute} species, consisting of precursors to reactions,

71

byproducts, nutrients, and other regulatory chemicals. A generic

72

solute will be denoted by $\mathrm{s}$. In the treatment that follows,

73

an arbitrary species will be denoted by $\iota$, where $\iota =

74

\mathrm{c,f,s}$.

75

76

The fundamental quantities of interest are mass concentrations,

77

$\rho_0^\iota(\bX,t)$. These are the masses of each species per unit

108

\underbrace{\frac{\mathrm{d}}{\mathrm{d}t} \int\limits_{\Omega_0}

109

\rho_0^\iota (\bX,t)\mathrm{d}V}_{\text{Rate of change of mass}} =

110

\underbrace{\int\limits_{\Omega_0} \Pi^\iota

111

(\bX,t)\mathrm{d}V}_{\text{Mass being created}}

112

-\underbrace{\int\limits_{\partial\Omega_0}\bM^\iota(\bX,t)\cdot\bN

111

(\bX,t)\mathrm{d}V}_{\text{Mass being created}} -

112

\underbrace{\int\limits_{\partial\Omega_0} \bM^\iota(\bX,t)\cdot\bN

113

\mathrm{d}A}_{\text{Mass leaving the domain}},

114

\label{globalbalanceofmass}

115

\end{equation}

193

\psfrag{I}{$\bn\cdot\bm^\iota$}

194

\psfrag{J}{$\bg$}

195

\psfrag{K}{$\bq^\iota$}

196

\psfrag{L}{$\bP\bN$}

197

\psfrag{O}{$\Bsigma\bn$}

196

\psfrag{L}{$\bP^\iota\bN$}

197

\psfrag{O}{$\Bsigma^\iota\bn$}

198

\includegraphics[width=0.8\textwidth]{images/elucidation/%

199

continuum-potato-momentum}

200

\caption{Interaction forces, traction and body loads on the tissue.}

303

The final term with the gradient of total species velocity identifies

304

the contribution of the flux to the balance of momentum. In practise,

305

when the relative magnitude of the fluid mobility (and hence flux) is

306

small, the final term on the right hand-side of

307

Equation~(\ref{localbalanceofmomentum}) is negligible, resulting in a

308

more classical form of the balance of momentum. Furthermore, in the

309

absence of significant acceleration of the tissue during growth, the

310

left hand-side can also be neglected, reducing

311

(\ref{localbalanceofmomentum}) to the quasi-static balance of linear

312

momentum.

306

small, the final term on the right hand-side of Equation~(\ref{local%

307

balanceofmomentum}) is negligible, resulting in a more classical

308

form of the balance of momentum. Furthermore, in the absence of

309

significant acceleration of the tissue during growth, the left

310

hand-side can also be neglected, reducing (\ref{localbalanceofmoment%

311

um}) to the quasi-static balance of linear momentum.

313

312

314

313

The interaction forces, $\bq^\iota$, satisfy a relation with the mass

315

314

sources, $\Pi^\iota$, that is elucidated by the following

358

357

\left(\bV+\bV^\iota\right) \mathrm{d}V = 0,\nonumber

359

358

\end{eqnarray}

360

359

361

\noindent which, upon localisation (recalling

362

Equation~(\ref{masssummationresult})) , leads to

360

\noindent which, upon localisation (recalling Equation~(\ref{masss%

361

ummationresult})), leads to

363

362

364

363

\begin{eqnarray}

365

364

\sum\limits_{\iota}\left(\rho^\iota_0\bq^\iota +

422

421

\end{eqnarray}

423

422

424

423

\noindent where $\Bepsilon$ is the permutation symbol, and

425

$\Bepsilon\colon\bA$ is written as $\epsilon_{ijk}A_{jk}$ in indicial

426

form, for any second-order tensor $\bA$. Using the mass balance

427

equation (\ref{localbalanceofmass}), and balance of linear momentum

424

$\Bepsilon\colon\bA$ is written as $\epsilon_{ijk} A_{jk} \be_i$ in

425

indicial form for any second-order tensor $\bA$; $\be_i$ being the

426

$i^{\mathrm{th}}$ basis vector. Using the mass balance equation

427

(\ref{localbalanceofmass}), and balance of linear momentum

428

(\ref{localbalanceofmomentum}), we have

429

430

\begin{displaymath}

571

\rho^\iota_0\frac{\partial e^\iota}{\partial t} &=&

572

\bP^\iota\colon\mathrm{GRAD}\left[\bV+\bV^\iota\right]-

573

\mathrm{DIV}\left[\bQ^\iota\right] + \rho_0^\iota r^\iota

574

+\rho^\iota_0\tilde{e}^\iota\nonumber\\

575

& & - \mathrm{GRAD}\left[e^\iota\right]\cdot\bM^\iota,

576

\;\forall\,\iota.

574

+\rho^\iota_0\tilde{e}^\iota\nonumber\\ & & -

575

\mathrm{GRAD}\left[e^\iota\right]\cdot\bM^\iota, \;\forall\,\iota.

577

576

\label{localbalanceofenergy}

578

577

\end{eqnarray}

579

578

631

630

localisation that,

632

631

633

632

\begin{eqnarray}

634

\sum\limits_{\iota}\left(\rho_0^\iota\bq^\iota\cdot(\bV+\bV^\iota)

635

+ \Pi^\iota\left(e^\iota +

636

\frac{1}{2}\Vert\bV+\bV^\iota\Vert^2\right) +

637

\rho^\iota_0\tilde{e}^\iota\right) =

638

0. \label{energysummationresult}

633

\sum\limits_{\iota}\left(\rho_0^\iota\bq^\iota\cdot(\bV+\bV^\iota) +

634

\Pi^\iota\left(e^\iota + \frac{1}{2}\Vert\bV+\bV^\iota\Vert^2\right) +

635

\rho^\iota_0\tilde{e}^\iota\right) = 0. \label{energysummationresult}

639

636

\end{eqnarray}

640

637

641

This result, relating the interaction energies to interaction

642

forces between species, their sources and relative velocities, is

643

identical to that obtained from classical mixture theory

644

\citep{TruesdellNoll:65}, ensuring consistency of the present

645

formulation with mixture theory.

638

This result, relating the interaction energies to interaction forces

639

between species, their sources and relative velocities, is identical

640

to that obtained from classical mixture theory \citep{TruesdellNoll%

641

:65}, ensuring consistency of the present formulation with mixture

642

theory.

646

643

647

644

\section{The kinematics of growth}

648

645

\label{kinematics-of-growth}

682

679

notion of the \emph{growth component} of the deformation

683

680

gradient. This observation has led to an active field of study within

684

681

the literature on biological growth \citep{Skalak:81, SkalakHoger:96,

685

Klischetal:2001, TaberHumphrey:2001, LubardaHoger:02,

686

AmbrosiMollica:2002}, and the treatment below follows in the same

687

vein.

682

Klischetal:2001, TaberHumphrey:2001, LubardaHoger:02, AmbrosiMol%

683

lica:2002}, and the treatment below follows in the same vein.

688

684

689

685

In the setting of finite strain kinematics, the total deformation

690

686

gradient, $\bF$, is decomposed into the growth component of the solid

691

collagen, $\bF^{\mathrm{g}^\mathrm{c}}$, a

692

\emph{geometrically-necessitated elastic component} accompanying

693

growth, $\widetilde{\bF}^{\mathrm{e}^\mathrm{c}}$ and an

694

\emph{additional elastic component due to external stress},

695

$\overline{\bF}^{\mathrm{e}^\mathrm{c}}$. Later, we will write

696

$\bF^{\mathrm{e}^\mathrm{c}} = \overline{\bF}^{\mathrm{e}^\mathrm{c}}

697

\widetilde{\bF}^{\mathrm{e}^\mathrm{c}}$. This elastic-growth

698

decomposition is visualised in

699

Figure~\ref{continuum-potato-growth-kinematics} and is elaborated upon

700

below.

687

collagen, $\bF^{\mathrm{g}^\mathrm{c}}$, a \emph{geometrically-nece%

688

ssitated elastic component} accompanying growth, $\widetilde{\bF} ^

689

{\mathrm{e}^\mathrm{c}}$ and an \emph{additional elastic component due

690

to external stress}, $\overline{\bF} ^ {\mathrm{e}^\mathrm{c}}$.

691

Later, we will write $\bF^{\mathrm{e}^\mathrm{c}} = \overline{\bF} ^

692

{\mathrm{e}^\mathrm{c}} \widetilde{\bF}^{\mathrm{e}^\mathrm{c}}$. This

693

elastic-growth decomposition is visualised in Figure~\ref{continuum%

694

-potato-growth-kinematics} and is elaborated upon below.

701

695

702

696

This split of the total deformation gradient is analogous to the

703

classical decomposition of multiplicative plasticity

704

\citep{Bilbyetal:1956,Lee:1969}. As explained in Section

705

\ref{role-of-current-mass-balance}, we assume that the fluid-filled

706

pores also deform with $\bF$, and that a component,

707

$\bF^{\mathrm{e}^\mathrm{f}}$, of this total deformation gradient

708

tensor, determines the fluid stress. We also assume a fluid growth

709

component, $\bF^{\mathrm{g}^\mathrm{f}}$, which is detailed below, and

710

that $\bF^{\mathrm{e}^\mathrm{f}}\bF^{\mathrm{g}^\mathrm{f}} =

711

\bF$. As with the solid collagen we admit $\bF^{\mathrm{e}^\mathrm{f}}

712

= \overline{\bF}^{\mathrm{e}^\mathrm{f}}

713

\widetilde{\bF}^{\mathrm{e}^\mathrm{f}}$, the sub-components carrying

714

the same interpretation as for the solid collagen. However, this last

715

decomposition is not explicitly used.

697

classical decomposition of multiplicative plasticity \citep{Bilb%

698

yetal:1956,Lee:1969}. As explained in Section \ref{role-of-cur%

699

rent-mass-balance}, we assume that the fluid-filled pores also

700

deform with $\bF$, and that a component, $\bF^{\mathrm{e} ^

701

\mathrm{f}}$, of this total deformation gradient tensor, determines

702

the fluid stress. We also assume a fluid growth component,

703

$\bF^{\mathrm{g}^\mathrm{f}}$, which is detailed below, and that

704

$\bF^{\mathrm{e}^\mathrm{f}}\bF^{\mathrm{g}^\mathrm{f}} = \bF$. As

705

with the solid collagen we admit $\bF^{\mathrm{e}^\mathrm{f}} =

706

\overline{\bF}^{\mathrm{e}^\mathrm{f}} \widetilde{\bF}^{\mathrm{e} ^

707

\mathrm{f}}$, the sub-components carrying the same interpretation as

708

for the solid collagen. However, this last decomposition is not

709

explicitly used.

716

710

717

711

Assuming that the volume changes associated with growth described

718

712

above are isotropic, a simple form for the growth part of the

786

780

It is worth emphasising that this argument holds for

787

781

$\bF^{\mathrm{g}^\mathrm{f}}$, which is the local stress-free state of

788

782

deformation of the fluid-containing pores at a point. The actual

789

deformation gradient, $\bF =

790

\bF^{\mathrm{e}^\mathrm{f}}\bF^{\mathrm{g}^\mathrm{f}}$, also depends

791

on the the elastic part, $\bF^{\mathrm{e}^\mathrm{f}}$, which is

792

determined by the constitutive response of the fluid. Under stress, an

793

incompressible fluid will have

794

$\mathrm{det}(\bF^{\mathrm{e}^\mathrm{f}}) = 1$, where

795

$\mathrm{det}(\bullet)$ denotes the determinant of a second-order

796

tensor. Therefore, a fluid-saturated tissue will swell with fluid

797

influx, \mbox{$\mathrm{det}(\bF) =

798

\mathrm{det}(\bF^{\mathrm{g}^\mathrm{f}}) > 1$}. A compressible

799

fluid may have

800

$\mathrm{det}(\bF^{\mathrm{e}^\mathrm{f}}) < 1$ allowing

801

$\mathrm{det}(\bF) < 1$ even with

802

$\mathrm{det}(\bF^{\mathrm{g}^\mathrm{f}}) >1$. But, even in this case,

803

in the stress-free state there will be swelling.

783

deformation gradient, $\bF = \bF^{\mathrm{e}^\mathrm{f}} \bF ^

784

{\mathrm{g}^\mathrm{f}}$, also depends on the elastic part,

785

$\bF^{\mathrm{e}^\mathrm{f}}$, which is determined by the constitutive

786

response of the fluid. Under stress, an incompressible fluid will have

787

$\mathrm{det}(\bF^{\mathrm{e}^\mathrm{f}}) = 1$, where $\mathrm{det}(%

788

\bullet)$ denotes the determinant of a second-order tensor. Therefore,

789

a fluid-saturated tissue will swell with fluid influx, \mbox{$\mathrm%

790

{det}(\bF) = \mathrm{det}(\bF^{\mathrm{g}^\mathrm{f}}) > 1$}. A

791

compressible fluid may have $\mathrm{det}(\bF^{\mathrm{e}^\mathrm{f}})

792

< 1$ allowing $\mathrm{det}(\bF) < 1$ even with $\mathrm{det}(\bF ^

793

{\mathrm{g}^\mathrm{f}}) >1$. But, even in this case, in the

794

stress-free state there will be swelling.

804

795

805

796

Thus, for the fluid phase, the isotropic swelling law can be extended

806

797

to the unsaturated case by introducing a degree of saturation,

811

802

intrinsic density in $\Omega_t$ and is given by $\tilde{\rho}^\iota =

812

803

\tilde{\rho}^\iota_0/\mathrm{det}(\bF)$. Note that the intrinsic

813

804

reference density, $\tilde{\rho}^\iota_0$, is a material

814

property. Upon solution of the mass balance equation

815

(\ref{current-mass-balance}) for $\rho^\iota$, the species volume

816

fractions, $\tilde{v}^\iota$, can be computed in a straightforward

817

fashion. The sum of these volume fractions is our required measure of

818

saturation defined in $\Omega_t$. We also recognise that for the

819

dilute solutions obtained with physiologically-relevant solute

820

concentrations, the saturation condition is very well approximated by

821

$\tilde{v}^\mathrm{f} + \tilde{v}^\mathrm{c} = 1$. So, we proceed to

822

redefine the fluid growth-induced component of the pore deformation

823

gradient tensor as follows:

805

property. Upon solution of the mass balance equation (\ref{current%

806

-mass-balance}) for $\rho^\iota$, the species volume fractions,

807

$\tilde{v}^\iota$, can be computed in a straightforward fashion. The

808

sum of these volume fractions is our required measure of saturation

809

defined in $\Omega_t$. We also recognise that for the dilute solutions

810

obtained with physiologically-relevant solute concentrations, the

811

saturation condition is very well approximated by $\tilde{v} ^

812

\mathrm{f} + \tilde{v}^\mathrm{c} = 1$. So, we proceed to redefine the

813

fluid growth-induced component of the pore deformation gradient tensor

814

as follows:

824

815

825

816

\begin{equation}

826

817

\bF^{\mathrm{g}^\mathrm{f}} = \left\{ \begin{array}{ll} \left(

841

832

non-physical. Saturation holds in the sense that $\tilde{v}^\mathrm{f}

842

833

+ \tilde{v}^\mathrm{c} = 1$. It has been common in soft tissue

843

834

literature to assume that, under normal physiological conditions, soft

844

tissues are fully saturated by the fluid and

845

\mbox{Equation~(\ref{isotropicgrowth})} is appropriate for $\iota =

846

\mathrm{f}$. However, this treatment of saturation and swelling

847

induced by the fluid phase is necessary background for Section

848

\ref{caviation-under-tension}, where we examine the response of the

849

fluid phase under tension. This approach also holds relevance for

850

partial drying, which \emph{ex vivo} or \emph{in vitro} tissue may be

851

subject to under certain laboratory conditions. It also significantly

852

expands the relevance of the formulation by making it applicable to

853

the mechanics of drained porous media other than biological tissue;

854

most prominently, soils.

835

tissues are fully saturated by the fluid and \mbox{Equation~(\ref%

836

{isotropicgrowth})} is appropriate for $\iota = \mathrm{f}$.

837

However, this treatment of saturation and swelling induced by the

838

fluid phase is necessary background for Section \ref{caviation-under%

839

-tension}, where we examine the response of the fluid phase under

840

tension. This approach also holds relevance for partial drying, which

841

\emph{ex vivo} or \emph{in vitro} tissue may be subject to under

842

certain laboratory conditions. It also significantly expands the

843

relevance of the formulation by making it applicable to the mechanics

844

of drained porous media other than biological tissue; most

845

prominently, soils.

855

846

856

847

\section{The entropy inequality and its restrictions on constitutive

857

848

relations}

880

871

\int\limits_{\Omega_0} \rho_0^\iota \eta^\iota

881

872

\mathrm{d}V}_{\text{Rate of change of entropy}} &\geq&

882

873

\underbrace{\sum\limits_{\iota}\int\limits_{\Omega_0}

883

\frac{\rho_0^\iota r^\iota}{\theta}

884

\mathrm{d}V}_{\text{Entropy added through heat supply}}\nonumber\\ &

885

&- \underbrace{\sum\limits_{\iota}\int\limits_{\partial \Omega_0}

874

\frac{\rho_0^\iota r^\iota}{\theta} \mathrm{d}V}_{\text{Entropy

875

added through heat supply}}\nonumber\\ & &-

876

\underbrace{\sum\limits_{\iota}\int\limits_{\partial \Omega_0}

886

877

\left(\bM^\iota \cdot \bN\eta^\iota + \frac{\bQ^\iota}{\theta} \cdot

887

878

\bN \right)\mathrm{d}A}_{\text{Entropy lost through mass and heat

888

879

flux}}.

896

887

897

888

\begin{eqnarray}

898

889

\sum\limits_{\iota}\left(\rho_0^\iota\frac{\partial\eta^\iota}{\partial

899

t} + \Pi^{\iota}\eta^{\iota}

900

\right) & & \geq

901

\sum\limits_{\iota}\left(\frac{\rho_0^\iota

902

r^\iota}{\theta} -\mathrm{GRAD}\left[\eta^\iota\right]\cdot\bM^\iota

903

\right) \nonumber \\ & & - \sum\limits_{\iota}\left(

890

t} + \Pi^{\iota}\eta^{\iota} \right) & & \geq

891

\sum\limits_{\iota}\left(\frac{\rho_0^\iota r^\iota}{\theta}

892

-\mathrm{GRAD}\left[\eta^\iota\right]\cdot\bM^\iota \right) \nonumber

893

\\ & & - \sum\limits_{\iota}\left(

904

894

\frac{\mathrm{DIV}\left[\bQ^\iota\right]}{\theta} -

905

895

\frac{\mathrm{GRAD}\left[\theta\right]\cdot\bQ^\iota}{\theta^2}\right).

906

896

\label{localentropyinequality}

950

940

the simplest possible (incorporating one field variable from each of

951

941

the different kinds of physics considered: Mechanics, heat transfer

952

942

and mass transport) and results in a restricted class of

953

constitutive relationships. As seen in

954

Section~\ref{eu-entropy-inequality}, mass-specific Helmholtz

955

free energies of species dependent upon other variables, such as

956

internal variables arising from mechanics, lead to a more general

957

class of constitutive relationships, such as viscoelastic materials

958

(Section~\ref{eu-viscoelastic-solid}).} $e^{\iota} =

959

\hat{e}^\iota(\bF^{\mathrm{e}^\iota}, \eta^{\iota}, \rho_0^\iota)$ and

960

substitute this into the Clausius-Duhem inequality. Upon applying the

961

chain rule of differentiation and regrouping some terms,

962

(\ref{clausiusduhemform}) becomes,

943

constitutive relationships. As seen in Section~\ref{eu-entropy%

944

-inequality}, mass-specific Helmholtz free energies of species

945

dependent upon other variables, such as internal variables arising

946

from mechanics, lead to a more general class of constitutive

947

relationships, such as viscoelastic materials (Section~\ref{eu%

948

-viscoelastic-solid}).} $e^{\iota} = \hat{e}^\iota(\bF^{\mathrm{e}

949

^ \iota}, \eta^{\iota}, \rho_0^\iota)$ and substitute this into the

950

Clausius-Duhem inequality. Upon applying the chain rule of

951

differentiation and regrouping some terms, (\ref{clausiusduhemform})

952

becomes,

963

953

964

\begin{eqnarray}

965

& &\sum\limits_{\iota}\left(\rho_{0}^{\iota}\frac{\partial

954

\begin{eqnarray} &

955

&\sum\limits_{\iota}\left(\rho_{0}^{\iota}\frac{\partial

966

956

e^\iota}{\partial\bF^{\mathrm{e}^\iota}} -

967

957

\bP^\iota\bF^{\mathrm{g}^{\iota\mathrm{T}}}\right)

968

958

\colon\dot{\bF}^{\mathrm{e}^\iota} +

975

965

\mathrm{DIV}\left[\bP^\iota\right] +

976

966

\mathrm{GRAD}\left[\bV+\bV^\iota\right]\bM^\iota\right)

977

967

\cdot(\bV^\iota+\bV)\nonumber\\ &+&

978

\sum\limits_{\iota}\left(\rho^\iota_0\bF^{-\mathrm{T}}

968

\sum\limits_{\iota}\left(\rho^\iota_0\bF^{-\mathrm{T}}

979

969

\left(\mathrm{GRAD}\left[e^\iota\right] -

980

\theta\ \mathrm{GRAD}\left[\eta^\iota\right]\right)

981

\right)\cdot\bV^\iota\nonumber\\

982

&+&\sum\limits_{\iota}\Pi^\iota\left(\psi^\iota

983

+ \frac{1}{2}\Vert\bV+\bV^\iota\Vert^2\right) +

984

\sum\limits_{\iota} \frac{\mathrm{GRAD}\left[\theta\right]

985

\cdot\bQ^\iota}{\theta}\nonumber\\ & & +

986

\sum\limits_{\iota}\rho^\iota_0\frac{\partial

987

e^\iota}{\partial\rho^\iota_0}\frac{\partial

988

\rho^\iota_0}{\partial t} -

989

\sum\limits_{\iota}\bP^\iota

990

\colon(\mathrm{GRAD}\left[\bV^\iota\right] +

991

\bF^{\mathrm{e}^\iota}\dot{\bF}^{\mathrm{g}^\iota}) \leq 0,

992

\label{reducedentropyinequality-1}

993

\end{eqnarray}

970

\theta\ \mathrm{GRAD}\left[\eta^\iota\right]\right)

971

\right)\cdot\bV^\iota\nonumber\\ &+&\sum\limits_{\iota}\Pi^\iota\left(\psi^\iota

972

+ \frac{1}{2}\Vert\bV+\bV^\iota\Vert^2\right) + \sum\limits_{\iota}

973

\frac{\mathrm{GRAD}\left[\theta\right]

974

\cdot\bQ^\iota}{\theta}\nonumber\\ & & +

975

\sum\limits_{\iota}\rho^\iota_0\frac{\partial

976

e^\iota}{\partial\rho^\iota_0}\frac{\partial

977

\rho^\iota_0}{\partial t} - \sum\limits_{\iota}\bP^\iota

978

\colon(\mathrm{GRAD}\left[\bV^\iota\right] +

979

\bF^{\mathrm{e}^\iota}\dot{\bF}^{\mathrm{g}^\iota}) \leq

980

0, \label{reducedentropyinequality-1}

981

\end{eqnarray}

994

982

995

983

\noindent which represents a fundamental restriction upon the physical

996

984

processes underlying biological growth. Any constitutive relationships

1025

1013

1026

1014

\begin{equation}

1027

1015

\begin{split}

1028

\rho^\iota_0\bV^\iota = & -\frac{\tilde{\bD}^\iota}{\rho^\iota_0} \left(

1029

\rho_0^\iota\bg - \mathrm{DIV}\left[\bP^\iota\right] +

1030

\mathrm{GRAD}\left[\bV\right]\bM^\iota\right)\\

1031

& -\frac{\tilde{\bD}^\iota}{\rho^\iota_0}\left(\rho^\iota_0

1032

\bF^{-\mathrm{T}}\left( \mathrm{GRAD}\left[ e^\iota\right]

1033

- \theta\ \mathrm{GRAD}\left[\eta^\iota\right]\right)\right),

1016

\rho^\iota_0\bV^\iota = & -\frac{\tilde{\bD}^\iota}{\rho^\iota_0}

1017

\left( \rho_0^\iota\bg - \mathrm{DIV}\left[\bP^\iota\right] +

1018

\mathrm{GRAD}\left[\bV\right]\bM^\iota\right)\\ &

1019

-\frac{\tilde{\bD}^\iota}{\rho^\iota_0}\left(\rho^\iota_0

1020

\bF^{-\mathrm{T}}\left( \mathrm{GRAD}\left[ e^\iota\right] -

1021

\theta\ \mathrm{GRAD}\left[\eta^\iota\right]\right)\right),

1034

1022

\label{constitutive-relation-flux}

1035

1023

\end{split}

1036

1024

\end{equation}

1060

1048

conduction.

1061

1049

1062

1050

With these constitutive relations (\ref{constitutive-relation%

1063

-hyperelasticity}--\ref{constitutive-relation-heat-conduction})

1051

-hyperelasticity}--\ref{constitutive-relation-heat-conduction})

1064

1052

ensuring that certain terms of the dissipation inequality vanish,

1065

1053

(\ref{reducedentropyinequality-1}) is further reduced to

1066

1054

1145

1133

\citet{LandLif} for a general formulation of statistical mechanics

1146

1134

models of long chain molecules. Fitting the WLC response function

1147

1135

derived by \citet{MarkoSiggia:95} to the collagen {\em fibril} data of

1148

\citet{Grahametal:2004} results in values of $A = 6$ nm and $L = 3480$

1149

nm. This is to be compared with $A=14.5$ nm and $L=309$ nm, reported

1150

by \citet{Sunetal:2002}, for a {\em single} collagen molecule. Taken

1151

together, these results demonstrate that the WLC analysis correctly

1152

predicts a collagen fibril to be longer and more compliant than its

1153

constituent molecule due to compliant intermolecular cross-links in

1154

fibrils.

1136

\citet{Grahametal:2004} results in values of $A = 6$~nm and $L =

1137

3480$~nm. This is to be compared with $A=14.5$~nm and $L=309$~nm,

1138

reported by \citet{Sunetal:2002}, for a {\em single} collagen

1139

molecule. Taken together, these results demonstrate that the WLC

1140

analysis correctly predicts a collagen fibril to be longer and more

1141

compliant than its constituent molecule due to compliant

1142

intermolecular cross-links in fibrils.

1155

1143

1156

1144

To model a collagen network structure, the WLC model has been embedded

1157

1145

as a single constituent chain of an eight-chain model

1167

1155

\begin{equation}

1168

1156

\begin{split}

1169

1157

\rho^\mathrm{c}_{0}\hat{e}^\mathrm{c}

1170

(\bF^{\mathrm{e}^\mathrm{c}},\rho^\mathrm{c}_{0}) &= \frac{N k

1171

\theta L}{4 A}\left(\frac{2 r^2}{L^2} + \frac{1}{1-r/L} -

1158

(\bF^{\mathrm{e}^\mathrm{c}},\rho^\mathrm{c}_{0}) &= \frac{N k \theta

1159

L}{4 A}\left(\frac{2 r^2}{L^2} + \frac{1}{1-r/L} -

1172

1160

\frac{r}{L}\right)\\ & +

1173

1161

\frac{\gamma}{\beta}({J^{\mathrm{e}^{\mathrm{c}^{-2\beta}}}} -1) +

1174

1162

\gamma{\bf 1}\colon(\bC^{\mathrm{e}^{\mathrm{c}}}-{\bf 1})\\ &-\frac{N

1206

1194

\sqrt{\bN_I\cdot\bC ^ {\mathrm{e}^{\mathrm{c}}}\bN_I}$, and $\bN_I,\,I

1207

1195

= 1,2,3$ are the unit vectors along the three unit cell axes,

1208

1196

respectively. Since the quantities $L$, $A$ and $r$ occur only as the

1209

ratios $L/A$ and $r/L$ in this model~(\ref{eight-chain-model}),

1197

ratios $L/A$ and $r/L$ in this model~(\ref{wlcm-energy}),

1210

1198

Table~\ref{parameters} lists the lengths $L,\ A,\ a,\ b,\ c$ used in

1211

1199

the computations in a non-dimensionalised manner.

1212

1200

1260

1248

\ref{saturation-and-tissue-swelling}, we have $\tilde{v}^\mathrm{f} +

1261

1249

\tilde{v}^\mathrm{c} = 1$ at $t = 0$, the saturation condition in

1262

1250

$\Omega_t$ when solutes are at low concentrations. At later times,

1263

Equation (\ref{saturation}) holds for

1264

$\bF^{\mathrm{g}^\mathrm{f}}$. If $\mathrm{det}(\bF(\bF ^

1265

{\mathrm{g}^\mathrm{f}})^{-1}) \le 1$ we set $\bF^{\mathrm{e} ^

1266

\mathrm{f}} = \bF(\bF^{\mathrm{g}^\mathrm{f}})^{-1}$ and

1267

$\bF^{\mathrm{v}} = {\bf 1}$ for no cavitation. Otherwise, since

1268

$\mathrm{det}(\bF(\bF^{\mathrm{g}^\mathrm{f}})^{-1}) > 1$, we specify

1269

$\bF^{\mathrm{e}^\mathrm{f}} = \mathrm{det}(\bF(\bF^{\mathrm{g} ^

1270

\mathrm{f}})^{-1})^{-1/3}\bF (\bF^{\mathrm{g}^\mathrm{f}})^{-1}$,

1271

thus restricting $\bF^{\mathrm{e}^\mathrm{f}}$ to be unimodular and

1272

allow cavitation by writing $\bF^{\mathrm{v}} =

1251

Equation (\ref{saturation}) holds for $\bF^{\mathrm{g}^\mathrm{f}}$.

1252

If $\mathrm{det}(\bF(\bF ^ {\mathrm{g}^\mathrm{f}})^{-1}) \le 1$ we

1253

set $\bF^{\mathrm{e} ^ \mathrm{f}} = \bF(\bF^{\mathrm{g}^\mathrm{f}})

1254

^ {-1}$ and $\bF^{\mathrm{v}} = {\bf 1}$ for no cavitation. Otherwise,

1255

since $\mathrm{det}(\bF(\bF^{\mathrm{g}^\mathrm{f}})^{-1}) > 1$, we

1256

specify $\bF^{\mathrm{e}^\mathrm{f}} = \mathrm{det}(\bF(\bF ^

1257

{\mathrm{g} ^ \mathrm{f}})^{-1})^{-1/3}\bF (\bF^{\mathrm{g} ^

1258

\mathrm{f}})^{-1}$, thus restricting $\bF^{\mathrm{e}^\mathrm{f}}$

1259

to be unimodular and allow cavitation by writing $\bF^{\mathrm{v}} =

1273

1260

\bF(\bF^{\mathrm{e}^\mathrm{f}} \bF^{\mathrm{g}^\mathrm{f}}) ^

1274

1261

{-1}$.\\ % Hack

1275

1262

1276

1277

1263

These conditional relations are summarised as:

1278

1264

1279

1265

\begin{equation}

1350

1336

potential.}

1351

1337

1352

1338

Experimentally determined transport coefficients (e.g. for mouse tail

1353

skin \citep{Swartzetal:99} and rabbit Achilles tendons

1354

\citep{Hanetal:2000}) are used for the fluid mobility values. Recall

1355

that the terms in the parenthesis on the right hand-side of

1356

\mbox{Equation (\ref{fluidflux})} sum to give the total driving force

1357

for transport. The first term is the contribution due to gravitational

1339

skin \citep{Swartzetal:99} and rabbit Achilles tendons \citep{Hanetal%

1340

:2000}) are used for the fluid mobility values. Recall that the

1341

terms in the parenthesis on the right hand-side of \mbox {Equation

1342

(\ref{fluidflux})} sum to give the total driving force for

1343

transport. The first term is the contribution due to gravitational

1358

1344

acceleration. In order to maintain physiological relevance, this term

1359

1345

has been neglected in the following treatment. The second term arises

1360

1346

from stress divergence. In the case of a non-uniform partial stress,

1516

1502

produced.

1517

1503

1518

1504

(\romannumeral 2) {\em Michaelis-Menten} enzyme kinetics (see, for

1519

e.g., \cite{Sengersetal:2004}), which involves a two-step reaction,

1505

e.g., \cite{Sengersetal:2004}), which involves a two-step reaction,

1520

1506

introduces collagen and solute source terms given by

1521

1507

1522

1508

\begin{equation}

1523

\Pi^\mathrm{s} =

1524

\frac{-(k_{\mathrm{max}}\rho^{\mathrm{s}})}

1525

{(\rho^{\mathrm{s}}_m+\rho^{\mathrm{s}})}

1526

\rho_{\mathrm{cell}}, \quad\Pi^\mathrm{c} = -\Pi^\mathrm{s},

1509

\Pi^\mathrm{s} = \frac{-(k_{\mathrm{max}}\rho^{\mathrm{s}})}

1510

{(\rho^{\mathrm{s}}_m+\rho^{\mathrm{s}})} \rho_{\mathrm{cell}},

1511

\quad\Pi^\mathrm{c} = -\Pi^\mathrm{s},

1527

1512

\label{enzyme-kinetics-source}

1528

1513

\end{equation}

1529

1514

1549

1534

mechanical influences on the strengths and stiffnesses of tendons by

1550

1535

comparing the stress-strain responses of unloaded control specimens

1551

1536

with those subjected to two different load cases (denoted {\em delay}

1552

and {\em no delay}) on the figure.

1537

and {\em no delay} in the figure).

1553

1538

1554

1539

An example of source terms of this form was originally proposed in the

1555

1540

context of bone growth \citep{HarriganHamilton:93}. I am not aware of

1620

1605

\begin{center}

1621

1606

\includegraphics[width=0.8\textwidth]{images/elucidation/%

1622

1607

concentration}

1623

\caption{Pore structure at the boundary deforming with the tissue.}

1608

\caption{Pore structure at the boundary deforming with the

1609

tissue.}

1624

1610

\label{current-conc-fluid-bc}

1625

1611

\end{center}

1626

1612

{If the pore structure at the boundary deforms with the tissue

1707

1693

1708

1694

From \mbox{Equation (\ref{current-mass-balance})}, the local form of

1709

1695

the balance of mass for the fluid species (assuming that the fluid

1710

species does not take part in reactions, i.e. $\Pi^\mathrm{f}=0$) in

1696

species does not take part in reactions, i.e. $\pi^\mathrm{f}=0$) in

1711

1697

the current configuration is

1712

1698

1713

1699

\begin{equation}

1714

\frac{\mathrm{d}\rho^f}{\mathrm{d}t} = - \mathrm{div}\left[\bm^f\right]

1715

- \rho^f \mathrm{div}\left[\bv\right].

1700

\frac{\mathrm{d}\rho^f}{\mathrm{d}t} = -

1701

\mathrm{div}\left[\bm^f\right] - \rho^f \mathrm{div}\left[\bv\right].

1716

1702

\label{fluidtransporteqn}

1717

1703

\end{equation}

1718

1704

1727

1713

\begin{equation}

1728

1714

\begin{split}

1729

1715

\rho_{0}^{\mathrm{f}}(\bX,0)

1730

&=:\rho_{0_{\mathrm{ini}}}^{\mathrm{f}}(\bX)\\

1731

&=\rho_{\mathrm{ini}}^{\mathrm{f}}(\bx\circ\Bvarphi)

1732

J(\bX, t)\\ &=\frac{\rho^{\mathrm{f}}

1733

(\bx\circ\Bvarphi,t)} {J^{f_\mathrm{g}}(\bX,t)}

1734

J(\bX,t)\\ &=\rho^{\mathrm{f}} (\bx\circ\Bvarphi,t)

1735

\cancelto{\approx 1\ \forall\ t}

1736

{J^{f_\mathrm{e}}}(\bX,t).\\

1716

&=:\rho_{0_{\mathrm{ini}}}^{\mathrm{f}}(\bX)\\ &=\rho_{\mathrm{ini}}

1717

^{\mathrm{f}} (\bx\circ\Bvarphi) J(\bX, t)\\ &=\frac{\rho^{\mathrm{f}}

1718

(\bx \circ \Bvarphi,t)} {J^{f_\mathrm{g}}(\bX,t)} J(\bX,t)

1719

\\ &=\rho^{\mathrm{f}} (\bx\circ\Bvarphi,t) \cancelto{\approx

1720

1\ \forall\ t} {J^{f_\mathrm{e}}}(\bX,t).\\

1737

1721

\label{incompderiv}

1738

1722

\end{split}

1739

1723

\end{equation}

1740

1724

1741

In (\ref{incompderiv}), $J := \mathrm{det}(\bF)$ and $J^{f_\mathrm{g}} :=

1742

\mathrm{det}(\bF^{\mathrm{g}^{\mathrm{f}}})$. The quantity

1743

$\rho_{\mathrm{ini}}^{\mathrm{f}}$ is defined by the right hand-sides

1744

of the first and second lines of (\ref{incompderiv}). To follow the

1725

In (\ref{incompderiv}), $J := \mathrm{det}(\bF)$ and $J^{f_\mathrm{g}}

1726

:= \mathrm{det}(\bF^{\mathrm{g}^{\mathrm{f}}})$. The quantity $\rho%

1727

_{\mathrm{ini}} ^ {\mathrm{f}}$ is defined by the right hand-sides of

1728

the first and second lines of (\ref{incompderiv}). To follow the

1745

1729

argument, consider, momentarily, a \emph{compressible} fluid. If the

1746

1730

current concentration, $\rho^\mathrm{f}$, changes due to elastic

1747

deformation of the fluid and by transport, then

1748

$\rho_{\mathrm{ini}}^{\mathrm{f}}$ as defined is not a

1749

physically-realised fluid concentration. It bears a purely

1750

mathematical relation to the current concentration, $\rho^\mathrm{f}$,

1751

since the latter quantity represents the effect of deformation of a

1752

tissue point as well as change in mass due to transport at that

1753

point.

1731

deformation of the fluid and by transport, then $\rho_{%

1732

\mathrm{ini}}^{\mathrm{f}}$ as defined is not a physically-realised

1733

fluid concentration. It bears a purely mathematical relation to the

1734

current concentration, $\rho^\mathrm{f}$, since the latter quantity

1735

represents the effect of deformation of a tissue point as well as

1736

change in mass due to transport at that point.

1754

1737

1755

1738

If the contribution due to mass change at a point is scaled out of

1756

1739

$\rho^\mathrm{f}$ the quotient is identical to the result of dividing

1757

1740

$\rho_{0_{\mathrm{ini}}}^{\mathrm{f}}$ by the deformation only. This

1758

1741

is expressed in the relation between the right hand-sides of the

1759

1742

second and third lines of (\ref{incompderiv}). The elastic component

1760

of fluid volume change in a pore is $J^{f_\mathrm{e}} :=

1761

\mathrm{det}(\bF^{\mathrm{e}^{\mathrm{f}}})$, which appears in the

1762

third line of (\ref{incompderiv}) via the preceding arguments. Clearly

1763

then, for a fluid demonstrating near incompressibility intrinsically

1764

(i.e., the true density is nearly constant), we have $J^{f_\mathrm{e}}

1765

\approx 1$ as indicated. Equation (\ref{incompderiv}) therefore shows

1766

that for a nearly incompressible fluid occupying the pores of a

1767

tissue, if we further assume that the pore structure deforms as the

1768

solid collagenous skeleton, $\rho_0^\mathrm{f}(\bX,0) \approx

1769

\rho^\mathrm{f}(\bx\circ\Bvarphi,t)$; i.e., the fluid concentration as

1743

of fluid volume change in a pore is $J^{f_\mathrm{e}} := \mathrm{det}

1744

(\bF^{\mathrm{e}^{\mathrm{f}}})$, which appears in the third line of

1745

(\ref{incompderiv}) via the preceding arguments. Clearly then, for a

1746

fluid demonstrating near incompressibility intrinsically (i.e., the

1747

true density is nearly constant), we have $J^{f_\mathrm{e}} \approx 1$

1748

as indicated. Equation (\ref{incompderiv}) therefore shows that for a

1749

nearly incompressible fluid occupying the pores of a tissue, if we

1750

further assume that the pore structure deforms as the solid

1751

collagenous skeleton, $\rho_0^\mathrm{f}(\bX,0) \approx \rho ^

1752

\mathrm{f} (\bx\circ\Bvarphi,t)$; i.e., the fluid concentration as

1770

1753

measured in the current configuration is approximately constant in

1771

1754

space and time. This allows us to write,

1772

1755

1773

1756

\vspace{-0.5cm} %Hack

1774

1757

1775

1758

\begin{equation}

1776

\frac{\partial} {\partial t}\left(

1777

\rho_{0_{\mathrm{ini}}}^{f}(\bX) \right) \equiv 0 \Rightarrow

1778

\frac{\partial} {\partial t}\left(\rho^{f} (\bx\circ\Bvarphi,t)

1779

\right)\Big\vert_{\bX} = 0,

1759

\frac{\partial} {\partial t}\left( \rho_{0_{\mathrm{ini}}}^{f}(\bX)

1760

\right) \equiv 0 \Rightarrow \frac{\partial} {\partial

1761

t}\left(\rho^{f} (\bx\circ\Bvarphi,t) \right)\Big\vert_{\bX} = 0,

1780

1762

\end{equation}

1781

1763

1782

1764

\noindent which is a hidden implication of our assumption of a

1808

1790

\begin{split}

1809

1791

\frac{\mathrm{d}\rho^\mathrm{s}}{\mathrm{d}t} &= \pi^\mathrm{s}-

1810

1792

\mathrm{div} \left[\widetilde{\bm^\mathrm{s}}+

1811

\frac{\rho^\mathrm{s}}{\rho^f}\bm^f\right] - \rho^\mathrm{s}

1812

\mathrm{div}[\bv]\\ &= \frac{\rho^\mathrm{s}}{\rho^f}\left(\cancelto{0}

1813

{-\mathrm{div}\left[\rho^f \bv^{f}\right] - \rho^f

1814

\mathrm{div}[\bv]}\right)\\ &\quad{} + \pi^\mathrm{s} -

1815

\mathrm{div}\left[\widetilde{\bm^\mathrm{s}}\right]

1816

-\bm^f\cdot\mathrm{grad}\left[\frac{ \rho^\mathrm{s}}{\rho^f}\right].\\

1793

\frac{\rho^\mathrm{s}}{\rho^f}\bm^f\right] - \rho^\mathrm{s}

1794

\mathrm{div}[\bv]\\ &=

1795

\frac{\rho^\mathrm{s}}{\rho^f}\left(\cancelto{0}

1796

{-\mathrm{div}\left[\rho^f \bv^{f}\right] - \rho^f

1797

\mathrm{div}[\bv]}\right)\\ &\quad{} + \pi^\mathrm{s} -

1798

\mathrm{div}\left[\widetilde{\bm^\mathrm{s}}\right]

1799

-\bm^f\cdot\mathrm{grad}\left[\frac{

1800

\rho^\mathrm{s}}{\rho^f}\right].\\

1817

1801

\end{split}

1818

1802

\end{equation}

1819

1803

1826

1810

\mathrm{div}\left[\widetilde{\bm^\mathrm{s}}\right] -

1827

1811

\frac{\bm^f\cdot\mathrm{grad}\left[\rho^\mathrm{s}\right]}{\rho^f} +

1828

1812

\frac{\rho^\mathrm{s} \bm^f \cdot \mathrm{grad}\left[\rho^f\right]}

1829

{\rho^{f^2}}.

1813

{\rho^{f^2}}.

1830

1814

\label{stdform}

1831

1815

\end{equation}

1832

1816

1866

1850

diffusivity, $\bm^{f}/\rho^{f}$ is the advective velocity, and

1867

1851

$\pi^\mathrm{s}$ is the volumetric source term. This form is well

1868

1852

suited for stabilisation schemes such as the streamline upwind

1869

Petrov-Galerkin (SUPG)

1870

method\footnote{Appendix~\ref{stabilisation-solute-transport}

1871

provides, in weak form, the SUPG-stabilised method for

1872

Equation~(\ref{morestdform}).} (see, for e.g., \cite{Paper6}), which

1873

limit spatial oscillations otherwise observed when the element

1874

P\'eclet number is large. Figure~\ref{stable-solution} shows the

1875

SUPG-stabilised solution for the simple advection-diffusion problem

1876

considered previously at an identical P\'eclet number.

1853

Petrov-Galerkin (SUPG) method\footnote{Appendix~\ref{stabilisation%

1854

-solute-transport} provides, in weak form, the SUPG-stabilised

1855

method for Equation~(\ref{morestdform}).} (see, for e.g.,

1856

\cite{Paper6}), which limit spatial oscillations otherwise observed

1857

when the element P\'eclet number is large. Figure~\ref{stable-solu%

1858

tion} shows the SUPG-stabilised solution for the simple advection%

1859

-diffusion problem considered previously at an identical P\'eclet

1860

number.

1877

1861

1878

1862

%

1879

1863