Coupling chemotaxis and growth poromechanics for the modelling of feather primordia patterning Nicol as A. Barnaﬁ Luis Miguel De Oliveira Vilaca Michel C. Milinkovitch

2025-05-06 0 0 3.63MB 25 页 10玖币

侵权投诉

Coupling chemotaxis and growth poromechanics

for the modelling of feather primordia patterning

Nicol´

as A. Barnaﬁ*

, Luis Miguel De Oliveira Vilaca†

, Michel C. Milinkovitch†

Ricardo Ruiz-Baier‡

October 18, 2022

Abstract

We propose a new mathematical model for the interaction of skin cell populations with ﬁbroblast

growth factor and bone morphogenetic protein, occurring within deformable porous media. The equations

for feather primordia pattering are based on the work by K.J. Painter et al. [J. Theoret. Biol., 437 (2018) 225–

238]. We perform a linear stability analysis to identify relevant parameters in the coupling mechanisms,

focusing in the regime of inﬁnitesimal strains. We also extend the model to the case of nonlinear poroelas-

ticity and include solid growth by means of Lee decompositions of the deformation gradient. We present

a few illustrative computational examples in 2D and 3D, and brieﬂy discuss the design of tailored efﬁcient

solvers.

1 Introduction

Chemotaxis models [24] describe the directed movement of cells in response to chemicals (attractants or

repellents), and can predict the formation of clustered structures. This process has been observed in a

variety of embryogenesis processes, such as gastrulation [48] and feather development [27]. In the latter, we

refer to the buds that then give origin to feathers as primordia, and their origin has been modelled through

Keller–Segel type chemotaxis models considering ﬁbroblast growth factor (FGF) and bone morphogenic

protein (BMP) [33, 35, 36].

Poroelasticity is a mixture model in which a solid phase coexists with (at least) one ﬂuid phase [12].

Biological organs and tissues are naturally porous at the tissue level, as they are composed, for example, of

both muscle and blood. This phase separation persists even up to the cellular level, as there is the cytoskele-

ton and the cytoplasm. For this reason, poroelastic models have become very pervasive in the modelling

of soft living tissue, such as in oedema formation [3, 28], cardiac perfusion [47, 2], lung characterisation [7],

and brain injury [46].

In the context of biologically-oriented problems, experiments have shown that the rheology of cyto-

plasm within living cells exhibits a poroelastic behaviour [31], and in turn, the composition of cells and the

extracellular matrix constitutes an overall poromechanical system. The presence of chemical solutes locally

