# The Avascular Tumour Growth in the Presence of Inhomogeneous Physical Parameters Imposed from a Finite Spherical Nutritive Environment

Authors: Foteini Kariotou; Panayiotis Vafeas
Journal: International Journal of Differential Equations
 (2012)

---

## Abstract

A well-known mathematical model of radially symmetric tumour growth is revisited in the present work. Under this aim, a cancerous spherical mass lying in a finite concentric nutritive surrounding is considered. The host spherical shell provides the tumor with vital nutrients, receives the debris of the necrotic cancer cells, and also transmits to the tumour the pressure imposed on its exterior boundary. We focus on studying the type of inhomogeneity that the nutrient supply and the pressure field imposed on the host exterior boundary, can exhibit in order for the spherical structure to be supported. It turns out that, if the imposed fields depart from being homogeneous, only a special type of interrelated inhomogeneity between nutrient and pressure can secure the spherical growth. The work includes an analytic derivation of the related boundary value problems based on physical conservation laws and their analytical treatment. Implementations in cases of special physical interest are examined, and also existing homogeneous results from the literature are fully recovered.


---

## Body

## 1. Introduction

Mathematical modelling of cancer tumour growth helps in understanding the mechanisms underlying the phenomenon in many ways. The main idea is that modelling physical hypotheses in mathematical terms, a biological phenomenon, such as tumour growth, is approached by a mathematical problem. Analytic and numerical or hybrid methods, applied to the problem studied, lead to conclusions on the solution of the problem. The conclusions are interpreted in biological terms that are subject to evaluation with respect to experimental evidence. Thus, mathematical modelling may potentially offer new perspective in the research directions of biological procedures.As cancer research and the related technology improved, an increasing amount of experimental data was produced, which made it hard to interpret. At the same time, mathematical models begun to develop at the aid of understanding the crucial parameters involved in tumour’s evolution [1, 2]. In the early years, mathematical models focused mainly on the avascular phase of the tumour’s development. This phase corresponds to the stages right after tumorigenesis and ends with a steady state, where the tumour’s volume gain, due to the new cancerous cells birth, balance the volume loss from the cells’ death and disintegration. In order for the tumour’s evolution to proceed, new phenomena, as angiogenesis, have to take place and the vascular phase may begin. In this phase, the tumour has developed a vascular net around it that provides the tumour cells with limitless nutrient supply and also it permits metastasis. As the present work focuses on the avascular phase on tumour evolution, we avoid presenting details on the proceeding phases, which can be found in many references, see, for example, [1, 2]. Since these next phases include the phenomena, appeared in the avascular phase and other more complicated ones, in depth studying of the avascular phase is well justified.Both deterministic, numerical and hybrid models have been developed towards the avascular growth study [2]. In the present work, we study a deterministic mathematical model that considers the tumour colony as a continuous medium, which evolves according to mass conservation law, Fick’s diffusion law, and fluid mechanics principles [3]. Deterministic models in tumour growth research have been studied since 1966 when Burton [4] modelled the tumour’s evolution with a Gompertzian curve. Though, the basis for mathematical modelling of avascular tumour growth has been established by the Greenspan works [5, 6] in the 1970s. Greenspan formulated the main procedures of the phenomenon, following existing experimental results, solving boundary value problems with reaction-diffusion partial differential equations in spherical domains. The main parameters that he investigated were the nutrient concentration in the tumour and its surrounding, as well as the concentration of a growth inhibitor factor. In addition, he introduced the idea of a surface tension that secures the solid structure of the tumour. Since then, his ideas were used widely and produced sophisticated models that investigate different aspects of avascular tumour growth. Indicatively, we refer to Adam’s papers [7, 8], to the series of works by Byrne and Chaplain [9–12], by Ward and King [13, 14], as well as by Mc Elwain and his collaborators [15, 16]. One can have a clear view on relevant references and analytic presentation of each contribution in the extensive reviews [1, 2, 11].In most works, a spherical tumour structure has been used and different forms of the physical parameters involved have been investigated. Heterogeneity with respect to different cell types has been taken into account [11], and the effect of inhomogeneity of the consumption rates and proliferation rates within the tumour colony has been investigated [15]. However, most works consider homogeneous nutrient supply from an infinite environment whereas this assumption ignores the realistic fact that tumours grow in finite organs with nonhomogeneous oxygen or glucose concentrations in the tumour’s microenvironment. Moreover, there is experimental evidence [17] that a geometrical confinement has significant impact on tumour growth, and to our knowledge there are only  few models that take this fact into account and stress out the importance of the tumour-host boundary interactions and the environmental heterogeneity [18, 19]. In the present work, we investigate the way that some key growth parameters interact, so as to permit a spherically symmetric tumour growth. The parameters that we consider are the finite growth environment, the nonhomogeneous pressure that it might induce, and the nonhomogeneous nutrient supply. Although a work on a spherical model is quite trivial after so much previous work on that geometric configuration, the point of view in our paper is to have a close look upon the detailed conditions that allow for such a configuration to occur and so to get the perspective of the conditions that might allow for other, less symmetric, geometries to appear, like the ellipsoidal one. Therefore, the present work should be viewed as a theoretical analysis of the basic mathematical model that provides proof that under the common assumptions the symmetric growth is inevitable and that radial growth is also possible under certain diversions from such assumptions. Moreover, concerning the cells’ motility that results to the tumour’s spatial expansion, we assume the rather logical approach that, additionally to the nutrient gradient and to the pressure gradient [6, 19], we further attribute the cells’ movement to the inhibitor gradient, a parameter that shows up to have a qualitative impact to the final results.In particular, in this work, we consider a spherical fully developed avascular tumour, meaning that it consists of a necrotic core, a quiescent layer, and a proliferating layer, which are under dynamical growth until the whole structure reaches a steady state. Then, the tumour may stay in a dormant state or it may further develop only by entering in the vascular phase, a case that is not examined in the present study. We further assume that the tumour grows in a finite concentric spherical host environment, modelling the healthy part of the organ in the interior of which the tumour develops. Both tissues, cancerous and healthy, are assumed to be incompressible fluids of different phase. Here, we need to clarify that the tumour’s avascular phase ends at a maximum tumour size, which is significantly smaller than a typical human host organ. Nevertheless, our study aims to reveal rather qualitative aspects than quantitative effects of the tumour growth. One such aspect is concerned with the confinement effects imposed by the growth environment. The effect of the stress, imposed by the surrounding, has been subject of investigation in several works [1, 8, 17, 18]. For example, the effects of a spatially restricted growth have been studied experimentally in [17], where a cancer cell culture was grown inside a very thin cylindrical tube. It turned out that the tube altered the shape and growth dynamics of the cell colony. Similar result was the outcome of the computational work [19]. Therefore, an analytical study that further investigates such dynamic growth inside a finite host environment is justified.We assume that the tumour receives the nutrient (in general oxygen or glycose) from the surrounding spherical layer with inhomogeneous concentration and also it is affected by a tissue oriented nonhomogeneous pressure field. The proof that either no kind of inhomogeneity or a special kind of inhomogeneity that relates the externally supplied nutrient field with the externally imposed pressure field can be supported by such model in order for its concentric spherical structure to be secured is one of the outcomes of the present work. In the first case, it turns out that a spherical tumour can be developed only under homogeneous nutrient supply and homogeneous pressure field imposed, which confirms all previous works in the field. Nevertheless, the model reveals another option that allows nonhomogeneity in a concentric spherical growth, provided that a special relation holds between the nutrient and the pressure fields, with respect to the cancer cells’ motility parameters. In any case, the concentration profile of both the nutrient and the inhibitor turns out to be strongly affected by the host tissue boundary. Similarly, the pressure distribution throughout the cancer tissue and the host tissue is affected by the host tissue boundary, but the evolution of the tumour boundary is independent of it.Before we proceed, we present briefly the basic assumptions we make within this work following the framework of [5, 6, 9, 12].We consider a spherical shell with radiusre, occupied by incompressible noncancerous tissue, which covers a growing spherical tumour colony. We assume that the tumour is avascular and fully developed, while it consists of three concentric spherical domains. The exterior spherical shell Ωp is occupied by proliferating cells, the subsequent shell Ωq is occupied by quiescent cells, and the interior necrotic sphere Ωn includes the dead cells and debris. The crucial physical parameters that govern the tumour growth at every fixed point r(t) are a nutritive substance supplied externally (usually considered as oxygen or glucose) with concentration σ(r), a growth inhibitor factor with concentration β(r), which is produced by the necrotic cell debris and the pressure field P(r) generated by the incompressible finite environment and the growing mass in its interior. The nutrient is supplied into the host tissue with constant concentration σ∞(re) and diffuses inwardly with diffusion constant kσ. Tumour cells consume nutrient and proliferate, as long as the nutrient concentration is maintained above the critical value σp*, which is fulfilled in the exterior shell Ωp. In the inner shell Ωq, where σ(r) lies in the interval 0<σn*<σ<σp*, the tumour cells cannot proliferate and live in a quiescent state. Finally, in Ωn, the concentration falls below σn* and the tumour cells die out of starvation. Moreover, Ωn also includes the dead cells due to apoptosis or programmed cell death, which are forwarded there by pressure gradients. The disintegration of the dead cells into simpler compounds produces a growth inhibitor substance in Ωn, which diffuses outwardly with diffusion constant kβ. This inhibitor factor prohibits cell proliferation, unless its concentration falls below a critical value βp*, which is assumed to hold in Ωp. As both the nutrient and the inhibitor diffuse inward to or outward from the tumour, respectively, their distribution is described with parabolic partial differential equations, as we will show in the sequel. Though, it is a common assumption [1, 3, 9] that the diffusion time scaling is significantly shorter than the time scaling one is really interested in, which is the tumour volume doubling time scale. Hence, as we demonstrate in Section 2.1.1, the elimination of the time derivative is justified concerning the nutrient and the inhibitor concentration, which are then considered to be in diffusive equilibrium state. The pressure field is generated from the pressure P∞(re), which is imposed by the host boundary and affects the interior pressure following the Young-Laplace law [5, 6] for two face incompressible fluids, together with the pressure imposed by cell proliferation and loss of necrotic mass, due to disintegration. The produced pressure differentials cause passive cell motion opposite to the direction of the pressure gradient. Additionally to the widely used Darcy’s law, which connects the velocity v(r) of the tumour cells at the exterior tumour surface with the pressure gradient [3, 19], we moreover assume that exterior tumour cells actively move towards the nutrient gradient [18, 20] and also we suggest that they move opposite to the growth inhibitor gradient. This assumption leads to the evolution equation of the tumour’s exterior boundary, which is a highly nonlinear ordinary differential equation and its solution is the goal of this work, the exterior tumour boundary as a function of time. Finally, we assume that, as the tumour grows, it maintains its concentric spherical structure, so all the interfaces of the model are time-dependent functions, connected between each other and the whole model is formulated as a free boundary problem.In Section2, we derive the corresponding partial differential equations that describe the nutrient concentration and the inhibitor concentration field, and the pressure field in each region, applying appropriate physical laws. Also, the corresponding boundary conditions are derived in each interface. Finally, the corresponding boundary value problems are stated in the spherical coordinate system. The solution of these boundary value problems is included in Section 3, where the restrictions on the kind of the inhomogeneity that can be supported with the spherical structure are derived. Also, the connection between the inner boundaries and the evolution equation of the exterior tumour boundary is provided. In Section 4, we demonstrate our results by presenting particular implementations regarding the exterior supply of the nutrient and the pressure imposed to the exterior boundary of the host surrounding of the tumour. Also, reduction of the results obtained in the present work to already existing corresponding expressions included in the original works [5, 6, 11] is provided in Section 4. Finally, an outline of our work follows in Section 5.


## 2. Statement of the Problem

Let(r,θ,ϕ) be the spherical coordinates of the point r∈ℝ3. Then, each distinct region in the tumour is defined as follows. The necrotic core is defined as
(2.1)Ωn={(r,θ,ϕ)∈ℝ3:0≤r<rn,0≤θ≤π,0≤ϕ<2π},
while the quiescent layer, the proliferating layer, and the exterior layer occupied by the noncancerous host tissue are defined as
(2.2)Ωq={(r,θ,ϕ)∈ℝ3:rn<r<rq,0≤θ≤π,0≤ϕ<2π},Ωp={(r,θ,ϕ)∈ℝ3:rq<r<rp,0≤θ≤π,0≤ϕ<2π},Ωe={(r,θ,ϕ)∈ℝ3:rp<r<re,0≤θ≤π,0≤ϕ<2π},
respectively. Each spherical region is bounded externally by the spherical surface Si:r=ri with i=n,q,p,e.Figure1 shows indicatively, with no reference to their real relative size, the aforementioned spherical regions under consideration along with their corresponding boundaries. In reality, the necrotic core is the dominant region of the tumour, while the exterior host layer is much larger than the included tumour regions.Figure 1
The spherical domains of the tumour growth model.


### 2.1. The Nutrient and the Inhibitor Concentration Problem

#### 2.1.1. Derivation of the Partial Differential Equations

According to the mass conservation law in each of the regionsΩi, i=n,q,p,e, the rate of change of the nutrient concentration in each region equals to the normal flow of nutrients Fin(r) inwards to Ωi, minus the flow outwards, Fout(r), and minus the amount of nutrient that is consumed in the interior of Ωi, that is,
(2.3)∫Ωi∂σi(r)∂tdv=∫∂Ωi(-n^)·Fin(r)ds-∫∂Ωin^·Fout(r)ds-∫ΩiΓidv,
where ∂Ωi is the boundary of the region Ωi, n^ is the outward unit normal vector to ∂Ωi, and Γi is the nutrient consumption rate, considered to be constant in Ωi. Since the nutrient flow is directed monotonically inwards in all spherical layers, the flow inwards takes place in the exterior boundary surface, denoted by Si+, and the flow outwards is realized from the inner boundary Si-. Considering also that for the spherical surface the normal direction is the radial one, (2.3) becomes
(2.4)∫Ωi∂σi(r)∂tdv=∫Si+(-r^)·Fin(r)ds-∫Si-(-r^)·Fout(r)ds-∫ΩiΓidvwithi=n,q,p,e.Assuming that the nutrient flow is only diffusive, it follows Fick’s law, which means that it is directed to regions of lower nutrient concentrationF(r)=-kσ∇σ(r). Then, relationship (2.4) yields
(2.5)∫Ωi∂σi(r)∂tdv=∫Si+r^·(kσ∇σi(r))ds-∫Si-r^·(kσ∇σi(r))ds-∫ΩiΓidv.Applying the divergence theorem inΩi in the right-hand side of (2.5), (2.5) becomes ∫Ωi(((∂σi(r))/∂t)-kσΔσi(r)+Γi)dv=0, which results in
(2.6)∂σi(r)∂t-kσΔσi(r)+Γi=0forr∈Ωiwithi=n,q,p,e.In particular, considering the physical assumptions made in the Introduction, the consumption rate in each region is defined asΓe=Γn=0, in Ωe and Ωn, respectively, where no consumption of nutrients takes place. In Ωq and Ωp, the consumption rate is assumed to depend on the vital state of the nutrients, which depends on the amount of the available nutrient, while it is defined as Γq=Γ(σn*/σ∞)and Γp=Γ(σp*/σ∞), respectively, where Γ is a proportional constant depending on the species of the host tissue cells, σn*, σp* are the nutrient critical values mentioned in the introduction, and σ∞ is a reference value corresponding to the nutrient concentration in the vicinity of the tumour’s host organ.We turn now to the inhibitor concentration. From the mathematical point of view, the inhibitor concentration problem is very similar to the nutrient concentration problem. The difference lies on the interpretation of the parameters and in their particular expressions. As mentioned in the introduction, the inhibitor is produced inside the necrotic core with production ratePn and is consumed in the quiescent and proliferating layer with consumption rate Ci with i=q,p, respectively. Therefore, the inhibitor flow is directed outwards and the corresponding partial differential equations read as follows
(2.7)∂βn(r)∂t-kβΔβn(r)-Pn=0forr∈Ωn,∂βi(r)∂t-kβΔβi(r)+Ci=0forr∈Ωiwithi=q,p,
and considering that the inhibitor is neither consumed nor produced in the host tissue, we obtain
(2.8)∂βe(r)∂t-kβΔβe(r)=0forr∈ℝ3-{Ωn∪Ωq∪Ωp}.Before we proceed, we will examine the significance of the terms with the time derivatives in the partial differential (2.6)–(2.8). Assuming the reference variables tg, rg and σ∞, we define the dimensionless variables t-=(t/tg),r-=(r/rg)and σ-i(r-)=σi(r)/σ∞for i=n,q,p,e. Hence, (2.6) becomes
(2.9)ε∂σ-i(r-)∂t--Δr-σ-i(r-)+rg2σ∞kσΓi=0withi=n,q,p,e,
where Δr-=rg2Δ, while the constant ε=rg2/(tgkσ) becomes much smaller than one, when the scaling rg=10-2 cm, tg=1day=86400sec, and kσ=10-6cm2/sec is used, which are a typical length scale, a typical time scale and a typical diffusion constant for the growth of an avascular tumour [9]. Then, ε=rg2/(tgkσ)=(10-4/(86400·10-6))=1,12·10-3, and the time derivative terms can be neglected. Under the same scaling, (2.7)-(2.8) can also be considered as elliptic partial differential equations. This means that when the problem is studied per day and the length is measured in μm, the diffusion of substances inwards to and outwards from the tumour is almost immediate and the tumour is always in a diffusive equilibrium state.


#### 2.1.2. Derivation of the Boundary Conditions

Let us consider an elementary test cylinder lying symmetrically across the tumour interfaceSi with i=n,q,p,e. Applying the mass conservation law in this cylinder, with elementary volume Vi and bases Si+, Si-, we follow similar arguments as in Section 2.1.1 to obtain
(2.10)∫Si+r^·∇σi+(r)ds-∫Si-r^·∇σi-(r)ds=∫ViΓikσdv,
where σi+ and σi- are the nutrient concentration at the exterior and the interior region to Si, respectively. Taking the limit at which the test volume is eliminated, both its bases coincide to the corresponding elementary section on Si and the volume integral of the right-hand side of (2.10) is eliminated, as Γi/kσis a continuous and bounded function. Applying then the integral mean value theorem, we obtain
(2.11)r^·∇σi+(r)-r^·∇σi-(r)=0forr∈Siwithi=n,q,p,e.Finally, we assume continuity conditions of the nutrient concentration on each interface and an external nutrient supply traced onSe as
(2.12)σ∞(re)=σ∞(re)g(θ,ϕ)forre∈Se.Similarly, we derive continuity boundary conditions for the inhibitor flow through the tumour’s interfaces and we assume that the inhibitor concentration is also continuous there. Finally, the asymptotic condition(2.13)limr→∞βe(r)=0
is supposed to hold uniformly over all directions.


### 2.2. The Pressure Field Problem

#### 2.2.1. Derivation of the Partial Differential Equations

Combining previous assumptions [6, 19], which attribute the cells movement either to the nutrient gradient or to the pressure gradient, as it is dictated by the Darcy law, we further assume that cells move towards the direction of the nutrient gradient and opposite to the direction of the pressure and inhibitor gradients. Denoting by vi(r) the velocity of the cell at the point r of the region Ωi with i=n,q,p,e, we have
(2.14)vi(r)=-μP∇Pi(r)+μσ∇σi(r)-μβ∇βi(r),
where μP, μσ, and μβ are constants of proportionality characterizing the motility of the cell. Applying the divergence operator on both sides of (2.14), we obtain
(2.15)∇·vi(r)=-μPΔPi(r)+μσΔσi(r)-μβΔβi(r).Here we recall that both the tumour tissue and the host tissue are assumed to be incompressible fluids with constant density. The mass conservation law implies that((∂ρt(r,t))/∂t)+∇·(ρt(r,t)vi(r,t))=Mi(r,t), where ρt(r,t) is the density of the tumour and Mi(r,t) is the mass per unit volume, per unit time that is produced or lost in the region Ωi for i=n,q,p. Since the density is considered to be constant, we derive that
(2.16)∇·vi(r,t)=Mi(r,t)ρt(r,t)≡Gi
in the tumour regions Ωi for i=n,q,p. We assume that Gi is not time dependent and it is invariable in every point inside Ωi. With similar arguments, we conclude that ∇·ve(r,t)=0 in Ωe, where no cell mass is neither produced nor lost.Substituting the partial differential (2.6) and (2.7)-(2.8) in (2.15) for each region, we conclude to the following:
(2.17)ΔPn(r)=-GnμP+μβμPPnkβ≡Fnforr∈Ωn,ΔPi(r)=-GiμP+μσμPΓikσ-μβμPCikβ≡Fiforr∈Ωiwithi=q,p,
while
(2.18)ΔPe(r)=0forr∈Ωe.


#### 2.2.2. Derivation of the Boundary Conditions

The tumour colony and its host surrounding are assumed to be incompressible fluids of different phase. So the Young-Laplace equation is considered on their interfaceSp, that is,
(2.19)Pp(r)-Pe(r)=aprpforr∈Sp,
where ap is a positive proportionality constant. Moreover, both the normal component and the tangential component of the pressure gradient are considered to be continuous across Sp in order to secure the tumour’s boundary compactness. Between the tumour’s inner interfaces we consider continuity conditions for both the pressure field and its normal derivatives, as the tumour’s distinct regions are considered to correspond to cells at different stage of their biological cycle but not to different fluids, characterized by different physical parameters. The tangential derivatives of the pressure field on the exterior surface Se are also considered to be continuous, so that it maintains its compactness. Finally, at Se, the pressure field’s trace is
(2.20)P∞(re)=P∞(re)h(θ,ϕ)forre∈Se.


### 2.3. The Evolution Equation and Relevant Topics

The aim of this work is to determine the evolution of the tumour’s exterior boundarySp under the assumption that this boundary evolves normally to itself. Considering that vi(rp)=(drp)/dt with i=n,q,p,e, then directly from (2.14), we easily get the following relationship:
(2.21)r^·drpdt=-μPr^·∇Pp(rp)+μσr^·∇σp(rp)-μβr^·∇βp(rp)
or
(2.22)drpdt=∂∂r(-μPPp(rp)+μσσp(rp)-μββp(rp)).We note here that relationship (2.22) is an ordinary differential equation with respect to the function rp(t). The uniqueness of its solution is secured by the initial condition rp(0)=Rp, where Rp is the initial radius of the Sp boundary.The right-hand side of (2.22) depends, as it will be shown in what follows, on the time-dependent boundaries rn(t), rq(t), rp(t), and re(t). Hence, (2.22) is solvable under constraints, which interrelate these boundaries and secure that (2.22) is dependent only on rp(t). These constraints are provided by the critical values of the nutrient and inhibitor concentrations. In particular, the critical nutrient value σn* determines if a cell dies out of starvation or not, so this value is met on the surface Sn, that is,
(2.23)σq(rn)=σn(rn)=σn*.Also, the nutrient valueσp* and the inhibitor value βp* determine if a cell proliferates or not, so these critical values are met on surface Sq, that is,
(2.24)σq(rq)=σp(rq)=σp*,(2.25)βq(rq)=βp(rq)=βp*.


## 2.1. The Nutrient and the Inhibitor Concentration Problem

### 2.1.1. Derivation of the Partial Differential Equations

According to the mass conservation law in each of the regionsΩi, i=n,q,p,e, the rate of change of the nutrient concentration in each region equals to the normal flow of nutrients Fin(r) inwards to Ωi, minus the flow outwards, Fout(r), and minus the amount of nutrient that is consumed in the interior of Ωi, that is,
(2.3)∫Ωi∂σi(r)∂tdv=∫∂Ωi(-n^)·Fin(r)ds-∫∂Ωin^·Fout(r)ds-∫ΩiΓidv,
where ∂Ωi is the boundary of the region Ωi, n^ is the outward unit normal vector to ∂Ωi, and Γi is the nutrient consumption rate, considered to be constant in Ωi. Since the nutrient flow is directed monotonically inwards in all spherical layers, the flow inwards takes place in the exterior boundary surface, denoted by Si+, and the flow outwards is realized from the inner boundary Si-. Considering also that for the spherical surface the normal direction is the radial one, (2.3) becomes
(2.4)∫Ωi∂σi(r)∂tdv=∫Si+(-r^)·Fin(r)ds-∫Si-(-r^)·Fout(r)ds-∫ΩiΓidvwithi=n,q,p,e.Assuming that the nutrient flow is only diffusive, it follows Fick’s law, which means that it is directed to regions of lower nutrient concentrationF(r)=-kσ∇σ(r). Then, relationship (2.4) yields
(2.5)∫Ωi∂σi(r)∂tdv=∫Si+r^·(kσ∇σi(r))ds-∫Si-r^·(kσ∇σi(r))ds-∫ΩiΓidv.Applying the divergence theorem inΩi in the right-hand side of (2.5), (2.5) becomes ∫Ωi(((∂σi(r))/∂t)-kσΔσi(r)+Γi)dv=0, which results in
(2.6)∂σi(r)∂t-kσΔσi(r)+Γi=0forr∈Ωiwithi=n,q,p,e.In particular, considering the physical assumptions made in the Introduction, the consumption rate in each region is defined asΓe=Γn=0, in Ωe and Ωn, respectively, where no consumption of nutrients takes place. In Ωq and Ωp, the consumption rate is assumed to depend on the vital state of the nutrients, which depends on the amount of the available nutrient, while it is defined as Γq=Γ(σn*/σ∞)and Γp=Γ(σp*/σ∞), respectively, where Γ is a proportional constant depending on the species of the host tissue cells, σn*, σp* are the nutrient critical values mentioned in the introduction, and σ∞ is a reference value corresponding to the nutrient concentration in the vicinity of the tumour’s host organ.We turn now to the inhibitor concentration. From the mathematical point of view, the inhibitor concentration problem is very similar to the nutrient concentration problem. The difference lies on the interpretation of the parameters and in their particular expressions. As mentioned in the introduction, the inhibitor is produced inside the necrotic core with production ratePn and is consumed in the quiescent and proliferating layer with consumption rate Ci with i=q,p, respectively. Therefore, the inhibitor flow is directed outwards and the corresponding partial differential equations read as follows
(2.7)∂βn(r)∂t-kβΔβn(r)-Pn=0forr∈Ωn,∂βi(r)∂t-kβΔβi(r)+Ci=0forr∈Ωiwithi=q,p,
and considering that the inhibitor is neither consumed nor produced in the host tissue, we obtain
(2.8)∂βe(r)∂t-kβΔβe(r)=0forr∈ℝ3-{Ωn∪Ωq∪Ωp}.Before we proceed, we will examine the significance of the terms with the time derivatives in the partial differential (2.6)–(2.8). Assuming the reference variables tg, rg and σ∞, we define the dimensionless variables t-=(t/tg),r-=(r/rg)and σ-i(r-)=σi(r)/σ∞for i=n,q,p,e. Hence, (2.6) becomes
(2.9)ε∂σ-i(r-)∂t--Δr-σ-i(r-)+rg2σ∞kσΓi=0withi=n,q,p,e,
where Δr-=rg2Δ, while the constant ε=rg2/(tgkσ) becomes much smaller than one, when the scaling rg=10-2 cm, tg=1day=86400sec, and kσ=10-6cm2/sec is used, which are a typical length scale, a typical time scale and a typical diffusion constant for the growth of an avascular tumour [9]. Then, ε=rg2/(tgkσ)=(10-4/(86400·10-6))=1,12·10-3, and the time derivative terms can be neglected. Under the same scaling, (2.7)-(2.8) can also be considered as elliptic partial differential equations. This means that when the problem is studied per day and the length is measured in μm, the diffusion of substances inwards to and outwards from the tumour is almost immediate and the tumour is always in a diffusive equilibrium state.


### 2.1.2. Derivation of the Boundary Conditions

Let us consider an elementary test cylinder lying symmetrically across the tumour interfaceSi with i=n,q,p,e. Applying the mass conservation law in this cylinder, with elementary volume Vi and bases Si+, Si-, we follow similar arguments as in Section 2.1.1 to obtain
(2.10)∫Si+r^·∇σi+(r)ds-∫Si-r^·∇σi-(r)ds=∫ViΓikσdv,
where σi+ and σi- are the nutrient concentration at the exterior and the interior region to Si, respectively. Taking the limit at which the test volume is eliminated, both its bases coincide to the corresponding elementary section on Si and the volume integral of the right-hand side of (2.10) is eliminated, as Γi/kσis a continuous and bounded function. Applying then the integral mean value theorem, we obtain
(2.11)r^·∇σi+(r)-r^·∇σi-(r)=0forr∈Siwithi=n,q,p,e.Finally, we assume continuity conditions of the nutrient concentration on each interface and an external nutrient supply traced onSe as
(2.12)σ∞(re)=σ∞(re)g(θ,ϕ)forre∈Se.Similarly, we derive continuity boundary conditions for the inhibitor flow through the tumour’s interfaces and we assume that the inhibitor concentration is also continuous there. Finally, the asymptotic condition(2.13)limr→∞βe(r)=0
is supposed to hold uniformly over all directions.


## 2.1.1. Derivation of the Partial Differential Equations

According to the mass conservation law in each of the regionsΩi, i=n,q,p,e, the rate of change of the nutrient concentration in each region equals to the normal flow of nutrients Fin(r) inwards to Ωi, minus the flow outwards, Fout(r), and minus the amount of nutrient that is consumed in the interior of Ωi, that is,
(2.3)∫Ωi∂σi(r)∂tdv=∫∂Ωi(-n^)·Fin(r)ds-∫∂Ωin^·Fout(r)ds-∫ΩiΓidv,
where ∂Ωi is the boundary of the region Ωi, n^ is the outward unit normal vector to ∂Ωi, and Γi is the nutrient consumption rate, considered to be constant in Ωi. Since the nutrient flow is directed monotonically inwards in all spherical layers, the flow inwards takes place in the exterior boundary surface, denoted by Si+, and the flow outwards is realized from the inner boundary Si-. Considering also that for the spherical surface the normal direction is the radial one, (2.3) becomes
(2.4)∫Ωi∂σi(r)∂tdv=∫Si+(-r^)·Fin(r)ds-∫Si-(-r^)·Fout(r)ds-∫ΩiΓidvwithi=n,q,p,e.Assuming that the nutrient flow is only diffusive, it follows Fick’s law, which means that it is directed to regions of lower nutrient concentrationF(r)=-kσ∇σ(r). Then, relationship (2.4) yields
(2.5)∫Ωi∂σi(r)∂tdv=∫Si+r^·(kσ∇σi(r))ds-∫Si-r^·(kσ∇σi(r))ds-∫ΩiΓidv.Applying the divergence theorem inΩi in the right-hand side of (2.5), (2.5) becomes ∫Ωi(((∂σi(r))/∂t)-kσΔσi(r)+Γi)dv=0, which results in
(2.6)∂σi(r)∂t-kσΔσi(r)+Γi=0forr∈Ωiwithi=n,q,p,e.In particular, considering the physical assumptions made in the Introduction, the consumption rate in each region is defined asΓe=Γn=0, in Ωe and Ωn, respectively, where no consumption of nutrients takes place. In Ωq and Ωp, the consumption rate is assumed to depend on the vital state of the nutrients, which depends on the amount of the available nutrient, while it is defined as Γq=Γ(σn*/σ∞)and Γp=Γ(σp*/σ∞), respectively, where Γ is a proportional constant depending on the species of the host tissue cells, σn*, σp* are the nutrient critical values mentioned in the introduction, and σ∞ is a reference value corresponding to the nutrient concentration in the vicinity of the tumour’s host organ.We turn now to the inhibitor concentration. From the mathematical point of view, the inhibitor concentration problem is very similar to the nutrient concentration problem. The difference lies on the interpretation of the parameters and in their particular expressions. As mentioned in the introduction, the inhibitor is produced inside the necrotic core with production ratePn and is consumed in the quiescent and proliferating layer with consumption rate Ci with i=q,p, respectively. Therefore, the inhibitor flow is directed outwards and the corresponding partial differential equations read as follows
(2.7)∂βn(r)∂t-kβΔβn(r)-Pn=0forr∈Ωn,∂βi(r)∂t-kβΔβi(r)+Ci=0forr∈Ωiwithi=q,p,
and considering that the inhibitor is neither consumed nor produced in the host tissue, we obtain
(2.8)∂βe(r)∂t-kβΔβe(r)=0forr∈ℝ3-{Ωn∪Ωq∪Ωp}.Before we proceed, we will examine the significance of the terms with the time derivatives in the partial differential (2.6)–(2.8). Assuming the reference variables tg, rg and σ∞, we define the dimensionless variables t-=(t/tg),r-=(r/rg)and σ-i(r-)=σi(r)/σ∞for i=n,q,p,e. Hence, (2.6) becomes
(2.9)ε∂σ-i(r-)∂t--Δr-σ-i(r-)+rg2σ∞kσΓi=0withi=n,q,p,e,
where Δr-=rg2Δ, while the constant ε=rg2/(tgkσ) becomes much smaller than one, when the scaling rg=10-2 cm, tg=1day=86400sec, and kσ=10-6cm2/sec is used, which are a typical length scale, a typical time scale and a typical diffusion constant for the growth of an avascular tumour [9]. Then, ε=rg2/(tgkσ)=(10-4/(86400·10-6))=1,12·10-3, and the time derivative terms can be neglected. Under the same scaling, (2.7)-(2.8) can also be considered as elliptic partial differential equations. This means that when the problem is studied per day and the length is measured in μm, the diffusion of substances inwards to and outwards from the tumour is almost immediate and the tumour is always in a diffusive equilibrium state.


## 2.1.2. Derivation of the Boundary Conditions

Let us consider an elementary test cylinder lying symmetrically across the tumour interfaceSi with i=n,q,p,e. Applying the mass conservation law in this cylinder, with elementary volume Vi and bases Si+, Si-, we follow similar arguments as in Section 2.1.1 to obtain
(2.10)∫Si+r^·∇σi+(r)ds-∫Si-r^·∇σi-(r)ds=∫ViΓikσdv,
where σi+ and σi- are the nutrient concentration at the exterior and the interior region to Si, respectively. Taking the limit at which the test volume is eliminated, both its bases coincide to the corresponding elementary section on Si and the volume integral of the right-hand side of (2.10) is eliminated, as Γi/kσis a continuous and bounded function. Applying then the integral mean value theorem, we obtain
(2.11)r^·∇σi+(r)-r^·∇σi-(r)=0forr∈Siwithi=n,q,p,e.Finally, we assume continuity conditions of the nutrient concentration on each interface and an external nutrient supply traced onSe as
(2.12)σ∞(re)=σ∞(re)g(θ,ϕ)forre∈Se.Similarly, we derive continuity boundary conditions for the inhibitor flow through the tumour’s interfaces and we assume that the inhibitor concentration is also continuous there. Finally, the asymptotic condition(2.13)limr→∞βe(r)=0
is supposed to hold uniformly over all directions.


## 2.2. The Pressure Field Problem

### 2.2.1. Derivation of the Partial Differential Equations

Combining previous assumptions [6, 19], which attribute the cells movement either to the nutrient gradient or to the pressure gradient, as it is dictated by the Darcy law, we further assume that cells move towards the direction of the nutrient gradient and opposite to the direction of the pressure and inhibitor gradients. Denoting by vi(r) the velocity of the cell at the point r of the region Ωi with i=n,q,p,e, we have
(2.14)vi(r)=-μP∇Pi(r)+μσ∇σi(r)-μβ∇βi(r),
where μP, μσ, and μβ are constants of proportionality characterizing the motility of the cell. Applying the divergence operator on both sides of (2.14), we obtain
(2.15)∇·vi(r)=-μPΔPi(r)+μσΔσi(r)-μβΔβi(r).Here we recall that both the tumour tissue and the host tissue are assumed to be incompressible fluids with constant density. The mass conservation law implies that((∂ρt(r,t))/∂t)+∇·(ρt(r,t)vi(r,t))=Mi(r,t), where ρt(r,t) is the density of the tumour and Mi(r,t) is the mass per unit volume, per unit time that is produced or lost in the region Ωi for i=n,q,p. Since the density is considered to be constant, we derive that
(2.16)∇·vi(r,t)=Mi(r,t)ρt(r,t)≡Gi
in the tumour regions Ωi for i=n,q,p. We assume that Gi is not time dependent and it is invariable in every point inside Ωi. With similar arguments, we conclude that ∇·ve(r,t)=0 in Ωe, where no cell mass is neither produced nor lost.Substituting the partial differential (2.6) and (2.7)-(2.8) in (2.15) for each region, we conclude to the following:
(2.17)ΔPn(r)=-GnμP+μβμPPnkβ≡Fnforr∈Ωn,ΔPi(r)=-GiμP+μσμPΓikσ-μβμPCikβ≡Fiforr∈Ωiwithi=q,p,
while
(2.18)ΔPe(r)=0forr∈Ωe.


### 2.2.2. Derivation of the Boundary Conditions

The tumour colony and its host surrounding are assumed to be incompressible fluids of different phase. So the Young-Laplace equation is considered on their interfaceSp, that is,
(2.19)Pp(r)-Pe(r)=aprpforr∈Sp,
where ap is a positive proportionality constant. Moreover, both the normal component and the tangential component of the pressure gradient are considered to be continuous across Sp in order to secure the tumour’s boundary compactness. Between the tumour’s inner interfaces we consider continuity conditions for both the pressure field and its normal derivatives, as the tumour’s distinct regions are considered to correspond to cells at different stage of their biological cycle but not to different fluids, characterized by different physical parameters. The tangential derivatives of the pressure field on the exterior surface Se are also considered to be continuous, so that it maintains its compactness. Finally, at Se, the pressure field’s trace is
(2.20)P∞(re)=P∞(re)h(θ,ϕ)forre∈Se.


## 2.2.1. Derivation of the Partial Differential Equations

Combining previous assumptions [6, 19], which attribute the cells movement either to the nutrient gradient or to the pressure gradient, as it is dictated by the Darcy law, we further assume that cells move towards the direction of the nutrient gradient and opposite to the direction of the pressure and inhibitor gradients. Denoting by vi(r) the velocity of the cell at the point r of the region Ωi with i=n,q,p,e, we have
(2.14)vi(r)=-μP∇Pi(r)+μσ∇σi(r)-μβ∇βi(r),
where μP, μσ, and μβ are constants of proportionality characterizing the motility of the cell. Applying the divergence operator on both sides of (2.14), we obtain
(2.15)∇·vi(r)=-μPΔPi(r)+μσΔσi(r)-μβΔβi(r).Here we recall that both the tumour tissue and the host tissue are assumed to be incompressible fluids with constant density. The mass conservation law implies that((∂ρt(r,t))/∂t)+∇·(ρt(r,t)vi(r,t))=Mi(r,t), where ρt(r,t) is the density of the tumour and Mi(r,t) is the mass per unit volume, per unit time that is produced or lost in the region Ωi for i=n,q,p. Since the density is considered to be constant, we derive that
(2.16)∇·vi(r,t)=Mi(r,t)ρt(r,t)≡Gi
in the tumour regions Ωi for i=n,q,p. We assume that Gi is not time dependent and it is invariable in every point inside Ωi. With similar arguments, we conclude that ∇·ve(r,t)=0 in Ωe, where no cell mass is neither produced nor lost.Substituting the partial differential (2.6) and (2.7)-(2.8) in (2.15) for each region, we conclude to the following:
(2.17)ΔPn(r)=-GnμP+μβμPPnkβ≡Fnforr∈Ωn,ΔPi(r)=-GiμP+μσμPΓikσ-μβμPCikβ≡Fiforr∈Ωiwithi=q,p,
while
(2.18)ΔPe(r)=0forr∈Ωe.


## 2.2.2. Derivation of the Boundary Conditions

The tumour colony and its host surrounding are assumed to be incompressible fluids of different phase. So the Young-Laplace equation is considered on their interfaceSp, that is,
(2.19)Pp(r)-Pe(r)=aprpforr∈Sp,
where ap is a positive proportionality constant. Moreover, both the normal component and the tangential component of the pressure gradient are considered to be continuous across Sp in order to secure the tumour’s boundary compactness. Between the tumour’s inner interfaces we consider continuity conditions for both the pressure field and its normal derivatives, as the tumour’s distinct regions are considered to correspond to cells at different stage of their biological cycle but not to different fluids, characterized by different physical parameters. The tangential derivatives of the pressure field on the exterior surface Se are also considered to be continuous, so that it maintains its compactness. Finally, at Se, the pressure field’s trace is
(2.20)P∞(re)=P∞(re)h(θ,ϕ)forre∈Se.


## 2.3. The Evolution Equation and Relevant Topics

The aim of this work is to determine the evolution of the tumour’s exterior boundarySp under the assumption that this boundary evolves normally to itself. Considering that vi(rp)=(drp)/dt with i=n,q,p,e, then directly from (2.14), we easily get the following relationship:
(2.21)r^·drpdt=-μPr^·∇Pp(rp)+μσr^·∇σp(rp)-μβr^·∇βp(rp)
or
(2.22)drpdt=∂∂r(-μPPp(rp)+μσσp(rp)-μββp(rp)).We note here that relationship (2.22) is an ordinary differential equation with respect to the function rp(t). The uniqueness of its solution is secured by the initial condition rp(0)=Rp, where Rp is the initial radius of the Sp boundary.The right-hand side of (2.22) depends, as it will be shown in what follows, on the time-dependent boundaries rn(t), rq(t), rp(t), and re(t). Hence, (2.22) is solvable under constraints, which interrelate these boundaries and secure that (2.22) is dependent only on rp(t). These constraints are provided by the critical values of the nutrient and inhibitor concentrations. In particular, the critical nutrient value σn* determines if a cell dies out of starvation or not, so this value is met on the surface Sn, that is,
(2.23)σq(rn)=σn(rn)=σn*.Also, the nutrient valueσp* and the inhibitor value βp* determine if a cell proliferates or not, so these critical values are met on surface Sq, that is,
(2.24)σq(rq)=σp(rq)=σp*,(2.25)βq(rq)=βp(rq)=βp*.


## 3. Model Analysis

In what follows the boundary value problems for the nutrient concentration field, for the inhibitor concentration, and for the pressure field are solved, using standard spectral analysis of the Laplace operator. The inhomogeneity imposed by the nutrient supply is reflected upon the conditions (2.23) and (2.24), which are considered to hold either globally on the whole surfaces Sn and Sq, respectively, or as an average. Both aspects lead to the same evolution equation but allow different distribution of the nutrient supply and of the pressure on Se. The latter approach is illustrated by an example in Section 4.1, where a special inhomogeneity of the nutrient supply is assumed.


### 3.1. The Nutrient Concentration Problem

With respect to Section2.1, the boundary value problem that determines the nutrient concentration field assumes the form
(3.1)Δσi(r)=0forr∈Ωiwithi=n,e,
while
(3.2)Δσi(r)=Γikσ≡γiforr∈Ωiwithi=q,p,
where the boundary conditions to be satisfied are
(3.3)σq(r)=σn(r)forr∈Sn,∂∂rσq(r)-∂∂rσn(r)=0forr∈Sn,
while
(3.4)σp(r)=σq(r)forr∈Sq,∂∂rσp(r)-∂∂rσq(r)=0forr∈Sq,σe(r)=σp(r)forr∈Sp,∂∂rσe(r)-∂∂rσp(r)=0forr∈Sp,
whilst at the outer boundary we obtain
(3.5)σ∞(re)=σ∞(re)g(r^e)=σ∞(re)g(θ,ϕ)forre∈Se.Before we proceed with the solution of the problem, it is worth mentioning the eigenfunctions we will use, which are the surface spherical harmonics, defined as(3.6)Ylm,s(θ,ϕ)=Plm(cosθ){cosmϕ,s=esinmϕ,s=o,
which satisfy the following orthogonality property:
(3.7)∫S2Ylm,s(θ,ϕ)Yl′m′,s′(θ,ϕ)ds(θ,ϕ)=1εm4π(2l+1)(l+m)!(l-m)!δll′δmm′δss′,whereεm={1,m=02,m≥1.Moreover, for notation convenience, we will replace the triple summation symbol by the following∑l=0+∞∑m=0l∑s=e,o⋯≡∑l,m,s⋯ in all the expansions.Within this aim, we expand the nutrient concentration in each domain, in spherical harmonics, demanding the field to be continuous at the center of the spherical core. In that way we, obtain(3.8)σn(r)=∑l,m,sal,σ(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,σq(r)=γq6r2+∑l,m,s(bl,σ(i)m,srl+bl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,σp(r)=γp6r2+∑l,m,s(cl,σ(i)m,srl+cl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,σe(r)=∑l,m,s(dl,σ(i)m,srl+dl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωe.In addition, the spherical expansion of the nutrient supply fieldσ∞(re) on Se reads as
(3.9)σ∞(re)=σ∞(re)g(θ,ϕ)=σ∞(re)∑l,m,sglm,sYlm,s(θ,ϕ),
where glm,s=((2l+1)/4π)((l-m)!/(l+m)!)εm∫02π∫0πg(θ,ϕ)Ylm,s(θ,ϕ)sinθdθdϕ. Applying boundary conditions (3.3)–(3.5) into expansions (3.8) with (3.9), one has to deal with a seventh-order linear algebraic system with respect to the seven sequences of unknown coefficients. Long but straightforward calculations lead to the nutrient concentration field, where σi(r,θ,ϕ) is defined for r∈Ωi with i=n,q,p,e, that is,
(3.10)σn(r,θ,ϕ)=γp-γq2rq2(1-23rqre)+γq2rn2(1-23rnre)-γp2rp2(1-23rpre)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.11)σq(r,θ,ϕ)=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1r-1re)+γq6r2+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.12)σp(r,θ,ϕ)=γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.13)σe(r,θ,ϕ)=[γq3rn3+γp-γq3rq3-γp3rp3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ).Substituting expression (3.10) into condition (2.23), we perform trivial calculations to obtain the following equation:
(3.14)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞(re)∑l,m,sglm,s(rnre)lYlm,s(θ,ϕ).Equation (3.14) can be considered in two ways, either as one that holds true pointwise on all (rn,θ,ϕ)∈Sn or as an equation that holds true at the average over the surface Sn.In the first case, we imply the orthogonality property ofYlm,s(θ,ϕ) in (3.14) and we conclude that all coefficients glm,s for l≥1 should vanish and the nutrient supply should be defined uniformly over all directions on Se, meaning that σ∞(re)=σ∞(re), where we have assumed for simplicity and without loss of generality that g00,e=1.In the second case, we take the mean value of both sides of (3.14) by integrating over Sn. Due to orthogonality of Ylm,s(θ,ϕ), the angular dependence will vanish and the result will read as
(3.15)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞(re),
which is also the result if the first argument is considered. The difference lies in that, considering (3.14) as an average equation, the external nutrient supply is allowed to be defined nonuniformly. Though the drawback in this case lies on the perturbation implied in the locus of the points that enjoy the critical nutrient concentration σn*, which is no longer a spherical surface.In what follows, the treatment of either case is indifferent. Nevertheless, we will come back to the different approaches in Section3.4, when the evolution equation is studied. Therefore, relationship (2.24) along with (3.12) gives in a similar way
(3.16)σp*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rq-1re)+γq6rq2+σ∞(re).Both (3.15) and (3.16) should hold throughout all time of the tumour’s evolution, which provides restrictions on the evolution of all boundaries rn(t), rq(t), rp(t) and re(t). In addition, from relations (3.15) and (3.16), we are provided with a relation depending only on rn(t) and rq(t), that is,
(3.17)σn*-σp*=γq3rn3(1rn-1rq)+γq6(rn2-rq2).


### 3.2. The Inhibitor Concentration Problem

The spherical expansions of the inhibitor concentration that solve the boundary value problem (2.7)-(2.8), with continuous boundary conditions and also with (2.13), assume the following forms:
(3.18)βn(r)=pn6r2+∑l,m,sal,β(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,βq(r)=cq6r2+∑l,m,s(bl,β(i)m,srl+bl,β(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,βp(r)=cp6r2+∑l,m,s(cl,β(i)m,srl+cl,β(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,βe(r)=∑l,m,sdl,β(e)m,sr-(l+1)Ylm,s(θ,ϕ)forr∈ℝ3-(Ωp∪Ωq∪Ωn),
where for simplicity we have used the notation -(Pn/kβ)=pn and (Ci/kβ)=ci with i=p,q.Following similar track of calculations as in Section3.1, we arrive at the following expression for the inhibitor concentration field:
(3.19)βn(r)=cq-pn2rn2+cp-cq2rq2-cp2rp2+pn6r2,(3.20)βq(r)=cp-cq2rq2-cp2rp2+cq-pn3rn3r+cq6r2,(3.21)βp(r)=-cp2rp2+[cq-pn3rn3+cp-cq3rq3]1r+cp6r2,(3.22)βe(r)=[cq-pn3rn3+cp-cq3rq3-cp3rp3]1r.The inhibitor concentration is, by definition, a decreasing function ofr in ℝ3-Ωn. The critical inhibitor concentration βp* appears on Sq, and its further decrease permits cell proliferation in Ωp. Hence, (2.25) together with (3.20) leads to
(3.23)βp*=(cp2-cq3)rq2-cp2rp2+cq-pn3rn3rq,
which provides another restriction between the boundaries rn(t), rq(t), and rp(t), which should hold for every t>0.


### 3.3. The Pressure Field

The pressure field satisfies the boundary value problem (2.17)–(2.20) completed with continuity conditions for both the field and its gradient in all inner interfaces. The corresponding spherical expansions read as
(3.24)Pn(r)=Fn6r2+∑l,m,sal,P(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,Pq(r)=Fq6r2+∑l,m,s(bl,P(i)m,srl+bl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,Pp(r)=Fp6r2+∑l,m,s(cl,P(i)m,srl+cl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,Pe(r)=∑l,m,s(dl,P(i)m,srl+dl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωe,
while the exterior pressure distribution P∞(re) on Se is expanded in spherical harmonics as
(3.25)P∞(re)=P∞(re)h(θ,ϕ)=P∞(re)∑l,m,shlm,sYlm,s(θ,ϕ),
where hlm,s=((2l+1)/4π)((l-m)!/(l+m)!))εm∫02π∫0πh(θ,ϕ)Ylm,s(θ,ϕ)sinθdθdϕ. The solution of the above boundary value problem is obtained in the form
(3.26)Pn(r)=Fq-Fn2rn2(1-23rnre)+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+aprp+Fn6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.27)Pq(r)=Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+aprp+Fq-Fn3rn3(1r-1re)+Fq6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.28)Pp(r)=-Fp2rp2(1-23rpre)+aprp+(Fq-Fn3rn3+Fp-Fq3rq3)(1r-1re)+Fp6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.29)Pe(r)=(Fq-Fn3rn3+Fp-Fq3rq3-Fp3rp3)(1r-1re)+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ).


### 3.4. The Evolution Equation

In the sequel, we collect all the results for the concentration fields and for the pressure field, as traced onSp, and we substitute them into the differential equation (2.22). Here, we have to recall the two different approaches that are possible for the nutrient concentration field, mentioned in Section 3.1.


#### 3.4.1. The Homogeneous Approach

Following the demand that the nutrient trace on the exterior boundarySe is homogeneous, the nutrient concentration (3.12) in the proliferating region is
(3.30)σp(r,θ,ϕ)=σ∞(re)+γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re).Substituting (3.30), (3.21), and (3.28) in (2.22) and performing algebraic manipulations the following evolution equation is obtained:
(3.31)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ).We note here that the whole model’s treatment is based on the assumption that as the tumour evolves, all its boundaries remain members of the concentric spherical coordinate system, so that all the Laplace and the Poisson partial differential equations can be treated analytically through the separation of variables method. This implies thatrp=rp(t), and consequently its derivative cannot depend on the direction r^≡(θ,ϕ). As a consequence, only the zeroth-order term of the spherical expansion can be included in (3.31), a fact that eliminates the effect of any pressure that is imposed to the host tissue by the exterior medium on the tumour’s evolution. Moreover, the host tissue’s boundary re=re(t) itself does not seem to be involved by the tumour’s evolution. Therefore, the evolution equation is written as follows:
(3.32)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3},
where the physical implementation was secured mathematically by imposing hlm,s=0 for l≥1, m=0,1,…,l, and s=e,o into (3.31) concluding to (3.32). This equation describes the transitory phase until it reaches its steady state. From (3.32), one can obtain the radii of the tumour’s steady-state boundaries by setting the right-hand side equal to zero.


#### 3.4.2. The Nonhomogeneous Approach

By considering (3.14) as an average condition, the nutrient concentration in the proliferation shell reads as
(3.33)σp(r,θ,ϕ)=γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ).Then, the corresponding evolution equation fordrp/dt becomes
(3.34)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ)+μσσ∞(re)rp∑l,m,sglm,sl(rpre)lYlm,s(θ,ϕ),
or
(3.35)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}+1rp∑l,m,sl(rpre)l(μσσ∞(re)glm,s-μpP∞(re)hlm,s)Ylm,s(θ,ϕ).The arguments stated in Section3.4.1 hold true in this case as well, provided that between the spherical coefficients of the nutrient exterior field and the pressure field, as they are dictated on Se, the following condition is considered:
(3.36)glm,s=μpP∞(re)μσσ∞(re)hlm,sforl=1,2,…,m=1,2,…,lands=e,o,
and hence
(3.37)σ∞(re)=σ∞(re)-μpμσP∞(re)h00,e+μpμσP∞(re)
or
(3.38)μσ(σ∞(re)-σ∞(re))=μp(P∞(re)-P∞(re)h00,e).Consequently, the tumour’s boundaries evolve as concentric spheres and the evolution (3.35), provided (3.36), coincides with (3.32). From there after, it is obvious that the nonhomogeneous and the homogeneous approaches are treated in the same way.We stress that the restrictions on the type of inhomogeneity exhibited by the exterior pressure are implied by the assumption (2.22) for the evolution equation, in the model presented in this work. Nevertheless, following a different approach on the evolution condition, one could investigate further the possibility of imposing an inhomogeneous pressure trace upon Se, which is under our current investigation. For this reason, we do not imply this restriction on the results of Section 3.3, since they hold true up to the point that they are imposed properly in Section 3.4 and so after.Equation (3.32) is a nonlinear ordinary differential equation, which includes the rates (rn/rp) and (rq/rp), accompanied by rp, which are all unknown functions of time. In order to make (3.32) solvable, we need the expressions that connect all the unknown rates with function rp=rp(t). This is accomplished via expressions (3.15)–(3.17) and (3.23). In particular, (3.17) can be written as
(3.39)σn*-σp*rp2=γq3(rnrp)3(rprn-rprq)+γq6((rnrp)2-(rqrp)2),
or using the following simple symbols for the rates under consideration, those are (rn/rp)≡x=x(t), (rq/rp)≡y=y(t), and (rp/re)≡z=z(t), where by definition 0<x,y,z<1, expression (3.39) yields
(3.40)-γq3x3+γq2x2y-γq6y3-σn*-σp*rp2y=0.The third-degree nonlinear algebraic equation (3.40) provides x as a function of y and rp. Similarly, (3.23) assumes the form
(3.41)cq-pn3x3+(cp2-cq3)y3-[cp2+βp*rp2]y=0.Equations (3.40) and (3.41) form a nonlinear algebraic system, which when it is solved, provide the rates (rn/rp)≡x(t) and (rq/rp)=y(t) with respect to rp(t) and the differential (3.32) is solvable, under the initial condition rp(0)=Rp, which corresponds to the first observed tumour radius.Continuing, the rate(rp/re) is also connected to rp=rp(t) via (3.15), or equivalently (3.16); it is worthwhile to specify this connection, as it provides the evolution of the host boundary and also it appears in the nutrient field, the inhibitor field, and the pressure field provided in Sections 3.1, 3.2, and 3.3, respectively. In particular, using (3.15), we arrive at
(3.42)σn*-σ∞(re)rp2=γp-γq2(rqrp)2(1-23rqrprpre)-γp2(1-23rpre)+γq3(rnrp)3(rprn-rpre)+γq6(rnrp)2
or, in terms of the x,y,z symbols,
(3.43)(-γq3x3-γp-γq3y3+γp3)z=-γq2x2-γp-γq2y2+γp2+σn*-σ∞(re)rp2,
and if one substitutes the (x,y) solution of the system (3.40) and (3.41), then (3.43) provides the rate (rp/re)≡z(t) as a function of rp(t).


## 3.1. The Nutrient Concentration Problem

With respect to Section2.1, the boundary value problem that determines the nutrient concentration field assumes the form
(3.1)Δσi(r)=0forr∈Ωiwithi=n,e,
while
(3.2)Δσi(r)=Γikσ≡γiforr∈Ωiwithi=q,p,
where the boundary conditions to be satisfied are
(3.3)σq(r)=σn(r)forr∈Sn,∂∂rσq(r)-∂∂rσn(r)=0forr∈Sn,
while
(3.4)σp(r)=σq(r)forr∈Sq,∂∂rσp(r)-∂∂rσq(r)=0forr∈Sq,σe(r)=σp(r)forr∈Sp,∂∂rσe(r)-∂∂rσp(r)=0forr∈Sp,
whilst at the outer boundary we obtain
(3.5)σ∞(re)=σ∞(re)g(r^e)=σ∞(re)g(θ,ϕ)forre∈Se.Before we proceed with the solution of the problem, it is worth mentioning the eigenfunctions we will use, which are the surface spherical harmonics, defined as(3.6)Ylm,s(θ,ϕ)=Plm(cosθ){cosmϕ,s=esinmϕ,s=o,
which satisfy the following orthogonality property:
(3.7)∫S2Ylm,s(θ,ϕ)Yl′m′,s′(θ,ϕ)ds(θ,ϕ)=1εm4π(2l+1)(l+m)!(l-m)!δll′δmm′δss′,whereεm={1,m=02,m≥1.Moreover, for notation convenience, we will replace the triple summation symbol by the following∑l=0+∞∑m=0l∑s=e,o⋯≡∑l,m,s⋯ in all the expansions.Within this aim, we expand the nutrient concentration in each domain, in spherical harmonics, demanding the field to be continuous at the center of the spherical core. In that way we, obtain(3.8)σn(r)=∑l,m,sal,σ(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,σq(r)=γq6r2+∑l,m,s(bl,σ(i)m,srl+bl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,σp(r)=γp6r2+∑l,m,s(cl,σ(i)m,srl+cl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,σe(r)=∑l,m,s(dl,σ(i)m,srl+dl,σ(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωe.In addition, the spherical expansion of the nutrient supply fieldσ∞(re) on Se reads as
(3.9)σ∞(re)=σ∞(re)g(θ,ϕ)=σ∞(re)∑l,m,sglm,sYlm,s(θ,ϕ),
where glm,s=((2l+1)/4π)((l-m)!/(l+m)!)εm∫02π∫0πg(θ,ϕ)Ylm,s(θ,ϕ)sinθdθdϕ. Applying boundary conditions (3.3)–(3.5) into expansions (3.8) with (3.9), one has to deal with a seventh-order linear algebraic system with respect to the seven sequences of unknown coefficients. Long but straightforward calculations lead to the nutrient concentration field, where σi(r,θ,ϕ) is defined for r∈Ωi with i=n,q,p,e, that is,
(3.10)σn(r,θ,ϕ)=γp-γq2rq2(1-23rqre)+γq2rn2(1-23rnre)-γp2rp2(1-23rpre)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.11)σq(r,θ,ϕ)=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1r-1re)+γq6r2+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.12)σp(r,θ,ϕ)=γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ),(3.13)σe(r,θ,ϕ)=[γq3rn3+γp-γq3rq3-γp3rp3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ).Substituting expression (3.10) into condition (2.23), we perform trivial calculations to obtain the following equation:
(3.14)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞(re)∑l,m,sglm,s(rnre)lYlm,s(θ,ϕ).Equation (3.14) can be considered in two ways, either as one that holds true pointwise on all (rn,θ,ϕ)∈Sn or as an equation that holds true at the average over the surface Sn.In the first case, we imply the orthogonality property ofYlm,s(θ,ϕ) in (3.14) and we conclude that all coefficients glm,s for l≥1 should vanish and the nutrient supply should be defined uniformly over all directions on Se, meaning that σ∞(re)=σ∞(re), where we have assumed for simplicity and without loss of generality that g00,e=1.In the second case, we take the mean value of both sides of (3.14) by integrating over Sn. Due to orthogonality of Ylm,s(θ,ϕ), the angular dependence will vanish and the result will read as
(3.15)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞(re),
which is also the result if the first argument is considered. The difference lies in that, considering (3.14) as an average equation, the external nutrient supply is allowed to be defined nonuniformly. Though the drawback in this case lies on the perturbation implied in the locus of the points that enjoy the critical nutrient concentration σn*, which is no longer a spherical surface.In what follows, the treatment of either case is indifferent. Nevertheless, we will come back to the different approaches in Section3.4, when the evolution equation is studied. Therefore, relationship (2.24) along with (3.12) gives in a similar way
(3.16)σp*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rq-1re)+γq6rq2+σ∞(re).Both (3.15) and (3.16) should hold throughout all time of the tumour’s evolution, which provides restrictions on the evolution of all boundaries rn(t), rq(t), rp(t) and re(t). In addition, from relations (3.15) and (3.16), we are provided with a relation depending only on rn(t) and rq(t), that is,
(3.17)σn*-σp*=γq3rn3(1rn-1rq)+γq6(rn2-rq2).


## 3.2. The Inhibitor Concentration Problem

The spherical expansions of the inhibitor concentration that solve the boundary value problem (2.7)-(2.8), with continuous boundary conditions and also with (2.13), assume the following forms:
(3.18)βn(r)=pn6r2+∑l,m,sal,β(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,βq(r)=cq6r2+∑l,m,s(bl,β(i)m,srl+bl,β(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,βp(r)=cp6r2+∑l,m,s(cl,β(i)m,srl+cl,β(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,βe(r)=∑l,m,sdl,β(e)m,sr-(l+1)Ylm,s(θ,ϕ)forr∈ℝ3-(Ωp∪Ωq∪Ωn),
where for simplicity we have used the notation -(Pn/kβ)=pn and (Ci/kβ)=ci with i=p,q.Following similar track of calculations as in Section3.1, we arrive at the following expression for the inhibitor concentration field:
(3.19)βn(r)=cq-pn2rn2+cp-cq2rq2-cp2rp2+pn6r2,(3.20)βq(r)=cp-cq2rq2-cp2rp2+cq-pn3rn3r+cq6r2,(3.21)βp(r)=-cp2rp2+[cq-pn3rn3+cp-cq3rq3]1r+cp6r2,(3.22)βe(r)=[cq-pn3rn3+cp-cq3rq3-cp3rp3]1r.The inhibitor concentration is, by definition, a decreasing function ofr in ℝ3-Ωn. The critical inhibitor concentration βp* appears on Sq, and its further decrease permits cell proliferation in Ωp. Hence, (2.25) together with (3.20) leads to
(3.23)βp*=(cp2-cq3)rq2-cp2rp2+cq-pn3rn3rq,
which provides another restriction between the boundaries rn(t), rq(t), and rp(t), which should hold for every t>0.


## 3.3. The Pressure Field

The pressure field satisfies the boundary value problem (2.17)–(2.20) completed with continuity conditions for both the field and its gradient in all inner interfaces. The corresponding spherical expansions read as
(3.24)Pn(r)=Fn6r2+∑l,m,sal,P(i)m,srlYlm,s(θ,ϕ)forr∈Ωn,Pq(r)=Fq6r2+∑l,m,s(bl,P(i)m,srl+bl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωq,Pp(r)=Fp6r2+∑l,m,s(cl,P(i)m,srl+cl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωp,Pe(r)=∑l,m,s(dl,P(i)m,srl+dl,P(e)m,sr-(l+1))Ylm,s(θ,ϕ)forr∈Ωe,
while the exterior pressure distribution P∞(re) on Se is expanded in spherical harmonics as
(3.25)P∞(re)=P∞(re)h(θ,ϕ)=P∞(re)∑l,m,shlm,sYlm,s(θ,ϕ),
where hlm,s=((2l+1)/4π)((l-m)!/(l+m)!))εm∫02π∫0πh(θ,ϕ)Ylm,s(θ,ϕ)sinθdθdϕ. The solution of the above boundary value problem is obtained in the form
(3.26)Pn(r)=Fq-Fn2rn2(1-23rnre)+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+aprp+Fn6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.27)Pq(r)=Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+aprp+Fq-Fn3rn3(1r-1re)+Fq6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.28)Pp(r)=-Fp2rp2(1-23rpre)+aprp+(Fq-Fn3rn3+Fp-Fq3rq3)(1r-1re)+Fp6r2+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ),(3.29)Pe(r)=(Fq-Fn3rn3+Fp-Fq3rq3-Fp3rp3)(1r-1re)+P∞(re)∑l,m,shlm,s(rre)lYlm,s(θ,ϕ).


## 3.4. The Evolution Equation

In the sequel, we collect all the results for the concentration fields and for the pressure field, as traced onSp, and we substitute them into the differential equation (2.22). Here, we have to recall the two different approaches that are possible for the nutrient concentration field, mentioned in Section 3.1.


### 3.4.1. The Homogeneous Approach

Following the demand that the nutrient trace on the exterior boundarySe is homogeneous, the nutrient concentration (3.12) in the proliferating region is
(3.30)σp(r,θ,ϕ)=σ∞(re)+γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re).Substituting (3.30), (3.21), and (3.28) in (2.22) and performing algebraic manipulations the following evolution equation is obtained:
(3.31)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ).We note here that the whole model’s treatment is based on the assumption that as the tumour evolves, all its boundaries remain members of the concentric spherical coordinate system, so that all the Laplace and the Poisson partial differential equations can be treated analytically through the separation of variables method. This implies thatrp=rp(t), and consequently its derivative cannot depend on the direction r^≡(θ,ϕ). As a consequence, only the zeroth-order term of the spherical expansion can be included in (3.31), a fact that eliminates the effect of any pressure that is imposed to the host tissue by the exterior medium on the tumour’s evolution. Moreover, the host tissue’s boundary re=re(t) itself does not seem to be involved by the tumour’s evolution. Therefore, the evolution equation is written as follows:
(3.32)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3},
where the physical implementation was secured mathematically by imposing hlm,s=0 for l≥1, m=0,1,…,l, and s=e,o into (3.31) concluding to (3.32). This equation describes the transitory phase until it reaches its steady state. From (3.32), one can obtain the radii of the tumour’s steady-state boundaries by setting the right-hand side equal to zero.


### 3.4.2. The Nonhomogeneous Approach

By considering (3.14) as an average condition, the nutrient concentration in the proliferation shell reads as
(3.33)σp(r,θ,ϕ)=γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ).Then, the corresponding evolution equation fordrp/dt becomes
(3.34)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ)+μσσ∞(re)rp∑l,m,sglm,sl(rpre)lYlm,s(θ,ϕ),
or
(3.35)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}+1rp∑l,m,sl(rpre)l(μσσ∞(re)glm,s-μpP∞(re)hlm,s)Ylm,s(θ,ϕ).The arguments stated in Section3.4.1 hold true in this case as well, provided that between the spherical coefficients of the nutrient exterior field and the pressure field, as they are dictated on Se, the following condition is considered:
(3.36)glm,s=μpP∞(re)μσσ∞(re)hlm,sforl=1,2,…,m=1,2,…,lands=e,o,
and hence
(3.37)σ∞(re)=σ∞(re)-μpμσP∞(re)h00,e+μpμσP∞(re)
or
(3.38)μσ(σ∞(re)-σ∞(re))=μp(P∞(re)-P∞(re)h00,e).Consequently, the tumour’s boundaries evolve as concentric spheres and the evolution (3.35), provided (3.36), coincides with (3.32). From there after, it is obvious that the nonhomogeneous and the homogeneous approaches are treated in the same way.We stress that the restrictions on the type of inhomogeneity exhibited by the exterior pressure are implied by the assumption (2.22) for the evolution equation, in the model presented in this work. Nevertheless, following a different approach on the evolution condition, one could investigate further the possibility of imposing an inhomogeneous pressure trace upon Se, which is under our current investigation. For this reason, we do not imply this restriction on the results of Section 3.3, since they hold true up to the point that they are imposed properly in Section 3.4 and so after.Equation (3.32) is a nonlinear ordinary differential equation, which includes the rates (rn/rp) and (rq/rp), accompanied by rp, which are all unknown functions of time. In order to make (3.32) solvable, we need the expressions that connect all the unknown rates with function rp=rp(t). This is accomplished via expressions (3.15)–(3.17) and (3.23). In particular, (3.17) can be written as
(3.39)σn*-σp*rp2=γq3(rnrp)3(rprn-rprq)+γq6((rnrp)2-(rqrp)2),
or using the following simple symbols for the rates under consideration, those are (rn/rp)≡x=x(t), (rq/rp)≡y=y(t), and (rp/re)≡z=z(t), where by definition 0<x,y,z<1, expression (3.39) yields
(3.40)-γq3x3+γq2x2y-γq6y3-σn*-σp*rp2y=0.The third-degree nonlinear algebraic equation (3.40) provides x as a function of y and rp. Similarly, (3.23) assumes the form
(3.41)cq-pn3x3+(cp2-cq3)y3-[cp2+βp*rp2]y=0.Equations (3.40) and (3.41) form a nonlinear algebraic system, which when it is solved, provide the rates (rn/rp)≡x(t) and (rq/rp)=y(t) with respect to rp(t) and the differential (3.32) is solvable, under the initial condition rp(0)=Rp, which corresponds to the first observed tumour radius.Continuing, the rate(rp/re) is also connected to rp=rp(t) via (3.15), or equivalently (3.16); it is worthwhile to specify this connection, as it provides the evolution of the host boundary and also it appears in the nutrient field, the inhibitor field, and the pressure field provided in Sections 3.1, 3.2, and 3.3, respectively. In particular, using (3.15), we arrive at
(3.42)σn*-σ∞(re)rp2=γp-γq2(rqrp)2(1-23rqrprpre)-γp2(1-23rpre)+γq3(rnrp)3(rprn-rpre)+γq6(rnrp)2
or, in terms of the x,y,z symbols,
(3.43)(-γq3x3-γp-γq3y3+γp3)z=-γq2x2-γp-γq2y2+γp2+σn*-σ∞(re)rp2,
and if one substitutes the (x,y) solution of the system (3.40) and (3.41), then (3.43) provides the rate (rp/re)≡z(t) as a function of rp(t).


## 3.4.1. The Homogeneous Approach

Following the demand that the nutrient trace on the exterior boundarySe is homogeneous, the nutrient concentration (3.12) in the proliferating region is
(3.30)σp(r,θ,ϕ)=σ∞(re)+γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re).Substituting (3.30), (3.21), and (3.28) in (2.22) and performing algebraic manipulations the following evolution equation is obtained:
(3.31)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ).We note here that the whole model’s treatment is based on the assumption that as the tumour evolves, all its boundaries remain members of the concentric spherical coordinate system, so that all the Laplace and the Poisson partial differential equations can be treated analytically through the separation of variables method. This implies thatrp=rp(t), and consequently its derivative cannot depend on the direction r^≡(θ,ϕ). As a consequence, only the zeroth-order term of the spherical expansion can be included in (3.31), a fact that eliminates the effect of any pressure that is imposed to the host tissue by the exterior medium on the tumour’s evolution. Moreover, the host tissue’s boundary re=re(t) itself does not seem to be involved by the tumour’s evolution. Therefore, the evolution equation is written as follows:
(3.32)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3},
where the physical implementation was secured mathematically by imposing hlm,s=0 for l≥1, m=0,1,…,l, and s=e,o into (3.31) concluding to (3.32). This equation describes the transitory phase until it reaches its steady state. From (3.32), one can obtain the radii of the tumour’s steady-state boundaries by setting the right-hand side equal to zero.


## 3.4.2. The Nonhomogeneous Approach

By considering (3.14) as an average condition, the nutrient concentration in the proliferation shell reads as
(3.33)σp(r,θ,ϕ)=γp6r2-γp2rp2(1-23rpre)+[γq3rn3+γp-γq3rq3](1r-1re)+σ∞(re)∑l,m,sglm,s(rre)lYlm,s(θ,ϕ).Then, the corresponding evolution equation fordrp/dt becomes
(3.34)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}-μpP∞(re)rp∑l,m,shlm,sl(rpre)lYlm,s(θ,ϕ)+μσσ∞(re)rp∑l,m,sglm,sl(rpre)lYlm,s(θ,ϕ),
or
(3.35)drpdt=rp{(rnrp)3(μPFq-Fn3+μβcq-pn3-μσγq3)+(rqrp)3(μPFp-Fq3+μβcp-cq3-μσγp-γq3)-μPFp3-μβcp3+μσγp3}+1rp∑l,m,sl(rpre)l(μσσ∞(re)glm,s-μpP∞(re)hlm,s)Ylm,s(θ,ϕ).The arguments stated in Section3.4.1 hold true in this case as well, provided that between the spherical coefficients of the nutrient exterior field and the pressure field, as they are dictated on Se, the following condition is considered:
(3.36)glm,s=μpP∞(re)μσσ∞(re)hlm,sforl=1,2,…,m=1,2,…,lands=e,o,
and hence
(3.37)σ∞(re)=σ∞(re)-μpμσP∞(re)h00,e+μpμσP∞(re)
or
(3.38)μσ(σ∞(re)-σ∞(re))=μp(P∞(re)-P∞(re)h00,e).Consequently, the tumour’s boundaries evolve as concentric spheres and the evolution (3.35), provided (3.36), coincides with (3.32). From there after, it is obvious that the nonhomogeneous and the homogeneous approaches are treated in the same way.We stress that the restrictions on the type of inhomogeneity exhibited by the exterior pressure are implied by the assumption (2.22) for the evolution equation, in the model presented in this work. Nevertheless, following a different approach on the evolution condition, one could investigate further the possibility of imposing an inhomogeneous pressure trace upon Se, which is under our current investigation. For this reason, we do not imply this restriction on the results of Section 3.3, since they hold true up to the point that they are imposed properly in Section 3.4 and so after.Equation (3.32) is a nonlinear ordinary differential equation, which includes the rates (rn/rp) and (rq/rp), accompanied by rp, which are all unknown functions of time. In order to make (3.32) solvable, we need the expressions that connect all the unknown rates with function rp=rp(t). This is accomplished via expressions (3.15)–(3.17) and (3.23). In particular, (3.17) can be written as
(3.39)σn*-σp*rp2=γq3(rnrp)3(rprn-rprq)+γq6((rnrp)2-(rqrp)2),
or using the following simple symbols for the rates under consideration, those are (rn/rp)≡x=x(t), (rq/rp)≡y=y(t), and (rp/re)≡z=z(t), where by definition 0<x,y,z<1, expression (3.39) yields
(3.40)-γq3x3+γq2x2y-γq6y3-σn*-σp*rp2y=0.The third-degree nonlinear algebraic equation (3.40) provides x as a function of y and rp. Similarly, (3.23) assumes the form
(3.41)cq-pn3x3+(cp2-cq3)y3-[cp2+βp*rp2]y=0.Equations (3.40) and (3.41) form a nonlinear algebraic system, which when it is solved, provide the rates (rn/rp)≡x(t) and (rq/rp)=y(t) with respect to rp(t) and the differential (3.32) is solvable, under the initial condition rp(0)=Rp, which corresponds to the first observed tumour radius.Continuing, the rate(rp/re) is also connected to rp=rp(t) via (3.15), or equivalently (3.16); it is worthwhile to specify this connection, as it provides the evolution of the host boundary and also it appears in the nutrient field, the inhibitor field, and the pressure field provided in Sections 3.1, 3.2, and 3.3, respectively. In particular, using (3.15), we arrive at
(3.42)σn*-σ∞(re)rp2=γp-γq2(rqrp)2(1-23rqrprpre)-γp2(1-23rpre)+γq3(rnrp)3(rprn-rpre)+γq6(rnrp)2
or, in terms of the x,y,z symbols,
(3.43)(-γq3x3-γp-γq3y3+γp3)z=-γq2x2-γp-γq2y2+γp2+σn*-σ∞(re)rp2,
and if one substitutes the (x,y) solution of the system (3.40) and (3.41), then (3.43) provides the rate (rp/re)≡z(t) as a function of rp(t).


## 4. Special Cases

In this section, we implement the results obtained in Section3 in three special cases that correspond to different physical approaches. At the end, we reduce the results to recover the corresponding existing results that appear in the relative literature, mainly in [5, 6].


### 4.1. Inhomogeneous Data Imposed onSe

For the purpose of illuminating the inhomogeneous case, we suppose that the pressure’s trace onSe implies an inhomogeneity in the form
(4.1)P∞(re)=P∞(re)(1±ae(rpre)sin2θ)forre∈Se,
where ae is a positive proportionality constant. Expression (4.1) in terms of spherical harmonics assumes the form
(4.2)P∞(re)=P∞(re)(1±2ae3(rpre)∓2ae3(rpre)P2(cosθ)),
which in terms of (3.25) corresponds to
(4.3)h00,e=(1±2ae3(rpre)),h20,e=54π2∫02π∫0π(∓2ae3(rpre)P2(cosθ))Y20,e(θ,ϕ)sinθdθdϕ=∓2ae3(rpre),
and h20,o=0, h2m,s=0 for m=1,2 and s=e,o, while hlm,s=0 for every l=1,3,…, m=0,1,…,l and s=e,o. Then, the pressure field (3.26)–(3.29) enjoys a corresponding modification. In order for the evolution of the tumour’s exterior boundary to maintain angular independence, a solution of the model exists, when the nutrient supply is also nonhomogeneous and its spherical coefficients are defined by (3.36). Then,
(4.4)σ∞(re)=σ∞(re)+σ∞(re)∑l,m,sl≥1μpP∞(re)μσσ∞(re)hlm,sYlm,s(θ,ϕ)
or, for the case under consideration,
(4.5)σ∞(re)=σ∞(re)+μpμσP∞(re)(∓2ae3(rpre)P2(cosθ)).The evolution equation is not affected by this inhomogeneity, and it is given by (3.32).


### 4.2. Homogeneous Pressure Distribution Implying the Young-Laplace Law onSe

In this Section we consider the host nutritive environment to lie within an infinite medium, considered to be an incompressible fluid characterized by physical parameters different from those of both the host tissue and the tumour. Let us suppose that this medium imposes in the host tissue a uniform pressureP∞. Then, in the interface Se, the Young-Laplace equation should hold, which yields the boundary condition
(4.6)Pe(re)-P∞=aereforre∈Se,
where ae is a positive proportionality constant. Condition (4.6) corresponds to the results of the problem under investigation, if the trace of the pressure field on Se is interpreted as
(4.7)Pe(re)=P∞(re)≡P∞+aere.In terms of expression (2.20), which yields P∞(re)=P∞+(ae/re), we have hlm,s=0 for l≥1, m=0,1,…,l, and s=e,o, while h00,e=1. Then, the pressure field (3.26)–(3.29) would be affected as follows:
(4.8)Pn(r)=P∞+aere+aprp+Fq-Fn2rn2(1-23rnre)+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+Fn6r2,Pq(r)=P∞+aere+aprp+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+Fq-Fn3rn3(1r-1re)+Fq6r2,Pp(r)=P∞+aere+aprp-Fp2rp2(1-23rpre)+(Fq-Fn3rn3+Fp-Fq3rq3)(1r-1re)+Fp6r2,Pe(r)=P∞+aere+(Fq-Fn3rn3+Fp-Fq3rq3-Fp3rp3)(1r-1re).Though, the exterior pressure should not affect the evolution equation at all, as it is explained in Section3.4.


### 4.3. Homogeneous Nutrient Concentration Implying the Gibbs-Thomson Relation onSe

Following a related idea presented in [11], one can assume the Gibbs-Thomson relation for connecting the nutrient concentration trace at the point re∈Se with the corresponding curvature of Se, that is,
(4.9)σe(re)-σ∞=ae,σreforre∈Se,
where ae,σ is a positive proportionality constant. This jump condition is justified if we consider two facts. Firstly, as in Section 4.1, we consider that the host nutritive environment lies within an infinite medium of constant nutrient concentration σ∞. Secondly, we assume that the noncancerous cells at the boundary Se of the host organ consume nutrients in order to provide the energy necessary for maintaining the bonds between them and the solidity of the whole structure and this energy is curvature related as dictated by the Gibbs-Thomson relation (4.9). This case corresponds to the results obtained in the present paper if we imply
(4.10)σ∞(re)=σ∞(re)≡σ∞+ae,σre
inside (3.9).The nutrient field in this case is immediately derived from relationships (3.10)–(3.13), by substitution of (4.10), in the form
(4.11)σn(r)=σ∞+ae,σre+γp-γq2rq2(1-23rqre)+γq2rn2(1-23rnre)-γp2rp2(1-23rpre),σq(r)=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1r-1re)+γq6r2+σ∞+ae,σre,σp(r)=-γp2rp2(1-23rpre)+(γq3rn3+γp-γq3rq3)(1r-1re)+γp6r2+σ∞+ae,σre,σe(r)=(γq3rn3+γp-γq3rq3-γp3rp3)(1r-1re)+σ∞+ae,σre.The corresponding modifications upon (3.15) and (3.16), related to the nutrient critical values, are
(4.12)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞+ae,σre,σp*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rq-1re)+γq6rq2+σ∞+ae,σre.Since (3.15) and (3.16) provide the restrictions that interconnect the tumour’s boundaries, the evolution equation is affected by the assumption (4.10) through its effect upon (4.12).The above results fully recover the corresponding nutrient results obtained in [11] for the case where γp=γq and re=rp.


### 4.4. Recovery of Existing Results

In order to let the results obtained in Section3 to recover existing results in the literature [5, 6], we take the limit of our results where re→rp, at the case where γp=γq and σe(re)=σ∞. Then, expressions (3.10)–(3.13) become
(4.13)σn(r,θ,ϕ)=γq2rn2(1-23rnrp)-γp6rp2+σ∞,σq(r,θ,ϕ)=σp(r,θ,ϕ)=-γp6rp2+γp3rn3(1r-1rp)+γp6r2+σ∞.Moreover, if the inhibitor is considered to vanish far away from the tumour and not on the tumour’s boundary and also ifcp=cq=0, then the results (3.19)–(3.22) for the inhibitor concentration coincide with the case examined in [5] and yield
(4.14)βn(r)=-pn2rn2+pn6r2,βq(r)=βp(r)=-pn3rn3r.Then, (3.15) and (3.23) reveal
(4.15)σ∞-σn*=-γp3(12(rp2-rn2)-rn3(1rn-1rp)),βp*=-pn3rn3rq,
respectively, while they provide the inner boundaries rn, rq with respect to the exterior boundary rp of the tumour.The corresponding results for the pressure field when no exterior pressure is imposed on the tumor are found in [6], where the interior of the tumor is characterized by Fn=Fq=Fp and P∞(re)=0. Then, the expressions (3.26)–(3.28) meet the corresponding result in [6] and become
(4.16)Pn(r)=Pq(r)=Pp(r)=Fp6(r2-rp2)+aprp.The evolution equation studied in [6] coincides with the one included in the present work with the parameter values μp=1, μσ=(λ/μ), and μβ=0. The results cannot be compared further with [5] or [6] as they have been obtained either due to the use of different evolution equation or because the nutrient satisfies different conditions on the tumour’s boundary.


## 4.1. Inhomogeneous Data Imposed onSe

For the purpose of illuminating the inhomogeneous case, we suppose that the pressure’s trace onSe implies an inhomogeneity in the form
(4.1)P∞(re)=P∞(re)(1±ae(rpre)sin2θ)forre∈Se,
where ae is a positive proportionality constant. Expression (4.1) in terms of spherical harmonics assumes the form
(4.2)P∞(re)=P∞(re)(1±2ae3(rpre)∓2ae3(rpre)P2(cosθ)),
which in terms of (3.25) corresponds to
(4.3)h00,e=(1±2ae3(rpre)),h20,e=54π2∫02π∫0π(∓2ae3(rpre)P2(cosθ))Y20,e(θ,ϕ)sinθdθdϕ=∓2ae3(rpre),
and h20,o=0, h2m,s=0 for m=1,2 and s=e,o, while hlm,s=0 for every l=1,3,…, m=0,1,…,l and s=e,o. Then, the pressure field (3.26)–(3.29) enjoys a corresponding modification. In order for the evolution of the tumour’s exterior boundary to maintain angular independence, a solution of the model exists, when the nutrient supply is also nonhomogeneous and its spherical coefficients are defined by (3.36). Then,
(4.4)σ∞(re)=σ∞(re)+σ∞(re)∑l,m,sl≥1μpP∞(re)μσσ∞(re)hlm,sYlm,s(θ,ϕ)
or, for the case under consideration,
(4.5)σ∞(re)=σ∞(re)+μpμσP∞(re)(∓2ae3(rpre)P2(cosθ)).The evolution equation is not affected by this inhomogeneity, and it is given by (3.32).


## 4.2. Homogeneous Pressure Distribution Implying the Young-Laplace Law onSe

In this Section we consider the host nutritive environment to lie within an infinite medium, considered to be an incompressible fluid characterized by physical parameters different from those of both the host tissue and the tumour. Let us suppose that this medium imposes in the host tissue a uniform pressureP∞. Then, in the interface Se, the Young-Laplace equation should hold, which yields the boundary condition
(4.6)Pe(re)-P∞=aereforre∈Se,
where ae is a positive proportionality constant. Condition (4.6) corresponds to the results of the problem under investigation, if the trace of the pressure field on Se is interpreted as
(4.7)Pe(re)=P∞(re)≡P∞+aere.In terms of expression (2.20), which yields P∞(re)=P∞+(ae/re), we have hlm,s=0 for l≥1, m=0,1,…,l, and s=e,o, while h00,e=1. Then, the pressure field (3.26)–(3.29) would be affected as follows:
(4.8)Pn(r)=P∞+aere+aprp+Fq-Fn2rn2(1-23rnre)+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+Fn6r2,Pq(r)=P∞+aere+aprp+Fp-Fq2rq2(1-23rqre)-Fp2rp2(1-23rpre)+Fq-Fn3rn3(1r-1re)+Fq6r2,Pp(r)=P∞+aere+aprp-Fp2rp2(1-23rpre)+(Fq-Fn3rn3+Fp-Fq3rq3)(1r-1re)+Fp6r2,Pe(r)=P∞+aere+(Fq-Fn3rn3+Fp-Fq3rq3-Fp3rp3)(1r-1re).Though, the exterior pressure should not affect the evolution equation at all, as it is explained in Section3.4.


## 4.3. Homogeneous Nutrient Concentration Implying the Gibbs-Thomson Relation onSe

Following a related idea presented in [11], one can assume the Gibbs-Thomson relation for connecting the nutrient concentration trace at the point re∈Se with the corresponding curvature of Se, that is,
(4.9)σe(re)-σ∞=ae,σreforre∈Se,
where ae,σ is a positive proportionality constant. This jump condition is justified if we consider two facts. Firstly, as in Section 4.1, we consider that the host nutritive environment lies within an infinite medium of constant nutrient concentration σ∞. Secondly, we assume that the noncancerous cells at the boundary Se of the host organ consume nutrients in order to provide the energy necessary for maintaining the bonds between them and the solidity of the whole structure and this energy is curvature related as dictated by the Gibbs-Thomson relation (4.9). This case corresponds to the results obtained in the present paper if we imply
(4.10)σ∞(re)=σ∞(re)≡σ∞+ae,σre
inside (3.9).The nutrient field in this case is immediately derived from relationships (3.10)–(3.13), by substitution of (4.10), in the form
(4.11)σn(r)=σ∞+ae,σre+γp-γq2rq2(1-23rqre)+γq2rn2(1-23rnre)-γp2rp2(1-23rpre),σq(r)=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1r-1re)+γq6r2+σ∞+ae,σre,σp(r)=-γp2rp2(1-23rpre)+(γq3rn3+γp-γq3rq3)(1r-1re)+γp6r2+σ∞+ae,σre,σe(r)=(γq3rn3+γp-γq3rq3-γp3rp3)(1r-1re)+σ∞+ae,σre.The corresponding modifications upon (3.15) and (3.16), related to the nutrient critical values, are
(4.12)σn*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rn-1re)+γq6rn2+σ∞+ae,σre,σp*=γp-γq2rq2(1-23rqre)-γp2rp2(1-23rpre)+γq3rn3(1rq-1re)+γq6rq2+σ∞+ae,σre.Since (3.15) and (3.16) provide the restrictions that interconnect the tumour’s boundaries, the evolution equation is affected by the assumption (4.10) through its effect upon (4.12).The above results fully recover the corresponding nutrient results obtained in [11] for the case where γp=γq and re=rp.


## 4.4. Recovery of Existing Results

In order to let the results obtained in Section3 to recover existing results in the literature [5, 6], we take the limit of our results where re→rp, at the case where γp=γq and σe(re)=σ∞. Then, expressions (3.10)–(3.13) become
(4.13)σn(r,θ,ϕ)=γq2rn2(1-23rnrp)-γp6rp2+σ∞,σq(r,θ,ϕ)=σp(r,θ,ϕ)=-γp6rp2+γp3rn3(1r-1rp)+γp6r2+σ∞.Moreover, if the inhibitor is considered to vanish far away from the tumour and not on the tumour’s boundary and also ifcp=cq=0, then the results (3.19)–(3.22) for the inhibitor concentration coincide with the case examined in [5] and yield
(4.14)βn(r)=-pn2rn2+pn6r2,βq(r)=βp(r)=-pn3rn3r.Then, (3.15) and (3.23) reveal
(4.15)σ∞-σn*=-γp3(12(rp2-rn2)-rn3(1rn-1rp)),βp*=-pn3rn3rq,
respectively, while they provide the inner boundaries rn, rq with respect to the exterior boundary rp of the tumour.The corresponding results for the pressure field when no exterior pressure is imposed on the tumor are found in [6], where the interior of the tumor is characterized by Fn=Fq=Fp and P∞(re)=0. Then, the expressions (3.26)–(3.28) meet the corresponding result in [6] and become
(4.16)Pn(r)=Pq(r)=Pp(r)=Fp6(r2-rp2)+aprp.The evolution equation studied in [6] coincides with the one included in the present work with the parameter values μp=1, μσ=(λ/μ), and μβ=0. The results cannot be compared further with [5] or [6] as they have been obtained either due to the use of different evolution equation or because the nutrient satisfies different conditions on the tumour’s boundary.


## 5. Discussion

In the present work, we have studied a continuous model of avascular tumour growth, under two basic assumptions. Firstly, it evolves maintaining a spherical multilayer structure, lying inside a finite concentric spherical host medium. Secondly, its evolution is regulated by the diffusion of an inhomogeneous nutrient field and of a growth inhibitory factor and by the pressure that results from the cell proliferation and disintegration, as well as from an externally imposed inhomogeneous pressure from the finite surrounding medium. All the structure interfaces evolve in time, but the time scale of the space expansion of the tumour is so much different than the diffusion time scales of the substances in the tumour that the equilibrium assumption is well justified.Hence, the model is formulated in three boundary value problems that hold true as the tumour evolves and provide the nutrient field, the internally produced inhibitor field and the pressure field throughout the tumour and the host surrounding. The model includes an assumption for the boundaries’ evolution, which is formulated as a nonlinear ordinary differential equation with respect to the tumour’s exterior boundary, and also it includes connection formulae between all the other boundaries with respect to the tumour’s exterior one.The same formulation and analysis can probably be adapted in different models with alternative interpretations, for example, the nutrient concentration can be considered as a drug substance supplied externally that inhibits the growth, in which case the evolution equation would be respectively modified and so on.In the present work we focus on the type of inhomogeneity that can be compatible with the specific model that we present. The inhomogeneities that we are focused on concern the data that are independent of the tumour’s development, that is, the nutrient supply and the pressure field imposed by the exterior to the host organ surrounding. It turns out that a concentric spherical multilayer development could be secured under homogeneous nutrient supply and pressure field, or by a nutrient supply and an exterior pressure that exhibit a special connection between them, expressed in (3.36). In particular, the assumption that the spherical interfaces that separate the tumour regions are characterized by the critical nutrient values, demands that the nutrient supply should be homogeneous. Moreover, the assumption that the interfaces evolve maintaining their concentric spherical shape, together with the specific assumption for the tumour’s evolution, impose, the direction independence for the pressure trace on the exterior host tissue surface. On the other hand, relaxing the first demand and allowing the interfaces to be perturbed, the nutrient supply can be inhomogeneous and consequently the pressure trace can exhibit a type of inhomogeneity that secures the angular independence of the evolution equation.These considerations are illuminated, by implementing the general results in special cases, in Section4. In particular, it is shown that an angular-dependent pressure field of the form (4.1) permits a concentric spherical tumour development only under a similarly nonhomogeneous nutrient distribution in the surrounding, given in (4.5). A homogeneous pressure of the form (4.7), that takes into account the exterior boundary’s curvature and is dictated by the Young-Laplace law of different fluids’ interface, is shown to have no effect on the evolution equation of the tumour’s boundary and also to be irrelevant of the form of the nutrient supply. On the other hand, a homogeneous nutrient field that takes into account the exterior spherical curvature, dictated by the Gibbs-Thomson relation, has an indirect effect on the tumour’s evolution, as it effects the connection formulae (4.12), between the tumour’s interfaces. Though, as the two latter cases refer to homogeneous imposed fields, the interrelation (3.36), in view of (3.9) and (3.25), does not exist, as it is expected. Finally, the reduction to the simple cases of both homogeneous fields in an infinite host surrounding exhibits the lack of interdependence between the imposed pressure and the nutrient supply and fully recovers the existing radial results from the literature, as mentioned in Section 4.4. The consideration of the model presented in Section 3 under shape perturbations is not examined in the present research, as the focus of this work is on studying the consistency of the spherical development with nonhomogeneous imposed data. However, the sensitivity analysis of the proposed model is within our future consideration.Clearly, it is a drawback of the model presented, that it incorporates a rather ideal concentric spherical structure and the results refer to such an idealistic configuration or simplified assumptions. Alternative evolution approaches for the same spherical structure in avascular tumour development, and also alternative geometrical structure of the development, which is much more applicable to cancer growing in humans, are under our current investigation.



---

*Source: 175434-2012-09-06.xml*