*Department of Mathematics, Universit`

a di Pavia, Via Ferrata 1, 27100 Pavia, Italy. E-mail: nicolas.barnafi@unipv.it

†Laboratory of Artiﬁcial & Natural Evolution (LANE), Department of Genetics and Evolution, University of Geneva,

4 Boulevard d’Yvoy, 1205 Geneva, Switzerland; and SIB Swiss Institute of Bioinformatics, Geneva, Switzerland. E-mail:

LuisMiguel.DeOliveiraVilaca@unige.ch, Michel.Milinkovitch@unige.ch

‡School of Mathematics, Monash University, 9 Rainforest Walk, Melbourne 3800 VIC, Australia; and Universidad Adventista de

Chile, Casilla 7-D, Chill´

an, Chile. E-mail: ricardo.ruizbaier@monash.edu

arXiv:2210.08308v1 [math.NA] 15 Oct 2022

modiﬁes morphoelastic properties and these processes can be homogenised to obtain macroscopic models

of poroelasticity coupled with advection-reaction-diffusion equations (see e.g. [11, 37]). With this biological

basis, the scope of our work is twofold: on one hand, we extend the existing pattern formation models for

primordia by considering their interaction with the intracellular space in the outset of growth. On the other

hand, we perform a thorough stability analysis to understand the coupling mechanisms in the model and

to investigate the conditions that give rise to pattern formation.

We have structured the remainder of this paper in the following manner. Section 2 describes the coupled

model for poro-mechano-chemical interactions, restricting the presentation to the regime of linear poroe-

lasticity. We give an adimensional form of the governing equations, making precise boundary and initial

conditions. In Section 3 we perform a linear stability analysis addressing pattern formation according to

Turing instabilities. We separate the discussion in some relevant cases and derive and portray patterning

spaces. Section 4 is devoted to extending the model to the case of nonlinear (ﬁnite-strain) poroelasticity and

material growth, stating also the coupling with chemotaxis in the undeformed conﬁguration. Some numer-

ical examples are given in Section 5 and we close with a summary and discussion of model extensions in

Section 6.

2 A coupled model of linear poroelasticity and chemotaxis

Let us consider a piece of soft material as a porous medium in Rd,d=2, 3, composed by a mixture of

incompressible grains and interstitial ﬂuid, whose description can be placed in the context of the classical

Biot consolidation problem (see e.g. [42]). In the absence of gravitational forces, of body loads, and of mass

sources or sinks, we seek for each time t∈(0, tﬁnal], the displacement of the porous skeleton, u(t):Ω→Rd,

and the pore pressure of the ﬂuid, p(t):Ω→R, such that

∂tC0p+αBW div u−1

ηdiv{κ∇p}=0 in Ω×(0, tﬁnal], (2.1a)

σ=σporoelast +σact in Ω×(0, tﬁnal], (2.1b)

ρ∂ttu−div σ=0in Ω×(0, tﬁnal], (2.1c)

where κ(x)is the hydraulic conductivity of the porous medium, ρis the density of the solid material, ηis

the constant viscosity of the interstitial ﬂuid, C0is the constrained speciﬁc storage coefﬁcient, αBW is the

Biot-Willis consolidation parameter, In (2.1b) we are supposing that the poromechanical deformations are

also actively inﬂuenced by microscopic tension generation. A very simple description is given in terms of

active stresses: we assume that the total Cauchy stress contains a passive and an active component, where

σporoelast =λ(div u)I+2µε(u)−pI, (2.2)

and σact is speciﬁed in (2.5), below. The tensor ε(u) = 1

2(∇u+∇u|)is that of inﬁnitesimal strains, I

denotes the second-order identity tensor, and µ,λare the Lam´

e constants (shear and dilation moduli) of the

solid structure. Equations (2.1a)-(2.1c) represent the conservation of mass, the constitutive relation, and the

conservation of linear momentum, respectively.

In addition, let us consider a modiﬁed Patlak–Keller–Segel model for the distribution of chemotactic

cell populations of mesenchymal cells, m, epithelium activation state, e, ﬁbroblast growth factor (FGF),

f, and bone morphogenetic protein (BMP), b. The base-line model has been developed in [36] and (after

being properly modiﬁed to account for the motion of the underlying deformable porous media) it can be

summarised as follows

∂tm+∂tu·∇m−divDm∇m−αmexp(−γm)∇f=0 in Ω×(0, tﬁnal], (2.3a)

∂te−[κ1w(x,t)h1(m) + κ2h2(m)](1−e) + [1−h1(m)](κ3+κ4b)e=0 in Ω×(0, tﬁnal], (2.3b)

∂tf+∂tu·∇f−div(Df∇f)−κFe+δFf+ξfdiv u=0 in Ω×(0, tﬁnal], (2.3c)

∂tb+∂tu·∇b−div(Db∇b)−κBh3(m)m+δBb=0 in Ω×(0, tﬁnal], (2.3d)

where Dm,Df,Dbare positive deﬁnite diffusion matrices, and the spatio-temporal and nonlinear coefﬁ-

cients (in this case, the priming wave responsible for the generation of spatial patterns and the degree of

clustering of mesenchymal cells, respectively) are deﬁned as

w(x,t) = ω1

2{1+tanh(ω2[t−x2/ω3])},hi(m) = mPi[KPi

i+mPi]−1, (2.4)

where κi,ωi,Pi,γ,δi,ξfare positive model constants. Note that the mechano-chemical feedback (the process

where mechanical forces modify the reaction-diffusion effects) is here assumed only through an additional

reaction term in the FGF equation (2.3c), depending linearly on volume change. Since ξf>0, this contri-

bution acts as a local sink of FGF during dilation.

On the other hand, we assume that the active stress component acts isotropically on the medium (see

e.g. [22]), and it depends nonlinearly on the concentration of mesenchymal cells, as proposed for instance

in [34]

σact =λ+2µ

3τm

1+ζm2I, (2.5)

with τa model constant to be speciﬁed later on, which can be positive (implying that the function modiﬁes

the motion of the medium as a local dilation) or negative (isotropically distributed compression).

In order to reduce the number of model parameters, we focus on a dimensionless counterpart of systems

(2.1a)-(2.1c) and (2.3a)-(2.3d), which can be derived using the following transformation, suggested in [36]

m=m∗

γ,e=e∗,f=κFf∗

δB

,b=κBb∗

γ,

u=sDb

δB

u∗,p=p∗,t=t∗

δB

,x=sDb

δB

x∗,

where hereafter the stars are dropped for notational convenience. The adimensional coupled system is then

equipped with appropriate initial data at rest

m(0) = m0,e(0) = e0=κ1w(x, 0)h1(m0) + κ2h2(m0)

κ1w(x, 0)h1(m0) + κ2h2(m0) + [1−h1(m0)](κ3+κ4b),

f(0) = e0

δF

,b(0) = h3(m0)m0,u(0) = 0,∂tu(0) = 0,p(0) = 0, (2.6)

deﬁned in Ω; and boundary conditions in the following manner

{Dm∇m−αmexp(−γm)∇f} ·n=Df∇f·n=∇b·n=0 on ∂Ω×(0, tﬁnal], (2.7a)

u=uΓand κ

η∇p·n=0 on Γ×(0, tﬁnal], (2.7b)

σn=tand p=p0on Σ×(0, tﬁnal], (2.7c)

where ndenotes the outward unit normal on the boundary, and ∂Ω=Γ∪Σis disjointly split into Γand Σ

where we prescribe clamped boundaries and zero ﬂuid normal ﬂuxes; and a given traction ttogether with

constant ﬂuid pressure p0, respectively.

3 Linear stability analysis and dispersion relation

3.1 Preliminaries.

We proceed to derive a linear stability analysis for the coupled problem (2.1a)-(2.5). As usual the analysis

is performed on an inﬁnite domain in Rd, with d=2, 3. The ﬁrst step consists in linearising the poro-

mechano-chemical system around a steady state, deﬁned in (2.6). The linearised dimensionless equations

are given by

∂tC0p+αBW div u−1

ηdiv{κ∇p}=0,

ρ∂ttu−divσporoelast +σlin

act=0,

∂tm−div{Dm∇m−αm0e−m0∇f}=0,

∂te−Am(m0,e0)m+Ae(m0,b0)e+Ab(m0,e0)b=0,

∂tf−div{Df∇f} −e+δFf+ξfdiv u=0,

∂tb−∆b−H3(m0)m+δBb=0,

(3.1)

where

σlin

act = λ+2µ

3τ1−ζm2

1+ζm2

02!mI,

and Am,Ae,Ab,H3are functions of the steady state values, which will be made precise in Deﬁnition

3.1, below. We look then for solutions of the form u,p,m,e,f,b∝eik·x+φt, where kis the wave vector (a

measure of spatial structure) and φis the linear growth factor. By substituting this ansatz on the system

(3.1), we get a system of linear equations for the vector w=(u,p,m,e,f,b)|, where the associated complex

eigenvalues φ, give information on occurrence of instability of the steady state (i.e., pattern formation),

when its real component is positive.1In order to specify such system, we collect in Proposition 3.1 some

useful preliminary relations. The proof is postponed to the Appendix.

Proposition 3.1 Let g,vbe sufﬁciently regular scalar and vector functions deﬁned by g =exp(ik·x+φt)and

v=v0exp(ik·x+φt)(with v0, a constant vector), respectively. Let us also set θ=div v. Then

∇f=ikf,∆f=i2fk·k=−k2f,∂tf=φf,∇v=iv⊗k, div v=iv·k,∂tv=φv, (3.2a)

div ε(v) = −(v·k)k+k2v

2,∂tε(v) = φε(v),div(θI) = −(v·k)k,∂tθ=φθ. (3.2b)

The sought system is deﬁned next, starting from (3.1).

Deﬁnition 3.1 Let w=(u,p,m,e,f,b)|∈Rd+5, d =2, 3 be the vector of independent variables. Then the

associated linear system is given by Mw=0d+5, where the matrix M(d+5)×(d+5)adopts the form

M=M11 M12

M21 M22,

with the blocks deﬁned as

M11 =





A1k1+B A1k2··· A1kdik1

A2k1A2k2+B··· A2kdik2

.....

Adk1A1k2··· Adkd+Bikd

iαBW φk1iαBW φk2··· iαBW φkdC







1We stress that the present analysis will only address Turing patterning (Hopf bifurcations or other forms of instability that relate

to analysing the imaginary part of φ, are not considered).

M12 =





−iλ+2µ

3τ1−ζm2

(1+ζm2

0)2k1000

−iλ+2µ

3τ1−ζm2

(1+ζm2

0)2kd000





,

M21 =





0··· 0

iξfk1··· iξfkd

0··· 0





,

M22 =





φ+Dmk20−αm0e−m0k20

−Amφ+Ae0Ab

0−1φ+Dfk2+δF0

−H30 0 φ+k2+1





,

and where the relevant entries are deﬁned as B =ρφ2+µk2, C =C0φ+κ

ηk2, and

Ai=(µ+λ)ki,for all i in {1, . . . , d},

Am=κ1ω1h0

1(m0) + κ2h0

2(m0)(1−e0) + h0

1(m0)(κ3+κ4b0)e0,

Ae=κ1ω1h1(m0) + κ2h2(m0) + (1−h1(m0))(κ3+κ4b0),

Ab= (1−h1(m0))κ4e0,H3=h3(m0) + h0

3(m0)m0.

We then proceed to obtain a dispersion relation associated with the characteristic polynomial P(φ) =

det(M)of the matrix described in Deﬁnition 3.1. We obtain

P(φ;k2) = B(φ;k2)d−1P1(φ;k2) + P2(φ;k2)P3(φ;k2), (3.3)

where B(φ;k2) = ρφ2+µk2is a polynomial with pure imaginary roots. Consequently, it does not have an

inﬂuence on the stability of the steady state of system (2.1a)-(2.5). The polynomials Pi,i=1, 2, 3, are given

in what follows.

Deﬁnition 3.2 The polynomials conforming the characteristic equation (3.3) are

P1(φ;k2) = b3φ3+b2φ2+b1φ+b0,

P2(φ;k2) = a4φ4+a3φ3+a2φ2+a1φ+a0, (3.4)

P3(φ;k2) = c3φ3+c2φ2+c1φ+c0,

where all coefﬁcients in (3.4) adopt the following forms

a0=DmDf(kon +koff)k6+(DmδB+DmDf)(kon +koff)−αm0e−m0Amk4

+DmδB(kon +koff) + αm0e−m0H3Ab−αm0e−m0Amk2,

a1=DmDfk6+(DmDf+Dm+Df)(kon +koff) + DmδB+DmD f k4

+(DmδB+δB+Dm+Df)(kon +koff) + DmδB−αm0e−m0Amk2+δB(kon +koff),

a2=DmDf+Dm+Dfk4+DmδB+ (Dm+Df+1)(kon +koff) + δ+Dm+Dfk2

+ (δB+1)(kon +koff) + δB,

a3=Dm+Df+1k2+δB+kon +koff +1, a4=1,

文档加载中……请稍候！
如果长时间未打开，您也可以点击刷新试试。

下载文档到电脑，查找使用更方便

10 玖币 0人已下载

立即下载

摘要：

CouplingchemotaxisandgrowthporomechanicsforthemodellingoffeatherprimordiapatterningNicol´asA.Barna*,LuisMiguelDeOliveiraVilaca,MichelC.MilinkovitchRicardoRuiz-BaierOctober18,2022AbstractWeproposeanewmathematicalmodelfortheinteractionofskincellpopulationswithbroblastgrowthfactorandbonemorphogene...

展开>> 收起<<

Coupling chemotaxis and growth poromechanics for the modelling of feather primordia patterning Nicol as A. Barnaﬁ Luis Miguel De Oliveira Vilaca Michel C. Milinkovitch.pdf

共25页,预览5页

还剩页未读，继续阅读

声明：本站为文档C2C交易模式，即用户上传的文档直接被用户下载，本站只是中间服务平台，本站所有文档下载所得的收益归上传人(含作者)所有。玖贝云文库仅提供信息存储空间，仅对用户上传内容的表现方式做保护处理，对上载内容本身不做任何修改或编辑。若文档所含内容侵犯了您的版权或隐私，请立即通知玖贝云文库，我们立即给予删除！

Coupling chemotaxis and growth poromechanics for the modelling of feather primordia patterning Nicol as A. Barnaﬁ Luis Miguel De Oliveira Vilaca Michel C. Milinkovitch

相关推荐

开通VIP享超值会员特权

作者详情

相关内容

热门标签

举报选择: