# Tumor Growth with Nutrients: Regularity and Stability

Matt Jacobs  Department of Mathematics, Purdue University, West Lafayette, IN
<EMAIL_ADDRESS>,  Inwon Kim  Department of Mathematics, UCLA, Los
angeles, CA<EMAIL_ADDRESS>and  Jiajun Tong  Beijing International Center
for Mathematical Research, Peking University, Beijing, China
<EMAIL_ADDRESS>
###### Abstract.

In this paper we study a tumor growth model with nutrients. The model presents
dynamic patch solutions due to the contact inhibition among the tumor cells.
We show that when the nutrients do not diffuse and the cells do not die, the
tumor density exhibits regularizing dynamics. In particular, we provide
contraction estimates, exponential rate of asymptotic convergence, and
boundary regularity of the tumor patch. These results are in sharp contrast to
the models either with nutrient diffusion or with death rate in tumor cells.

I. Kim is partially supported by NSF Grant DMS-1900804. J. Tong is partially
supported by the Peking University Start-up Fund.

##  1\. Introduction

A model system that appears in literature describing tumor growth with
nutrients is

(P) $\left\\{\begin{array}[]{l}\rho_{t}-\nabla\cdot(\rho\nabla
p)=(n-b)\rho,\quad\rho\leq 1,\quad p\in P_{\infty}(\rho),\\\ \\\ n_{t}-D\Delta
n=-\rho n,\quad n\to c>0\mbox{ as }|x|\to\infty\end{array}\right.$

set in $Q:=\mathbb{R}^{d}\times[0,\infty)$, with $b$, $D$ and $c$ being non-
negative constants (see e.g. [MRCS14, PQV14, PTV14, DP21]). It is equipped
with the initial conditions

(1.1) $\rho(x,0)=\rho_{0}(x),\quad n(x,0)=n_{0}(x).$

Here $\rho$ and $n$ respectively denote the density of tumor cells and the
nutrients. The cells grow by consuming the nutrients which are supplied from
the external environment, while they die at a constant rate $b$ in the
meantime. The pressure $p\geq 0$ can be understood as the Lagrange multiplier
for the constraint $\rho\leq 1$ that represents the contact inhibition in
cells. In (P), $P_{\infty}(\rho)$ denotes the Hele-Shaw graph

$P_{\infty}(\rho)=\begin{cases}0,&\mbox{if }\rho\in[0,1),\\\
[0,+\infty),&\mbox{if }\rho=1.\end{cases}$

It is well-known that the resulting solution features time-evolving patches of
congested cell region where $\rho$ equals $1$. We call this set a tumor patch
for later reference.

This model, while relatively simple, presents a complex phenomena. One can
view the system as a singular limit of a reaction-diffusion system, where one
first takes $p=\rho^{m-1}$ and then sends $m$ to $+\infty$; see [PQV14, DP21].
With finite $m$, the reaction-diffusion system has been actively studied in
the literature: see [KMM+97, Kit97, Mim04].

The well-posedness of (P) is not hard to achieve: see Section 2 for more
discussions. On the other hand, qualitative behavior of the solutions of (P)
is much less understood. In particular, the growth of the tumor patch appears
to generate fingering phenomena, as observed by numerical experiments [Kit97,
MRCS14, PTV14] even when the patch is almost radial and when $n$ is initially
a constant. Such “dendric growth” is well-known to persist in models of
bacterial growth that aggressively consumes nutrients [BJST+94]. When $D>0$,
the formation of dendritic patterns is conjectured to occur because there are
more nutrients available near the tips of dendritic fingers compared to
valleys [MRCS14] (in the valleys there are more surrounding bacteria to
consume the nutrients). As a result, dendritic tips grow faster than the
valleys, leading to the amplification of instabilities. When $D=0$ and $b>0$,
numerical experiments still observe the fingering phenomena, possibly due to
the movement of tumor cells toward the necrotic core, where cells decay from
the maximal density (c.f. Figure 1). The scale of the aforementioned
instabilities in terms of $b$ and $D$ remains to be understood. While we do
not pursue regularity analysis in such general cases, a variational scheme is
introduced to approximate (P) and yield well-posedness in general settings:
see Theorem 2.2. Although there are many other possible ways to approximate
the equation and obtain the well-posedness, the variational scheme we
introduce is particularly numerically efficient and preserves many desirable
properties of the true system (c.f. the discussion in Section 2.1).

Figure 1. Numerical simulations from [JL22] of the tumor growth system with
$D=0$ and $b=0.4$ and $n_{0}=2$. The first image represents the initial patch
density. Brighter colored pixels indicate larger density values. The evolution
shows dendritic growth at the boundaries once cells begin to die in the
interior.

The main goal of this paper is to study the singular case $b=D=0$, with
particular focus on the dynamics and regularity of the tumor patch evolution.
This case is particularly interesting, as numerical experiments in [MRCS14]
suggest that, when $b$ is fixed to be zero, the dendritic behavior becomes
more and more branched and irregular as the diffusion parameter $D$ becomes
smaller and smaller. Surprisingly, we find that once the diffusion parameter
is set to zero, there is no dendritic growth whatsoever and the evolution is
regularizing. Roughly speaking, we will show that there is no rough growth of
the tumor patch other than those caused by topological changes (Theorem 2.7).
Moreover, when $n$ is initially a constant, we can further show that the
dynamics of the tumor patch can be understood in terms of a single-variable
nutrient-free system (Theorem 2.8). We hope our findings serve as the first
step to understanding the complex behavior of the system (P) with general
values of $b$ and $D$. In particular, reconciling our results with the
numerical experiments in [MRCS14] would be a very interesting future direction
of study.

Our results are based on the following rather unexpected comparison principle
when $b=D=0$ (Proposition 4.1), for two different solutions
$(\rho^{i},p^{i},n^{i})$ $(i=0,1)$ of (P):

$\mbox{ Let }\eta^{i}:=n^{i}_{0}-n^{i}.\mbox{ If }\rho^{1}\leq\rho^{0}\mbox{
and }n^{1}\leq n^{0}\mbox{ at }t=0\mbox{, then }\rho^{1}\leq\rho^{0}\mbox{ and
}\eta^{1}\leq\eta^{0}\mbox{ for all }t\geq 0.$

This comparison property is somewhat unintuitive, as the smaller density
should have more nutrients available for growth, raising the posibility that
the ordering could be violated as the system evolves. As such, observe that
this is not a standard comparison principle. While the initial ordering is
with the nutrient variable, the ordering at later times are associated with
the $\eta$-variable. To see why this should hold, we introduce the variable
$w(x,t):=\int_{0}^{t}p(x,s)\,ds$. Integrating the $\rho$-equation of (P) in
time and noticing that $n\rho=-n_{t}=\eta_{t}$ when $D=0$, we obtain

(1.2) $\rho-\Delta w=\rho_{0}+\int_{0}^{t}(n\rho)\,ds=\rho_{0}+\eta,$

and

(1.3) $\eta_{t}=(n_{0}-\eta)\rho.$

The time integrated system reveals that the total growth of the density at
time $t$ only depends on $\eta$, the total amount of nutrients consumed by the
time $t$, rather than the amount of available nutrients at time $t$. Thus, by
working with this version of the problem, the possibility of a comparison
property becomes much more evident. Indeed, a large part of our analysis,
including the proof of the above comparison principle, will be based on this
time-integrated system (1.2)-(1.3).

Heuristically speaking, since $\rho$ only grows and the set $\\{\rho=1\\}$
expands in time, the pressure variable is positive in the growing parts of
$\\{\rho=1\\}$. More precisely, the set $\\{\rho=1\\}$ can be decomposed as
the union $\\{w>0\\}\cup\\{\rho_{0}=1,n_{0}=0\\}$ (c.f. Lemma 4.6). On the
other hand, $\rho=\rho_{0}+\eta$ in $\\{w=0\\}$. Thus, $(w,\eta)$ satisfies

$\left\\{\begin{array}[]{l}(1-\eta-\rho_{0})\chi_{\\{w(\cdot,t)>0\\}}-\Delta
w(\cdot,t)=0\quad\hbox{in }\mathbb{R}^{d},\\\ \\\
\eta_{t}=(n_{0}-\eta)\rho.\end{array}\right.$

It is then not difficult to see that the aforementioned comparison principle
holds for this system. Since the above system is in the form of an obstacle
problem, let us discuss the free boundary regularity of the set $\\{w>0\\}$.
The standard theory for the obstacle problem yields that, as long as $\eta$ is
less than $1$ near the boundary of $\\{w>0\\}$ and is $C^{\alpha}$, the
boundary of $\\{w>0\\}$ has $C^{1,\alpha}$-regularity, away from cusp-type
singular points [Bla01, Caf98]. Thus, when we study the boundary regularity of
the tumor patches, an important step in our analysis is to ensure the
regularity of $n$, the non-degeneracy of $w$, and the nonexistence of cusp
points on the patch boundary.

We will see in Theorem 4.10 and Corollary 5.15 that, if $n_{0}(x)\geq 1$, the
pressure dominates the evolution and there is a generic (local-in-time)
regularity of the patch boundary after some finite time, whereas for
$n_{0}(x)<1$, the pressure diminishes exponentially fast as the nutrient
vanishes, and the tumor patch approaches an asymptotic profile. The finite-
time regularization result for $n_{0}\geq 1$ features similarity to those of
porous medium equations [CVW87] and the Hele-Shaw flow [Kim06]. In the case of
$n_{0}<1$, the large-time regularity of the tumor patch remains open in
general, and it may be possible that fractal structures persist as the set
approaches its asymptotic profile. Nonetheless, if the asymptotic tumor patch
is sufficiently far away from the convex hull of its initial position, we can
show that the patch evolution turns smooth within finite time.

Our last main result focuses on the problem when the initial nutrient is
constant, namely when $n_{0}(x)\equiv c$. In this case, very surprisingly, one
can characterize the behavior of the tumor patch solution in (P) through what
we call _master dynamics_. There are two of them, stated in Proposition 6.6
and Proposition 6.7 respectively. Let us take the first one as an example.
Consider the following problem that concerns the density evolution only:

(HS) $\partial_{t}\rho_{*}-\nabla\cdot(\rho_{*}\nabla
p_{*})=\rho_{0},\quad\rho_{*}\leq 1,\quad p_{*}\in
P_{\infty}(\rho_{*}),\quad\rho_{*}|_{t=0}=\rho_{0}.$

$\rho_{0}$ is the initial data in (1.1) which is assumed to be a patch here.
(HS) is reminiscent of the classic Hele-Shaw flow, whose regularity and long-
time convergence property is relatively well-understood in various settings
[EJ81, Kim06, CJK07]. We will show in Proposition 6.6 that, once $n_{0}=c>0$
is given, the corresponding $\rho$-evolution in (P) is simply a re-scaled (in
time) version of the $\rho_{*}$-evolution in (HS), and the re-scaling depends
on $n_{0}$ explicitly. The $n$-evolution in (P) can be readily represented as
well. In other words, once we understand the density evolution in (HS), which
is parameter-free and much simpler than (P), we can fully characterize all the
$(\rho,n)$-evolutions in (P) corresponding to all different values of
$n_{0}>0$. This is why we call $\\{\rho_{*}(\cdot,t)\\}_{t\geq 0}$ the master
dynamics. The second master dynamics is proposed in a similar spirit; see
Proposition 6.7. Given that the value of $n_{0}$ has a non-trivial impact on
the patch solution dynamics in (P), and that $n$ does not stay as a constant
over time, this connection is far from being apparent. We should emphasize
that neither of the master dynamics can be obtained by a mere change of space
and time variables in (P), and even the relation between the two master
dynamics is also highly non-trivial. Their proof relies on a very special
property of the system (P) when $b=D=0$, $n_{0}(x)$ is a constant, and $\rho$
is a patch solution (see Lemma 6.1): for any harmonic function, its average on
the set $\\{\rho=1\\}$ is time-invariant.

The rest of the paper is organized as follows. In Section 2, we state our main
results and discuss their implications. In Section 3, we introduce a discrete-
in-time variational scheme to approximate (P) in the general setting of
$b,D\geq 0$ and prove well-posedness. Sections 4-6 are devoted to the case
$b=D=0$. In Section 4, we first study the elliptic equation (1.2) (see a more
precise formulation in (4.2)), and then prove general properties of the time-
integrated system when the initial nutrient is bounded. Section 5 is focused
on free boundary regularity of the tumor patch solutions under suitable
geometric assumptions. Finally, in Section 6, we prove the master dynamics
when $n_{0}$ is constant. Then we apply that to characterize long-time
behavior and uniform boundary regularity of the patch solution $\rho$.

##  2\. Summary of Main Results

###  2.1. Well-posedness in the general setting

We first define weak solutions of (P) as follows.

######  Definition 2.1.

Let $n_{0}\in L^{\infty}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ and
$\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ such that
$\rho_{0}\in[0,1]$ almost everywhere. Fix $T>0$ and denote
$Q_{T}:=\mathbb{R}^{d}\times[0,T]$. Non-negative functions $\rho\in
C([0,T];L^{1}(\mathbb{R}^{d}))\cap L^{\infty}([0,T];BV(\mathbb{R}^{d}))$,
$p\in L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, and $n\in L^{\infty}(Q_{T})\cap
L^{\infty}([0,T];BV(\mathbb{R}^{d}))$ on $Q_{T}$ are said to form a weak
solution of (P) in $Q_{T}$, if they satisfy:

  1. (i)

$\rho\in[0,1]$ and $p(1-\rho)=0$ in $Q_{T}$;

  2. (ii)

For any $\psi\in H^{1}(Q_{T})$ that vanishes at $t=T$, we have

$\int_{0}^{T}\int_{\mathbb{R}^{d}}\nabla\psi\cdot\nabla
p-\rho\partial_{t}\psi\,dx\,dt=\int_{\mathbb{R}^{d}}\psi(x,0)\rho_{0}\,dx+\int_{0}^{T}\int_{\mathbb{R}^{d}}\psi(n-b)\rho\,dx\,dt;$

  3. (iii)

In addition,

$\partial_{t}n-D\Delta n=-\rho n\hbox{ in }\mathcal{D}^{\prime}(Q_{T}),\quad
n(\cdot,0)=n_{0}.$

Here and in what follows, by saying a function $f\in BV(\mathbb{R}^{d})$, we
mean that the total variation of $f$ is finite, and yet $f$ may not be
integrable on $\mathbb{R}^{d}$.

Using an adaptation of the minimizing movements scheme introduced in [JKT21],
we establish the following well-posedness result.

######  Theorem 2.2 (Well-posedness, a summary of Propositions 3.6, 3.7, and
3.13 and Remark 3.12).

Let $n_{0}\in L^{\infty}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$. Let
$\rho_{0}\in BV(\mathbb{R}^{d})$ be compactly supported, such that
$\rho_{0}\in[0,1]$ almost everywhere. Then for given $b,D\geq 0$, and any
$T>0$, there exists a unique weak solution $(\rho,p,n)$ of (P) in
$Q_{T}=\mathbb{R}^{d}\times[0,T]$ in the sense of Definition 2.1.

Moreover, $\rho$ and $p$ are compactly supported in $Q_{T}$, and satisfy the
complementarity relation in the distribution sense:

$p(\Delta p+n-b)=0\hbox{ in }\mathcal{D}^{\prime}(Q_{T}).$

If, additionally, $\rho_{0}(x)\in\\{0,1\\}$ almost everywhere in space, then
for every time $t\in[0,T]$ we have $\rho\in\\{0,1\\}$ almost everywhere in
space.

As we mentioned above, we obtain this result via a variational approximation
scheme, which is a simplified version of the more general scheme introduced in
[JKT21]. Although there are many different ways to establish the well-
posedness of the system (such as the degenerate diffusion approach considered
in [PQV14] and [GKM22]), we emphasize our variational scheme, as it has a very
efficient numerical implementation via the Back-and-Forth method [JL20,
JLL21]. For instance, the images displayed in Figure 1 are computed on a high-
resolution $1024\times 1024$ grid. Carrying out simulations of this size has
been out of reach for previous methods, owing to the difficult nonlinearities
in the equation. The efficiency of the scheme will be studied further in the
upcoming paper [JL22]. In addition to the favorable numerics, the scheme also
preserves many desirable properties of the true system, such as the patch-
preserving property and various important estimates that will prove useful in
our analysis of the scheme.

###  2.2. Contraction and stability estimates

The rest of the results focus on the case $b=D=0$. In this setting, it will be
useful to introduce the quantity $\eta:=n_{0}-n$ which represents the amount
of nutrient that has been consumed. It solves

$\partial_{t}\eta=(n_{0}-\eta)\rho,\quad\eta(x,0)=0,$

which can be derived from (P), and which is equivalent to the $n$-equation.
Therefore, in what follows, in the case of $b=D=0$, we will use $(\rho,p,n)$
and $(\rho,p,\eta)$ interchangeably as the solution of (P). In this case, the
system (P) enjoys a lot of nice properties, especially for patch solutions,
i.e., $\rho_{0}$ and $\rho$ take the value $0$ or $1$ almost everywhere in
space at all times.

We first state the following $L^{1}$-stability estimate for patch solutions,
from which the comparison principle can be readily derived.

######  Theorem 2.3 ($L^{1}$-contraction, a simplified version of Theorem
4.7).

Suppose that $(\rho^{i},p^{i},\eta^{i})$ $(i=0,1)$ are weak solutions of (P),
starting from initial datum $(\rho^{i}_{0},n^{i}_{0})$ respectively. If
$\rho_{0}^{i}$ are patches and $n^{0}_{0}\in L^{\infty}(\mathbb{R}^{d})$, then

$\lVert(\rho^{1}(\cdot,t)-\rho^{0}(\cdot,t))_{+}\rVert_{L^{1}}\leq
N(t)\lVert(n^{1}_{0}-n^{0}_{0})_{+}\rVert_{L^{1}}+M(t)\lVert(\rho_{0}^{1}-\rho_{0}^{0})_{+}\rVert_{L^{1}}.$

Here $L^{r}=L^{r}(\mathbb{R}^{d})$,

(2.1) $N(t)=\begin{cases}\frac{e^{(\lVert
n^{0}_{0}\rVert_{L^{\infty}}-1)t}-1}{\lVert
n^{0}_{0}\rVert_{L^{\infty}}-1}&\textup{if}\;\lVert
n^{0}_{0}\rVert_{L^{\infty}}\neq 1,\\\ t&\textup{otherwise},\end{cases}$

and

(2.2) $M(t)=\begin{cases}\frac{\lVert n^{0}_{0}\rVert_{L^{\infty}}e^{(\lVert
n^{0}_{0}\rVert_{L^{\infty}}-1)t}-1}{\lVert
n^{0}_{0}\rVert_{L^{\infty}}-1}&\textup{if}\;\lVert
n^{0}_{0}\rVert_{L^{\infty}}\neq 1,\\\ t+1&\textup{otherwise}.\end{cases}$

Based on the $L^{1}$-contraction and symmetries of the Laplacian operator, we
can derive the following $BV$-estimates for $\rho$. Note that here we consider
a somewhat non-standard norm that we call the $BV_{A}$-norm, where $A$ is an
antisymmetric matrix. The $BV_{A}$-norm is the $L^{1}$-norm of the product
$\nabla\rho\cdot Ax$. For antisymmetric matrices $A$, $\nabla\rho\cdot Ax$
picks up a non-radial component of the derivative of $A$, with a stronger
weight as one moves further from the origin. By bounding this norm, we obtain
a stronger control (compared to the vanilla $BV$-norm) on the behavior of the
non-radial components of the boundary variation as the tumor grows.

######  Theorem 2.4 ($BV$-estimates, Proposition 4.9).

Given an antisymmetric matrix $A\in\mathbb{R}^{d\times d}$ and some $g\in
L^{1}(\mathbb{R}^{d})$, define

$\lVert g\rVert_{BV_{A}(\mathbb{R}^{d})}:=\sup_{\varphi\in
C^{\infty}_{c}(\mathbb{R}^{d}),\lVert\varphi\rVert_{L^{\infty}(\mathbb{R}^{d})}\leq
1}\int_{\mathbb{R}^{d}}g(x)Ax\cdot\nabla\varphi\,dx.$

Let $(\rho,p,\eta)$ be a weak solution of (P) with initial data
$(\rho_{0},n_{0})$. If $\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap
BV(\mathbb{R}^{d})$ is a patch, and $n_{0}\in L^{\infty}(\mathbb{R}^{d})\cap
BV(\mathbb{R}^{d})$, then, with $N(t)$ and $M(t)$ given above in (2.1) and
(2.2),

$\displaystyle\lVert\rho(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\leq$
$\displaystyle\;N(t)\lVert
n_{0}\rVert_{BV(\mathbb{R}^{d})}+M(t)\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})},$
$\displaystyle\lVert\rho(\cdot,t)\rVert_{BV_{A}(\mathbb{R}^{d})}\leq$
$\displaystyle\;N(t)\lVert
n_{0}\rVert_{BV_{A}(\mathbb{R}^{d})}+M(t)\lVert\rho_{0}\rVert_{BV_{A}(\mathbb{R}^{d})}.$

From the comparison principle as well as analysis on radial solutions, one can
conclude the following result on the long-time behavior of the solution when
$\|n_{0}\|_{L^{\infty}}<1$. In this case, the tumor eventually stops growing
and approaches a stationary solution.

######  Theorem 2.5 (Theorem 4.10).

Suppose that $\lVert n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}<1$ and
$\rho_{0}(x)\in\\{0,1\\}$. Let $(\rho,p,\eta)$ be as given above. Then

$\displaystyle\lVert\rho(\cdot,t)-\rho_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq$
$\displaystyle\;\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t},$
$\displaystyle\lVert\eta(\cdot,t)-n_{0}\rho_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq$
$\displaystyle\;\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t},$

where $\rho_{\infty}$ solves the elliptic equation

$(1-n_{0})\rho_{\infty}-\Delta w_{\infty}=\rho_{0},\quad
w_{\infty}(1-\rho_{\infty})=0.$

###  2.3. Geometry and regularity of the patch boundary

Under suitable assumptions on the initial patch data $\rho_{0}$ and the
initial nutrient $n_{0}$, we can study boundary regularity of the growing part
of the set $\\{\rho(\cdot,t)=1\\}$. We will use reflection invariance of the
problem and comparison principle, using reflection-based geometry of the sets,
to achieve these. Such an argument was used first in [FK14] and later in
[KK20, KKP21] to obtain regularity results for interface motions with
reflection invariance.

The following version of the comparison principle plays an important role.

######  Theorem 2.6 (Reflection comparison, Proposition 5.4).

Suppose that $(\rho,p,\eta)$ is a solution of (P) starting from the initial
data $(\rho^{0},n^{0})$. Given a hyperplane $H$, let
$\rho_{H},\rho^{0}_{H},\eta_{H}$, and $n^{0}_{H}$ denote the reflections of
$\rho,\rho^{0},\eta$, and $n^{0}$ about the hyperplane $H$ respectively. Let
$H^{+}$ be one of the half spaces generated by $H$. If
$\rho^{0}_{H}\leq\rho^{0}$ and $n^{0}_{H}\leq n^{0}$ almost everywhere in
$H^{+}$, then $\rho_{H}\leq\rho$ and $\eta_{H}\leq\eta$ almost everywhere in
$H^{+}$.

We say a set $S\subset\mathbb{R}^{d}$ satisfies $r$-reflection if it contains
the ball $B_{r}(0)=\\{|x|<r\\}$, and if for any hyperplane $H$ not
intersecting with $B_{r}(0)$, the reflected image of $S$ with respect to $H$
is a subset of $S$ when restricted on the side of $H$ that contains the
origin; see Definition 5.1. It is a notion that is stronger than being star-
shaped, tailored to work with the reflection comparison above. With this
concept, we obtain the following results based on the above theorem.

######  Theorem 2.7 (Corollary 5.15).

Suppose $n_{0}$ is $C^{1}$, and its super-level sets satisfy $r$-reflection
for some $r>0$. Let $\Omega_{0}$ be an open bounded set in $\mathbb{R}^{d}$
contained in $B_{r}(0)$, and $\rho_{0}=\chi_{\Omega_{0}}$. Let $T(R)$ be
defined as in Corollary 5.5. Then the followings hold for any $0<\alpha<1$ and
a dimensional constant $c_{d}$:

  1. (a)

If $n_{0}(x)\geq 1$ on $\mathbb{R}^{d}$, then for any $R>c_{d}r$,
$\partial\\{w(\cdot,t)>0\\}$ is uniformly $C^{1,\alpha}$ in a unit
neighborhood for any finite time range within $[T(R),\infty)$.

  2. (b)

If $n_{0}(x)<1$ on $\mathbb{R}^{d}$, the same holds in for any finite time
range within $[T(R),\infty)$ if $R>c_{d}r$ satisfies
$B_{R}(0)\subset\\{w_{\infty}>0\\}$. Here $w_{\infty}$ is defined in Theorem
2.5 above (or equivalently, Theorem 4.10).

In both cases of the above theorem, since we can start with
$\rho_{0}=\chi_{\Omega_{0}}$ where $\Omega_{0}$ is quite arbitrary, the
evolution of the set $\Omega_{t}:=\\{\rho(\cdot,t)=1\\}$ may go through
topological singularities such as merging of the free boundaries. However, the
above results state that, if the initial nutrient $n_{0}$ is “well-prepared”
outside of $B_{r}$, then after a finite time there is no further topological
changes in the evolution, and $\Omega_{t}$ evolves with smooth boundary
outside of $B_{r}$.

###  2.4. The case of constant $n_{0}$

When $n_{0}$ is constant, even stronger characterization can be provided for
the evolution of patch solutions $\rho$. In particular, the $\rho$-evolution
in (P) coincides with density evolutions in some nutrient-free and parameter-
free Hele-Shaw-type systems, up to explicit rescaling in space and time. Such
very surprising relation cannot be derived from trivial change of variables,
but it crucially relies on special properties of the system (P) under the
given assumptions.

######  Theorem 2.8 (Two master dynamics, Lemma 6.1, Proposition 6.6 and
Proposition 6.7).

Let $(\rho,p,n)$ solve (P) with $\rho_{0}(x)\in\\{0,1\\}$ almost everywhere
and $n_{0}(x)=n_{0}>0$ being constant. Then $n$ can be represented in terms of
$\rho$ using (6.2). Let $m(t)$ be explicitly defined by (6.4). Then the total
mass of the tumor satisfies

$\int_{\mathbb{R}^{d}}\rho(x,t)\,dx=m(t)\int_{\mathbb{R}^{d}}\rho_{0}(x)\,dx.$

Moreover,

  1. (1)

Let $(\rho_{*},p_{*})$ be a weak solution of

(HS) $\partial_{t}\rho_{*}-\nabla\cdot(\rho_{*}\nabla
p_{*})=\rho_{0},\quad\rho_{*}\leq 1,\quad p_{*}\in
P_{\infty}(\rho_{*}),\quad\rho_{*}|_{t=0}=\rho_{0}.$

Then $\rho(x,t)=\rho_{*}(x,m(t)-1)$ for all $t\geq 0$.

  2. (2)

Let $(\rho_{\dagger},p_{\dagger})$ be a weak solution of

(HS’)
$\partial_{t}\rho_{\dagger}-\nabla\cdot\Big{[}\rho_{\dagger}\nabla\big{(}p_{\dagger}+V(x)\big{)}\Big{]}=0,\quad\rho_{\dagger}\leq
1,\quad p_{\dagger}\in
P_{\infty}(\rho_{\dagger}),\quad\rho_{\dagger}|_{t=0}=\rho_{0},$

where $V(x)=\frac{|x|^{2}}{2d}$. Then
$\rho(x,t)=\rho_{\dagger}(m(t)^{-\frac{1}{d}}x,\,\ln m(t))$ for all $t\geq 0$.

Here the notion of weak solutions is provided in Definition 3.9 below.

Thanks to the master dynamics, we can characterize the long-time behavior of
the patch solutions when $n_{0}>0$ is constant in $\mathbb{R}^{d}$ (c.f.
Theorem 2.5) under a suitable rescaling.

######  Theorem 2.9 (Proposition 6.8).

Suppose $n_{0}>0$ is constant in $\mathbb{R}^{d}$. Let $m(t)$ be defined in
(6.4). Assume $\Omega_{0}$ to be a bounded open set, such that
$B_{r_{1}}(0)\subset\Omega_{0}\subset B_{r_{2}}(0)$ for some $r_{1},r_{2}>0$.
Let $r_{\infty}>0$ be defined such that $|B_{r_{\infty}}(0)|=|\Omega_{0}|$.

Let $(\rho,p,n)$ solve (P) with $\rho_{0}=\chi_{\Omega_{0}}$. Then there
exists a constant $C>0$ only depending on $r_{1}$ and $r_{2}$, but not on
$n_{0}$, such that

$W_{2}\Big{(}\rho\big{(}m(t)^{\frac{1}{d}}x,t\big{)},\chi_{B_{r_{\infty}}(0)}(x)\Big{)}\leq
Cm(t)^{-\frac{1}{d}}$

for all $t\geq 0$. Here $W_{2}$ denotes the 2-Wasserstein distance.

Lastly, we present a regularity result for the rescaled solutions when
$n_{0}>1$. It is possible to use viscosity solutions approach to study (HS)
and (HS’) as pressure-driven free boundary problems. For instance, for (HS) we
have $\rho=\chi_{\\{p>0\\}}$ and the set $\\{p>0\\}$ evolves according to

$\left\\{\begin{array}[]{lll}-\Delta p=\chi_{\Omega_{0}}&\hbox{ in
}\\{p>0\\},\\\ V_{n}=|\nabla p|&\hbox{ on
}\partial\\{p>0\\}.\end{array}\right.$

In the case of a classic Hele-Shaw problem, the interior pressure equation is
replaced by $-\Delta p=0$ on a perforated domain, with the value of $p$ being
prescribed along the fixed inner boundary. In such case, the free boundary
regularity has been studied in [CJK07]. While we expect parallel results to
hold for our problem, it seems not straightforward to verify this.

######  Theorem 2.10 (Theorem 6.9).

Fix $n_{0}>1$, and let $\rho_{0}$, $r_{1}$, $r_{2}$, and $m(t)$ be as in above
theorem. Then there is $\alpha\in(0,1)$ and $T>0$ which depends on $r_{1}$,
$r_{2}$, and $n_{0}$, such that the followings hold for all $t\geq T$.

  1. (a)

The rescaled set $\tilde{\Omega}_{t}:=m(t)^{-1/d}\Omega_{t}$ has uniformly
$C^{1,\alpha}$-boundary;

  2. (b)

The rescaled nutrient variable $\tilde{n}(x,t):=n(m^{1/d}(t)x,t)$ is uniformly
bounded in $C^{\alpha}(\\{|x|\geq 2m^{-1/d}(t)r_{2}\\})$.

##  3\. Discrete-in-time Scheme and Well-posedness

In this section, we explicitly construct solutions for the PDE (P) under the
assumption that $n_{0}\in L^{\infty}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$
and $\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ along with
$\rho_{0}\in[0,1]$ almost everywhere. We shall use the following discrete-in-
time scheme introduced in [JKT21]. Given a time step $\tau>0$, an initial
density $\rho^{0,\tau}=\rho_{0}$, and an initial nutrient density
$n^{0,\tau}=n_{0}$, we iterate that

(3.1) $\displaystyle\rho^{k+1,\tau}=$
$\displaystyle\;\operatorname*{argmin}_{\rho\leq
1}\frac{1}{2\tau}W_{2}^{2}\Big{(}\rho,\rho^{k,\tau}\big{(}1+\tau(n^{k,\tau}-b)\big{)}\Big{)},$
$\displaystyle n^{k+1,\tau}=$ $\displaystyle\;e^{\tau
D\Delta}\big{(}n^{k,\tau}(1-\tau\rho^{k+1,\tau})\big{)},$

where $e^{\tau D\Delta}$ is the heat kernel. The optimal pressure variable can
be recovered by considering the dual problem to (3.1),

(3.2) $p^{k+1,\tau}=\operatorname*{argmax}_{p\geq
0}\int_{\mathbb{R}^{d}}p^{c_{\tau}}\rho^{k,\tau}(1+\tau(n^{k,\tau}-b))-p\,dx,$

where

$p^{c_{\tau}}(x)=\inf_{y\in\mathbb{R}^{d}}p(y)+\frac{1}{2\tau}|y-x|^{2},$

is the $c$-transform. We additionally define $p^{0,\tau}=0$.

We now import the following three crucial lemmas from [JKT21] with minor
adaptations. The first lemma gives a link between the primal and dual
variables; the second lemma establishes an energy dissipation property for the
scheme; and the final lemma establishes some useful properties enjoyed by the
discrete density.

######  Lemma 3.1 ([JKT21]).

The optimal primal and dual variables $\\{\rho^{k,\tau}\\}_{k}$ and
$\\{p^{k,\tau}\\}_{k}$ are linked through the following relations

(3.3) $p^{k+1,\tau}(1-\rho^{k+1,\tau})=0,\quad\big{(}id+\tau\nabla
p^{k+1,\tau}\big{)}_{\\#}\rho^{k+1,\tau}=\rho^{k,\tau}(1+\tau(n^{k,\tau}-b)).$

######  Lemma 3.2 ([JKT21]).

Each step of the scheme enjoys the energy dissipation property

(3.4) $\frac{1}{2}\lVert\nabla
p^{k+1,\tau}\rVert_{L^{2}(\mathbb{R}^{d})}^{2}\leq\int_{\mathbb{R}^{d}}\rho^{k,\tau}(n^{k,\tau}-b)p^{k+1,\tau}\,dx.$

######  Lemma 3.3 ([JKT21]).

For almost every $x\in\mathbb{R}^{d}$, we have $\rho^{k,\tau}(1-\tau
b)\leq\rho^{k+1,\tau}$. Furthermore, if $\rho_{0}(x)\in\\{0,1\\}$ almost
everywhere and $b=0$, then $\rho^{k,\tau}(x)\in\\{0,1\\}$ almost everywhere.

Now we are ready to introduce piecewise-constant-in-time interpolants: for
$t\in[k\tau,(k+1)\tau)$, define

$\rho^{\tau}(x,t):=\rho^{k,\tau}(x),\quad n^{\tau}(x,t):=n^{k,\tau}(x),\quad
p^{\tau}(x,t):=p^{k,\tau}(x).$

In addition, for $t\leq 0$, define $\rho^{\tau}(x,t)=\rho_{0}(x)$ and
$n^{\tau}(x,t)=n_{0}(x)$. Unlike in [JKT21], our growth rate here is
independent of the pressure. Hence, we will need the following estimates to
obtain compactness for the interpolants.

######  Lemma 3.4.

Let $\rho^{\tau},p^{\tau},n^{\tau}$ be the discrete interpolants defined
above. For any $t\geq 0$, we have

$\lVert\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq e^{t\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}\lVert\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}=:B(t),$

and

$\displaystyle\lVert p^{\tau}\rVert_{L^{2}(\mathbb{R}^{d}\times[0,t])}\leq$
$\displaystyle\;CB(t)^{\frac{d+4}{2d}}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}^{\frac{1}{2}},$ $\displaystyle\lVert
p^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}\leq$
$\displaystyle\;CB(t)^{\frac{d+2}{d}},$ $\displaystyle\lVert\nabla
p^{\tau}\rVert_{L^{2}(\mathbb{R}^{d}\times[0,t])}\leq$
$\displaystyle\;CB(t)^{\frac{d+2}{2d}}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}^{\frac{1}{2}}.$

Here $C>0$ is a universal constant only depending on $d$. Moreover, when
$\tau\ll 1$ such that $\tau b<1$,

$\lVert\nabla\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert\nabla
n^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq e^{(2\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}+1)t}\Big{(}\lVert\nabla\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert\nabla
n_{0}\rVert_{L^{1}(\mathbb{R}^{d})}\Big{)}.$

###### Proof.

It is clear that $\|n^{\tau}(\cdot,t)\|_{L^{\infty}}$ is non-increasing with
respect to time and thus

$\lVert\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\big{(}1+\tau\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\big{)}\lVert\rho^{\tau}(\cdot,t-\tau)\rVert_{L^{1}(\mathbb{R}^{d})}.$

Iterating and using the fact that $(1+\tau\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})})^{k}\leq e^{\tau k\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}$ for any $k$, the first result
follows.

For the second result, we can use the duality relation
$p^{\tau}(1-\rho^{\tau})=0$ to obtain

$\lVert
p^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\lVert\rho^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}\lVert
p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}\leq\lVert\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}^{\frac{1}{2}}\lVert
p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}.$

Next, the Gagliardo-Nirenberg inequality implies that

$\lVert p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}\leq C\lVert
p^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}^{\frac{2}{d+2}}\lVert\nabla
p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}^{\frac{d}{d+2}},$

where $C>0$ is a universal constant only depending on $d$. Thus,

$\lVert p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}\leq
C\lVert\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}^{\frac{1}{d}}\lVert\nabla
p^{\tau}(\cdot,t)\rVert_{L^{2}(\mathbb{R}^{d})}.$

Integrating in time and combining this with the energy dissipation inequality
(3.4), we see that

$\lVert p^{\tau}\rVert_{L^{2}(\mathbb{R}^{d}\times[0,t])}^{2}\leq
CB(t)^{\frac{2}{d}}\lVert\nabla
p^{\tau}\rVert_{L^{2}(\mathbb{R}^{d}\times[0,t])}^{2}\leq
CB(t)^{\frac{2}{d}}\lVert
p^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}.$

Here we used the fact that

$\lVert\nabla p^{k+1,\tau}\rVert_{L^{2}(\mathbb{R}^{d})}^{2}\leq
2\int_{\mathbb{R}^{d}}\rho^{k,\tau}n^{k,\tau}p^{k+1,\tau}\,dx\leq 2\lVert
p^{k+1,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}.$

Reusing the estimate on the pressure $L^{1}$-norm from above and noticing that

$\int_{0}^{t}B(s)\,ds\leq\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}^{-1}B(t),$

we get

$\lVert p^{\tau}\rVert_{L^{2}(\mathbb{R}^{d}\times[0,t])}\leq
CB(t)^{\frac{d+4}{2d}}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}^{\frac{1}{2}}.$

Combining this with our previous work, we get all of the estimates except for
the last one.

To prove the final estimate, we can use the $BV$-bound from [DPMSV16] to
obtain that, when $\tau b<1$,

$\begin{split}\lVert\nabla\rho^{k+1,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}\leq&\;\big{\|}\nabla\big{(}\rho^{k,\tau}(1+\tau(n^{k,\tau}-b))\big{)}\big{\|}_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\lVert\nabla\rho^{k,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}+\tau\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\nabla\rho^{k,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}+\tau\lVert\nabla
n^{k,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}.\end{split}$

It is then straightforward to see that the interpolants satisfy

$\lVert\nabla\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\lVert\nabla\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\nabla\rho^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}+\lVert\nabla
n^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}.$

From the discrete scheme, we also have

$\begin{split}\lVert\nabla
n^{k+1,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}\leq&\;\big{\|}\nabla e^{\tau
D\Delta}\big{(}n^{k,\tau}(1-\tau\rho^{k,\tau})\big{)}\big{\|}_{L^{1}(\mathbb{R}^{d})}\leq\big{\|}\nabla\big{(}n^{k,\tau}(1-\tau\rho^{k,\tau})\big{)}\big{\|}_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\lVert\nabla n^{k,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}+\tau\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\nabla\rho^{k,\tau}\rVert_{L^{1}(\mathbb{R}^{d})}.\end{split}$

Thus, the interpolants satisfy

$\lVert\nabla n^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\lVert\nabla
n_{0}\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\nabla\rho^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}.$

Summing the two estimates together, we see that

$\begin{split}&\;\lVert\nabla\rho^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert\nabla
n^{\tau}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\lVert\nabla\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert\nabla
n_{0}\rVert_{L^{1}(\mathbb{R}^{d})}\\\ &\;+\big{(}2\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}+1\big{)}\Big{(}\lVert\nabla\rho^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}+\lVert\nabla
n^{\tau}\rVert_{L^{1}(\mathbb{R}^{d}\times[0,t])}\Big{)}.\end{split}$

The final claimed inequality now follows by Gronwall’s inequality. ∎

In addition to the above pressure estimates, we have the following control on
the time derivatives of the density.

######  Lemma 3.5.

Let $\rho^{\tau}$ be the discrete density interpolant defined above, with
$\tau b<1$. For any $0\leq t_{0}<t_{1}\leq T$, we have

$\lVert\rho^{\tau}(\cdot,t_{1})-\rho^{\tau}(\cdot,t_{0})\rVert_{L^{1}(\mathbb{R}^{d})}\leq(t_{1}-t_{0}+\tau)\Big{(}2b+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\Big{)}\lVert\rho^{\tau}\rVert_{L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))}.$

###### Proof.

Using the fact that for any $k\geq 0$ we have $\rho^{k+1,\tau}\geq(1-\tau
b)\rho^{k,\tau}$, we can estimate

$\begin{split}&\;\lVert\rho^{\tau}(\cdot,t_{1})-\rho^{\tau}(\cdot,t_{0})\rVert_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\int_{\mathbb{R}^{d}}\rho^{\tau}(x,t_{1})+\rho^{\tau}(x,t_{0})-2(1-\tau
b)^{\lceil\frac{1}{\tau}(t_{1}-t_{0})\rceil}\rho^{\tau}(x,t_{0})\,dx\\\
\leq&\;2\left[1-(1-\tau
b)^{\lceil\frac{1}{\tau}(t_{1}-t_{0})\rceil}\right]\lVert\rho^{\tau}(\cdot,t_{0})\rVert_{L^{1}(\mathbb{R}^{d})}+\tau\sum_{k=\lfloor
t_{0}/\tau\rfloor}^{\lfloor
t_{1}/\tau\rfloor-1}\int_{\mathbb{R}^{d}}\rho^{k,\tau}n^{k,\tau}\,dx\\\
\leq&\;(t_{1}-t_{0}+\tau)\Big{(}2b+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\Big{)}\lVert\rho^{\tau}\rVert_{L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))}.\end{split}$

∎

Using the estimates above, we can prove that the discrete interpolants
converge to a weak solution of the continuum PDE as we send $\tau\to 0$.

######  Proposition 3.6.

Assume $n_{0}\in L^{\infty}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ and
$\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ along with
$\rho_{0}\in[0,1]$ almost everywhere. Take an arbitrary finite $T>0$ and
denote $Q_{T}=\mathbb{R}^{d}\times[0,T]$. The family
$\\{\rho^{\tau}\\}_{\tau>0}$ is strongly $L^{1}(Q_{T})$ precompact,
$\\{p^{\tau}\\}_{\tau>0}$ is weakly precompact in
$L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, and $\\{n^{\tau}\\}_{\tau>0}$ is
precompact in the weak-$*$ topology of $L^{\infty}(Q_{T})$. Let $(\rho,p,n)$
be a limit point in the above-mentioned topology. Then $(\rho,p,n)$ is a weak
solution of the tumor growth PDE (P). More precisely, we have $\rho\in
L^{1}(Q_{T})$, $p\in L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, and $n\in
L^{\infty}(Q_{T})$, satisfying that for any $\psi\in H^{1}(Q_{T})$ such that
$\psi(\cdot,T)=0$ almost everywhere, it holds

(3.5) $\int_{0}^{T}\int_{\mathbb{R}^{d}}\nabla\psi\cdot\nabla
p-\rho\partial_{t}\psi\,dx\,dt=\int_{\mathbb{R}^{d}}\psi(x,0)\rho_{0}(x)\,dx+\int_{0}^{T}\int_{\mathbb{R}^{d}}\psi\rho(n-b)\,dx\,dt,$

and

$\partial_{t}n-D\Delta n=-\rho n\mbox{ in }\mathcal{D}^{\prime}(Q_{T}),\quad
n(x,0)=n_{0}(x).$

Moreover, $\rho\in C^{0,1}([0,T];L^{1}(\mathbb{R}^{d}))\cap
L^{\infty}([0,T];BV(\mathbb{R}^{d}))$, $n\in
L^{\infty}([0,T];BV(\mathbb{R}^{d}))$, and $\rho\in[0,1]$, $p\geq 0$ with
$p(1-\rho)=0$ a.e. on $Q_{T}$. If $\rho_{0}\in\\{0,1\\}$ and $b=0$, then for
all $t\in[0,T]$ we have $\rho\in\\{0,1\\}$ a.e. in $\mathbb{R}^{d}$.

###### Proof.

From the estimates in Lemma 3.4, $\\{p^{\tau}\\}_{\tau>0}$ is weakly
precompact in $L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$. Thanks to Lemma 3.5, for
any $0\leq t_{0}\leq t_{1}\leq T$,

(3.6) $\limsup_{\tau\to
0}\lVert\rho^{\tau}(\cdot,t_{1})-\rho^{\tau}(\cdot,t_{0})\rVert_{L^{1}(\mathbb{R}^{d})}\leq(t_{1}-t_{0})(2b+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})})B(T).$

Combining this with Lemma 3.4, we know that $\\{\rho^{\tau}\\}_{\tau>0}$ is
uniformly bounded and equicontinuous in $L^{1}(Q_{T})$ (in space-time). Thus,
by the Riesz-Fréchet-Kolmogorov compactness theorem,
$\\{\rho^{\tau}\\}_{\tau>0}$ is strongly $L^{1}(Q_{T})$ precompact.

With these compactness properties, now we turn to considering the PDE. Given a
smooth test function $\psi$ that vanishes at time $T$, we can use (3.3) to
obtain

$\begin{split}&\;\int_{\mathbb{R}^{d}}\frac{\rho^{\tau}(x,t)-\rho^{\tau}(x,t-\tau)}{\tau}\psi(x,t)\,dx\\\
=&\;\int_{\mathbb{R}^{d}}\rho^{\tau}(x,t)\frac{\psi(x,t)-\psi(x+\tau\nabla
p^{\tau}(x,t),t)}{\tau}+(n^{\tau}(x,t-\tau)-b)\rho^{\tau}(x,t-\tau)\psi(x,t)\,dx.\end{split}$

Integrating both sides in time along $[\tau,T]$, we get

$\begin{split}&\;\int_{\tau}^{T-\tau}\int_{\mathbb{R}^{d}}\rho^{\tau}(x,t)\frac{\psi(x,t)-\psi(x,t+\tau)}{\tau}\,dx\,dt\\\
&\;+\frac{1}{\tau}\int_{0}^{\tau}\int_{\mathbb{R}^{d}}\rho^{\tau}(x,t+T-\tau)\psi(x,t+T-\tau)-\rho^{\tau}(x,t)\psi(x,t+\tau)\,dx\,dt\\\
=&\;\int_{\tau}^{T}\int_{\mathbb{R}^{d}}\rho^{\tau}(x,t)\frac{\psi(x,t)-\psi(x+\tau\nabla
p^{\tau}(x,t),t)}{\tau}+(n^{\tau}(x,t-\tau)-b)\rho^{\tau}(x,t-\tau)\psi(x,t)\,dx\,dt.\end{split}$

Thanks to the smoothness of $\psi$ and the estimates from above, the previous
line is equivalent to

(3.7)
$\int_{0}^{T}\int_{\mathbb{R}^{d}}-\rho^{\tau}\partial_{t}\psi\,dx\,dt-\int_{\mathbb{R}^{d}}\rho_{0}\psi(x,0)\,dx=\epsilon_{\tau}+\int_{0}^{T}\int_{\mathbb{R}^{d}}-\rho^{\tau}\nabla\psi\cdot\nabla
p^{\tau}+(n^{\tau}-b)\rho^{\tau}\psi\,dx\,dt.$

Here

$\begin{split}|\epsilon_{\tau}|\leq&\;C\sqrt{\tau}\|\nabla\psi\|_{L^{\infty}(Q_{T})}\lVert\nabla
p^{\tau}\rVert_{L^{2}(Q_{T})}\lVert\rho^{\tau}\rVert_{L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))}^{1/2}\\\
&\;+C\tau\left[\lVert D^{2}\psi\rVert_{L^{\infty}(Q_{T})}\lVert\nabla
p^{\tau}\rVert_{L^{2}(Q_{T})}^{2}\right.\\\
&\;\qquad\left.+\big{(}\lVert\partial_{t}\psi\rVert_{L^{\infty}(Q_{T})}+\lVert\partial_{t}^{2}\psi\rVert_{L^{\infty}(Q_{T})}\big{)}\lVert\rho^{\tau}\rVert_{L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))}\right],\end{split}$

where $C$ is a universal constant depending on $\|n_{0}\|_{L^{\infty}}$, $b$,
and $T$. Clearly, $\lim_{\tau\to 0}|\epsilon_{\tau}|=0$. Since $p^{\tau}\in
L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$ and $p^{\tau}(1-\rho^{\tau})=0$, it
follows that $\rho^{\tau}\nabla p^{\tau}=\nabla p^{\tau}$.

Now we can claim that, there exists $\rho\in L^{1}(Q_{T})$, $p\in
L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, $n\in L^{\infty}(Q_{T})$, and a sequence
$\\{\tau_{j}\\}_{j}$ converging to $0$, such that
$\rho^{\tau_{j}}\to\rho\in[0,1]$ strongly in $L^{1}(Q_{T})$,
$p^{\tau_{j}}\rightharpoonup p\geq 0$ weakly in
$L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, and $n^{\tau_{j}}$ converges weak-$*$ in
$L^{\infty}(Q_{T})$ to $n$. It is then clear that we can pass to the limit in
(3.7) to obtain (3.5). One may relax the regularity of $\psi$ to find (3.5)
actually holds for all $\psi\in H^{1}(Q_{T})$ satisfying $\psi(\cdot,T)=0$
almost everywhere.

The strong convergence from $\rho^{\tau_{j}}$ to $\rho$ in $L^{1}(Q_{T})$ also
implies

$\lim_{j\to+\infty}\|\rho^{\tau_{j}}(\cdot,t)-\rho(\cdot,t)\|_{L^{1}(\mathbb{R}^{d})}=0$

for almost every $t\in[0,T]$. Hence, by (3.6), up to modifying $\rho(x,t)$ for
those $t$ on a measure-zero set in $[0,T]$ if necessary, we can make $\rho$ be
Lipschitz continuous in $L^{1}(\mathbb{R}^{d})$ with respect to time. That
$\rho,n\in L^{\infty}([0,T];BV(\mathbb{R}^{d}))$ follows from Lemma 3.4.

Thanks to the strong $L^{1}$-convergence of $\\{\rho^{\tau_{j}}\\}_{j}$, we
can readily verify that the $n$-equation and $p(1-\rho)=0$ hold in the sense
of distribution (we omit the details). Given the regularity of $p$ and $\rho$,
$p(1-\rho)=0$ almost everywhere in $Q_{T}$. Hence, we may suitably modify $p$
on a measure-zero set of $Q_{T}$, if necessary, to achieve $p(1-\rho)=0$
everywhere on $Q_{T}$.

By the strong $L^{1}$-convergence of $\\{\rho^{\tau_{j}}\\}_{j}$, Lemma 3.3,
and the fact that $\rho$ is Lipschitz in $L^{1}(\mathbb{R}^{d})$ with respect
to time, we find that for all $t\in[0,T]$, $\rho\in\\{0,1\\}$ almost
everywhere if $\rho_{0}\in\\{0,1\\}$ almost everywhere and $b=0$. ∎

Then we show that our solution satisfies the so-called complementarity
condition.

######  Proposition 3.7.

If $(\rho,p,n)$ is a weak solution of the tumor growth PDE (P) (in the sense
of Definition 2.1), then for any bounded open set
$\Omega\subset\mathbb{R}^{d}$ and $\psi\in L^{2}([0,T];H^{1}(\Omega)\cap
L^{1}(\mathbb{R}^{d}))$, such that $\psi(1-\rho)=0$ and
$\psi|_{\partial\Omega\times[0,T]}=0$, we have

$\int_{0}^{T}\int_{\Omega}\nabla\psi\cdot\nabla p-\psi(n-b)\,dx\,dt=0.$

######  Remark 3.8.

Let us note that the space of $\psi$ satisfying the above conditions is
nontrivial. For instance, given a smooth function
$\eta:\mathbb{R}^{d}\times[0,T]\to\mathbb{R}$ such that $\eta$ is compactly
supported inside $\Omega\times[0,T]$, the choice $\psi:=p\eta$ satisfies all
of the above conditions.

###### Proof.

We begin by assuming that there exists $\delta>0$ such that $\psi(x,t)=0$ for
all $(x,t)\in[T-\delta,T]\times\Omega$. Fix $\epsilon\in(0,\delta)$ and let

$\psi^{\epsilon}(x,t)=\frac{1}{\epsilon}\int_{t}^{\min(T,t+\epsilon)}\psi(x,s)\,ds.$

Since $\psi^{\epsilon}\in H^{1}(\Omega\times[0,T])$, when extended by zero to
the whole $Q_{T}$, it is a valid test function for the weak formulation (see
(3.5)). Hence, we have

$\int_{\Omega\times[0,T]}-\rho\partial_{t}\psi^{\epsilon}+\nabla\psi^{\epsilon}\cdot\nabla
p-\psi^{\epsilon}\rho(n-b)\,dx\,dt=\int_{\Omega}\psi^{\epsilon}(x,0)\rho_{0}(x)\,dx.$

Note that

$\rho\partial_{t}\psi^{\epsilon}=\rho(x,t)\frac{\psi(x,\min(T,t+\epsilon))-\psi(x,t)}{\epsilon}\leq\frac{\psi(x,\min(T,t+\epsilon))-\psi(x,t)}{\epsilon}=\partial_{t}\psi^{\epsilon},$

where we use $\psi(1-\rho)=0$ and the non-negativity of $\psi$ to justify the
inequality. Therefore,

$\int_{\Omega\times[0,T]}\nabla\psi^{\epsilon}\cdot\nabla
p-\psi^{\epsilon}\rho(n-b)\,dx\,dt\leq\int_{\Omega}\psi^{\epsilon}(x,0)\big{(}\rho_{0}(x)-1\big{)}\,dx\leq
0.$

Sending $\epsilon\to 0$ and once again using $\psi(1-\rho)=0$, we see that

$\int_{\Omega\times[0,T]}\nabla\psi\cdot\nabla p-\psi(n-b)\,dx\,dt\leq 0,$

giving us one side of the equation.

To obtain the other direction, we instead smooth backwards in time and define

$\psi_{\epsilon}(x,t)=\frac{1}{\epsilon}\int_{\max(0,t-\epsilon)}^{t}\psi(x,s)\,ds.$

Note that $\psi_{\epsilon}(x,T)=0$ since $\psi$ vanishes on $[T-\delta,T]$. An
analogous argument shows that
$\rho\partial_{t}\psi_{\epsilon}\geq\partial_{t}\psi_{\epsilon}$, and thus

$\int_{\Omega\times[0,T]}\nabla\psi_{\epsilon}\cdot\nabla
p-\psi_{\epsilon}(n-b)\,dx\,dt\geq\int_{\Omega}\psi_{\epsilon}(x,0)\big{(}\rho_{0}(x)-1\big{)}\,dx\geq
0.$

Note that $\psi_{\epsilon}(x,0)\leq 0$. Sending $\epsilon\to 0$ and combining
with our previous work, we obtain

$\int_{\Omega\times[0,T]}\nabla\psi\cdot\nabla p-\psi(n-b)\,dx\,dt=0,$

for all non-negative $\psi$ such that $\psi(1-\rho)=0,$
$\psi|_{\partial\Omega\times[0,T]}=0$, and $\psi(x,t)=0$ on
$\Omega\times[T-\delta,T]$. Since the equation is linear in $\psi$ and does
not include time derivatives, we can drop the non-negativity assumption and
then take limits to drop the assumption that $\psi$ vanishes on
$\Omega\times[T-\delta,T]$. ∎

We now proceed to show the uniqueness of weak solutions. Similar results have
been established in the literature using the Hilbert duality principle [PQV14,
GKM22]. However, they require the nutrient variable to be at least
$L^{1}_{t}H^{1}_{x}$, which does not hold in our case when $D=0$. Instead, we
proceed with $L^{1}$-contraction approach to provide a unified proof of
uniqueness for all $D\geq 0$.

Let us focus on the $\rho$-equation for the moment. Consider the model problem

(3.8) $\rho_{t}-\nabla\cdot(\rho\nabla p)=f,\quad\rho\leq 1,\quad p\in
P_{\infty}(\rho),\quad\rho(x,0)=\rho_{0}.$

Similar to Definition 2.1, we introduce the notion of its weak solutions.

######  Definition 3.9.

Let $\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$ such that
$\rho_{0}\in[0,1]$ almost everywhere. Fix $T>0$ and denote
$Q_{T}=\mathbb{R}^{d}\times[0,T]$. Assume $f\in
L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))\cap L^{\infty}(Q_{T})$. Non-negative
functions $\rho\in C([0,T];L^{1}(\mathbb{R}^{d}))$ and $p\in
L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$, which are defined on
$\mathbb{R}^{d}\times[0,T]$, are said to form a weak solution of (3.8), if
they satisfy:

  1. (i)

$\rho\in[0,1]$ and $p(1-\rho)=0$ in $Q_{T}$;

  2. (ii)

For any $\psi\in C_{0}^{\infty}(Q_{T})$ that vanishes at $t=T$, we have

$\int_{0}^{T}\int_{\mathbb{R}^{d}}\nabla\psi\cdot\nabla
p-\rho\partial_{t}\psi\,dx\,dt=\int_{\mathbb{R}^{d}}\psi(x,0)\rho_{0}\,dx+\int_{0}^{T}\int_{\mathbb{R}^{d}}\psi
f\,dx\,dt.$

######  Remark 3.10.

Under the assumptions that $\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap
BV(\mathbb{R}^{d})$, and $f\in L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))\cap
L^{\infty}(Q_{T})\cap L^{\infty}([0,T];BV(\mathbb{R}^{d}))$, it is not
difficult to show existence of weak solutions of (3.8). In fact, one may use
still the discrete scheme (c.f. (3.1) and (3.2))

$\displaystyle\rho^{k+1,\tau}=$
$\displaystyle\;\operatorname*{argmin}_{\rho\leq
1}\frac{1}{2\tau}W_{2}^{2}\Big{(}\rho,\rho^{k,\tau}+\tau f^{k,\tau}\Big{)},$
$\displaystyle p^{k+1,\tau}=$ $\displaystyle\;\operatorname*{argmax}_{p\geq
0}\int_{\mathbb{R}^{d}}p^{c_{\tau}}\big{(}\rho^{k,\tau}+\tau
f^{k,\tau}\big{)}-p\,dx,$

where $f^{k,\tau}$ can be defined as, e.g., time-average of $f$ on each small
time interval of size $\tau$. Then arguing as before, we can prove the
existence. We omit the details.

We can show that the weak solution of (3.8) satisfies the comparison
principle. See Appendix A for the proof, which follows the Hilbert duality
argument in [PQV14] with minor modifications.

######  Lemma 3.11.

For $i=0,1$, assume $\rho_{0}^{i}\in L^{1}(\mathbb{R}^{d})\cap
BV(\mathbb{R}^{d})$ such that $\rho_{0}^{i}\in[0,1]$ almost everywhere. Take a
finite $T>0$ and denote $Q_{T}=\mathbb{R}^{d}\times[0,T]$. Let $f^{i}\in
L^{\infty}([0,T];L^{1}(\mathbb{R}^{d}))\cap L^{\infty}(Q_{T})$. Let
$(\rho^{i},p^{i})$ $(i=0,1)$ be weak solutions of

$\rho^{i}_{t}-\nabla\cdot(\rho^{i}\nabla p^{i})=f^{i},\quad\rho^{i}\leq
1,\quad p^{i}\in P_{\infty}(\rho^{i}),\quad\rho^{i}(x,0)=\rho^{i}_{0},$

respectively. If $\rho_{0}^{0}\leq\rho_{0}^{1}$ almost everywhere in
$\mathbb{R}^{d}$ and $f^{0}\leq f^{1}$ almost everywhere in $Q_{T}$, then
$\rho^{0}\leq\rho^{1}$ almost everywhere in $Q_{T}$.

######  Remark 3.12.

In a similar spirit, one can show that if $(\rho,p,n)$ is a weak solution of
(P) on $Q_{T}$ with $\rho_{0}$ being compactly supported, then $\rho$ and $p$
are compactly supported in $Q_{T}$. Indeed, assuming
$\rho_{0}\leq\chi_{B_{r_{0}}}$ for some $r_{0}>0$, we can prove
$\rho\leq\tilde{\rho}$ and $p\leq\tilde{p}$, where

$\tilde{\rho}(x,t):=\chi_{B_{r(t)}}(x),\;\tilde{p}(x,t):=\frac{\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}}{2d}\left(r(t)^{2}-|x|^{2}\right)_{+},\mbox{
and
}r(t):=r_{0}\exp\left(\frac{1}{d}\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}t\right).$

Note that $(\tilde{\rho},\tilde{p})$ is chosen to be a weak solution (in the
sense of Definition 3.9) of

$\partial_{t}\tilde{\rho}-\nabla(\tilde{\rho}\nabla\tilde{p})=\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}\tilde{\rho},\quad\tilde{\rho}\leq
1,\quad\tilde{p}\in
P_{\infty}(\tilde{\rho}),\quad\tilde{\rho}(x,0)=\chi_{B_{r_{0}}}(x).$

We skip the details.

Now we present uniqueness of the weak solution for (P).

######  Proposition 3.13.

Under the assumption of Proposition 3.6, the weak solution of (P) is unique.

###### Proof.

Fix $T>0$. Suppose $(\rho,p,n)$ and $(\tilde{\rho},\tilde{p},\tilde{n})$ are
two weak solutions. We first show an $L^{1}$-contraction principle regarding
$\rho$ and $\tilde{\rho}$.

Denote

$f=(n-b)\rho,\quad\tilde{f}=(\tilde{n}-b)\tilde{\rho},\quad
f_{*}=\max\big{\\{}f,\tilde{f}\big{\\}}.$

Since $f,\tilde{f}\in L^{\infty}([0,T];BV(\mathbb{R}^{d}))$ (c.f. Proposition
3.6), we have $f_{*}\in L^{\infty}([0,T];BV(\mathbb{R}^{d}))$. Let
$(\rho_{*},p_{*})$ be a weak solution (see Remark 3.10) of

$(\rho_{*})_{t}-\nabla\cdot(\rho_{*}\nabla p_{*})=f_{*},\quad\rho_{*}\leq
1,\quad p_{*}\in P_{\infty}(\rho_{*}),\quad\rho_{*}(x,0)=\rho_{0}.$

Then by Lemma 3.11, $\rho,\tilde{\rho}\leq\rho_{*}$ almost everywhere in
$Q_{T}$.

Take an arbitrary $t_{*}\in(0,T]$. Let $\zeta:\mathbb{R}^{d}\to\mathbb{R}$ be
a smooth non-negative function such that $\zeta=1$ on the unit ball, and
$\zeta=0$ outside the ball of radius 2. For any $R>0$, set
$\zeta_{R}(x):=\zeta(x/R)$. With $\delta<t_{*}$, let $\eta_{\delta}(t)\in
C^{\infty}([0,T])$, such that it is non-increasing, $\eta_{\delta}\equiv 1$ on
$[0,t_{*}-\delta]$, and $\eta_{\delta}\equiv 0$ on $[t_{*},T]$. Then we take
$\psi(x,t)=\zeta_{R}(x)\eta_{\delta}(t)\in H^{1}(Q_{T})$ in Definition 2.1 and
Definition 3.9, and take a difference of them to obtain that

$\int_{0}^{T}-\partial_{t}\eta_{\delta}(t)\big{\|}\big{(}\rho(\cdot,t)-\rho_{*}(\cdot,t)\big{)}\zeta_{R}\big{\|}_{L^{1}(\mathbb{R}^{d})}\,dt=\int_{0}^{T}\eta_{\delta}(t)\big{\|}\big{(}f(\cdot,t)-f_{*}(\cdot,t)\big{)}\zeta_{R}\big{\|}_{L^{1}(\mathbb{R}^{d})}\,dt+\epsilon_{R},$

where $\epsilon_{R}$ is an error term with $|\epsilon_{R}|\lesssim
R^{-2}(\lVert p\rVert_{L^{1}(Q_{T})}+\lVert p_{*}\rVert_{L^{1}(Q_{T})})$ and
we used the facts that $\rho\leq\rho_{*}$ and $f\leq f_{*}$. Here one can
derive an estimate for $\lVert p_{*}\rVert_{L^{1}(Q_{T})}$ as in Lemma 3.4.
Sending $R\to\infty$ and then $\delta\to 0$, we can use the time-continuity of
$\rho$ and $\rho_{*}$ (see Definition 2.1 and Definition 3.9), to obtain

$\|\rho_{*}(\cdot,t_{*})-\rho(\cdot,t_{*})\|_{L^{1}(\mathbb{R}^{d})}=\|f_{*}-f\|_{L^{1}(\mathbb{R}^{d}\times[0,t_{*}])}.$

Similarly,

$\|\rho_{*}(\cdot,t_{*})-\tilde{\rho}(\cdot,t_{*})\|_{L^{1}(\mathbb{R}^{d})}=\|f_{*}-\tilde{f}\|_{L^{1}(\mathbb{R}^{d}\times[0,t_{*}])}.$

As a result,

(3.9)
$\begin{split}\|\tilde{\rho}(\cdot,t_{*})-\rho(\cdot,t_{*})\|_{L^{1}(\mathbb{R}^{d})}\leq&\;\|f_{*}-f\|_{L^{1}(\mathbb{R}^{d}\times[0,t_{*}])}+\|f_{*}-\tilde{f}\|_{L^{1}(\mathbb{R}^{d}\times[0,t_{*}])}\\\
=&\;\|f-\tilde{f}\|_{L^{1}(\mathbb{R}^{d}\times[0,t_{*}])},\end{split}$

where in the last equality, we used the definition of $f_{*}$.

Now by the definition of $f$ and $\tilde{f}$, for any $t\in[0,T]$,

$\begin{split}&\;\|\tilde{\rho}(\cdot,t)-\rho(\cdot,t)\|_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\int_{0}^{t}(\|n_{0}\|_{L^{\infty}}+b)\|\tilde{\rho}(\cdot,\tau)-\rho(\cdot,\tau)\|_{L^{1}(\mathbb{R}^{d})}+\|\tilde{n}(\cdot,\tau)-n(\cdot,\tau)\|_{L^{1}(\mathbb{R}^{d})}\,d\tau.\end{split}$

Using the Duhamel’s formula for the heat equation, for every $t\in[0,T]$ and
almost every $x\in\mathbb{R}^{d}$, we have

$\tilde{n}(\cdot,t)-n(\cdot,t)=-\int_{0}^{t}e^{(t-\tau)D\Delta}\big{(}\tilde{n}(\cdot,\tau)\tilde{\rho}(\cdot,\tau)-n(\cdot,\tau)\rho(\cdot,\tau)\big{)}\,d\tau.$

Note that this formula is still valid even when $D=0$. Since the heat kernel
is a contraction on $L^{1}$, it follows that

$\begin{split}\lVert\tilde{n}(\cdot,t)-n(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq&\;\int_{0}^{t}\lVert\tilde{n}(\cdot,\tau)\tilde{\rho}(\cdot,\tau)-n(\cdot,\tau)\rho(\cdot,\tau)\rVert_{L^{1}(\mathbb{R}^{d})}\,d\tau\\\
\leq&\;\int_{0}^{t}\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\tilde{\rho}(\cdot,\tau)-\rho(\cdot,\tau)\rVert_{L^{1}(\mathbb{R}^{d})}+\lVert\tilde{n}(\cdot,\tau)-n(\cdot,\tau)\rVert_{L^{1}(\mathbb{R}^{d})}\,d\tau.\end{split}$

Thus,

$\begin{split}&\;\|\tilde{\rho}(\cdot,t)-\rho(\cdot,t)\|_{L^{1}(\mathbb{R}^{d})}+\lVert\tilde{n}(\cdot,t)-n(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\\\
\leq&\;\int_{0}^{t}(2\|n_{0}\|_{L^{\infty}}+b)\|\tilde{\rho}(\cdot,\tau)-\rho(\cdot,\tau)\|_{L^{1}(\mathbb{R}^{d})}+2\|\tilde{n}(\cdot,\tau)-n(\cdot,\tau)\|_{L^{1}(\mathbb{R}^{d})}\,d\tau.\end{split}$

It now follows from Gronwall’s inequality that
$\|\tilde{\rho}(\cdot,t)-\rho(\cdot,t)\|_{L^{1}(\mathbb{R}^{d})}+\lVert\tilde{n}(\cdot,t)-n(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}=0$
for all $t\in[0,T]$, and thus $\rho(\cdot,t)=\tilde{\rho}(\cdot,t)$ and
$n(\cdot,t)=\tilde{n}(\cdot,t)$ for all $t\in[0,T]$. Lastly, thanks to the
weak formulation (3.5), for any $\psi\in L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$

$\int_{0}^{T}\int_{\mathbb{R}^{d}}\nabla\psi\cdot\nabla\big{(}p-\tilde{p}\big{)}\,dx\,dt=0.$

Choosing $\psi$ that approximates $(p-\tilde{p})$, we conclude that
$p=\tilde{p}$. ∎

##  4\. The Case of Zero Nutrient Diffusion and No Death

In the rest of the paper, we will always assume $b=D=0$ in (P). This section
aims at studying qualitative properties of the solution. Recall that in this
case, $\eta:=n_{0}-n$ solves

(4.1) $\partial_{t}\eta=(n_{0}-\eta)\rho,\quad\eta(x,0)=0.$

Given $(\rho,p,n)$ as a weak solution of (P), it is straightforward to show
that $\eta$ is continuous in $L^{1}(\mathbb{R}^{d})$ with respect to time.

The assumption $b=0$ implies the tumor region is always expanding. This
enables us to study the $\rho$-equation of (P) through the so-called Baiocchi
transform [BCMP73]. Indeed, if we define $w(x,t)=\int_{0}^{t}p(x,s)\,ds$, we
see that for any $t\geq 0$,

$w(x,t)(1-\rho(x,t))=\int_{0}^{t}p(x,s)(1-\rho(x,t))\,ds\leq\int_{0}^{t}p(x,s)(1-\rho(x,s))\,ds=0$

everywhere on $\mathbb{R}^{d}$. Hence, by integrating the $\rho$-equation of
(P) over the time interval $[0,t]$ and using the time continuity of $\rho$ and
$\eta$ in $L^{1}(\mathbb{R}^{d})$, we get the elliptic equation

(4.2) $\rho-\Delta w=\rho_{0}+\eta,\quad w(1-\rho)=0,\quad\rho\in[0,1],\quad
w\geq 0$

along with the boundary condition that $w$ vanishes at infinity.

The main advantage of this formulation is that $w$ enjoys much better
regularity as compared to $p$ (clearly $w\in W^{2,r}(\mathbb{R}^{d})$ for any
$r<\infty$) while still describing the active tumor region. Section 4.1 will
be focused on proving comparison/contraction properties of this elliptic
equation.

On the other hand, that $D=0$ gives rise to many nice properties of the model.
For example, it gives a direct pointwise dependence of the nutrient variables
$n$ and $\eta$ on the density $\rho$

(4.3)
$\eta(x,t)=n_{0}(x)\left[1-\exp\left(-\int_{0}^{t}\rho(x,s)\,ds\right)\right].$

This allows us to extend the comparison/contraction results for
$\rho$-equation (4.2) to the full PDE (P). See Section 4.2.

###  4.1. The elliptic formulation

Let us start from the elliptic equation (4.2). We first establish the
following comparison principles and stability estimates.

######  Proposition 4.1.

Let $\Omega\subset\mathbb{R}^{d}$ be an open set with a smooth (possibly
empty) boundary. Suppose that for $i\in\\{0,1\\}$ we have non-negative
functions $\rho^{i},f^{i}\in C([0,T];L^{1}(\Omega))\cap
L^{\infty}(\Omega\times[0,T])$ and $w^{i}\in L^{\infty}([0,T];H^{2}(\Omega))$
such that for all $t\in[0,T]$

$\rho^{i}-\Delta w^{i}=f^{i},\quad
w^{i}(1-\rho^{i})=0,\quad\rho^{i}\in[0,1],\quad w^{i}\geq 0$

on $\Omega$. If $w^{0}(x,t)\leq w^{1}(x,t)$ for all
$(x,t)\in\partial\Omega\times[0,T]$, and $f^{0}(x,t)\leq f^{1}(x,t)$ for every
$t\in[0,T]$ and almost every $x\in\Omega$, then for all $t\in[0,T]$ and almost
every $x\in\Omega$ we have $\rho^{0}(x,t)\leq\rho^{1}(x,t)$ and
$w^{0}(x,t)\leq w^{1}(x,t)$.

######  Remark 4.2.

Let us emphasize that this comparison property does not immediately imply that
solutions to (P) satisfy a comparison principle, since we have not shown that
the nutrient variables stay ordered along solutions to the PDE (however, see
Theorem 4.7).

###### Proof.

Taking the difference of the two equations and integrating against
$(w^{0}-w^{1})_{+}$ on $\Omega\times[0,T]$, we obtain

$\int_{0}^{T}\int_{\Omega}(\rho^{0}-\rho^{1})(w^{0}-w^{1})_{+}+|\nabla(w^{0}-w^{1})_{+}|^{2}\,dx\,dt=\int_{0}^{T}\int_{\Omega}(f^{0}-f^{1})(w^{0}-w^{1})_{+}\,dx\,dt,$

where we have used the fact that $(w_{0}-w_{1})_{+}$ vanishes on
$\partial\Omega\times[0,T]$. The right hand side of the equation is clearly
non-positive. On the other hand, if $w^{0}>w^{1}$ then we must have $w^{0}>0$.
Therefore
$(\rho^{0}-\rho^{1})(w^{0}-w^{1})_{+}=(1-\rho^{1})(w^{0}-w^{1})_{+}\geq 0.$
Thus, the left hand side of the equation is non-negative. Hence, both sides
must be equal to zero, which allows us to conclude that
$(w^{0}-w^{1})_{+}=c(t)$ in $\Omega\times[0,T]$ for some non-negative constant
$c$ that depends only on time. Since both $w^{0}$ and $w^{1}$ approach zero at
infinity and $w^{0}\leq w^{1}$ on $\partial\Omega\times[0,T]$, it follows that
$c=0$. Thus, $w^{0}\leq w^{1}$ almost everywhere on $\Omega\times[0,T]$.

It is clear from the equation and the fact $w^{0}\leq w^{1}$ a.e. that
$\\{(x,t)\in\Omega\times[0,T]:\rho^{0}(x,t)>\rho^{1}(x,t)\\}\subset\\{(x,t)\in\Omega\times[0,T]:w^{1}(x,t)=0,w^{0}(x,t)=0\\}$
up to a set of measure zero. Since the equation implies that $w^{i}\in
L^{\infty}([0,T];W^{2,r}(\mathbb{R}^{d}))$ for any $r<\infty$, it follows that
$\Delta w^{i}=0$ a.e. on the set where $w^{i}=0$. Thus, $\rho^{i}=f^{i}$
almost everywhere on
$\\{(x,t)\in\Omega\times[0,T]:w^{1}(x,t)=0,w^{0}(x,t)=0\\}$. This allows us to
conclude that $\rho^{0}\leq\rho^{1}$ almost everywhere on $\Omega\times[0,T]$.
∎

A similar argument gives us the following $L^{1}$-stability property. It is a
strengthening of Lemma 3.11 and the $L^{1}$-contraction principle (3.9) in the
proof of Proposition 3.13.

######  Proposition 4.3.

Let $\Omega\subset\mathbb{R}^{d}$ be an open set with a smooth (possibly
empty) boundary. Suppose that for $i\in\\{0,1\\}$, $\rho_{0}^{i}(x)\in[0,1]$
are integrable on $\Omega$, $\rho^{i},\eta^{i}\in C([0,T];L^{1}(\Omega))\cap
L^{\infty}(\Omega\times[0,T])$ and $w^{i}\in L^{\infty}([0,T];H^{2}(\Omega))$,
such that

$\rho^{i}-\Delta w^{i}=\rho_{0}^{i}+\eta^{i},\quad
w^{i}(1-\rho^{i})=0,\quad\rho^{i}\in[0,1],\quad w^{i}\geq 0$

on $\Omega$. If $w^{0}(x,t)=w^{1}(x,t)$ for all
$(x,t)\in\partial\Omega\times[0,T]$, then for all $t\in[0,T]$

$\big{\|}(\rho^{1}(\cdot,t)-\rho^{0}(\cdot,t))_{+}\big{\|}_{L^{1}(\Omega)}\leq\left\|\Big{(}\rho^{1}_{0}(\cdot)-\rho^{0}_{0}(\cdot)+\eta^{1}(\cdot,t)-\eta^{0}(\cdot,t)\Big{)}_{+}\right\|_{L^{1}(\Omega)}.$

###### Proof.

Given small $\delta>0$, let $f_{\delta}:\mathbb{R}\to\mathbb{R}$ be a smooth
increasing function such that $f_{\delta}(a)=0$ for all $a\leq 0$ and
$f_{\delta}(a)=1$ for all $a\geq\delta$. Choose some $0\leq t_{1}<t_{2}\leq
T$. Taking the difference of the two equations and integrating against
$f_{\delta}(w^{1}-w^{0})$ on $\Omega\times[t_{1},t_{2}]$, we have

$\begin{split}&\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}-\rho^{0})f_{\delta}(w^{1}-w^{0})+f_{\delta}^{\prime}(w^{1}-w^{0})|\nabla(w^{1}-w^{0})|^{2}\,dx\,dt\\\
&=\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}_{0}-\rho^{0}_{0}+\eta^{1}-\eta^{0})f_{\delta}(w^{1}-w^{0})\,dx\,dt.\end{split}$

Sending $\delta\to 0$, we see that

$\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}-\rho^{0})\operatorname{sgn}_{+}(w^{1}-w^{0})\,dx\,dt\leq\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}_{0}-\rho^{0}_{0}+\eta^{1}-\eta^{0})\operatorname{sgn}_{+}(w^{1}-w^{0})\,dx\,dt.$

Let $E=\\{(x,t)\in\Omega\times[t_{1},t_{2}]:w^{1}(x,t)\leq
w^{0}(x,t),\,\rho^{0}(x,t)<\rho^{1}(x,t)\\}$. Then $w^{1},w^{0}=0$ almost
everywhere on $E$. Thus the regularity of $w^{i}$ implies that
$\rho^{i}=\rho^{i}_{0}+\eta^{i}$ almost everywhere on $E$. Now we use the
above inequality to derive that

$\begin{split}&\;\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}-\rho^{0})_{+}\,dx\,dt\\\
=&\;\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}-\rho^{0})\operatorname{sgn}_{+}(w^{1}-w^{0})\,dx\,dt+\int_{E}(\rho^{1}-\rho^{0})_{+}\,dx\,dt\\\
\leq&\;\int_{\Omega\times[t_{1},t_{2}]}(\rho^{1}_{0}-\rho^{0}_{0}+\eta^{1}-\eta^{0})\operatorname{sgn}_{+}(w^{1}-w^{0})\,dx\,dt+\int_{E}(\rho^{1}_{0}-\rho^{0}_{0}+\eta^{1}-\eta^{0})_{+}\,dx\,dt\\\
\leq&\;\big{\|}(\rho^{1}_{0}-\rho^{0}_{0}+\eta^{1}-\eta^{0})_{+}\big{\|}_{L^{1}(\Omega\times[t_{1},t_{2}])}.\end{split}$

Then the result follows since $t_{1}$ and $t_{2}$ are arbitrary and $\rho^{i}$
and $\eta^{i}$ are continuous in $L^{1}(\Omega)$ with respect to time. ∎

Thanks to the symmetries of the Laplacian, the above stability result implies
the following derivative estimates.

######  Proposition 4.4.

Suppose $(\rho,\rho_{0},w,\eta)$ solve

$\rho-\Delta w=\rho_{0}+\eta,\quad w(1-\rho)=0,\quad\rho\in[0,1],\quad w\geq
0$

on $\mathbb{R}^{d}\times[0,T]$, where $\rho_{0}\in[0,1]$ is integrable on
$\mathbb{R}^{d}$, $\eta\geq 0$, $\rho,\eta\in
C([0,T];L^{1}(\mathbb{R}^{d}))\cap L^{\infty}(\mathbb{R}^{d}\times[0,T])$ and
$w\in L^{\infty}([0,T];H^{2}(\mathbb{R}^{d}))$. Then for any antisymmetric
matrix $A\in\mathbb{R}^{d\times d}$, vector $v\in\mathbb{R}^{d}$, and scalar
$\lambda\geq 0$, we have for all $t\in[0,T]$,

$\lVert\nabla\rho(x,t)\cdot(Ax+\lambda
x+v)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\lVert\nabla(\rho_{0}(x,t)+\eta(x,t))\cdot(Ax+\lambda
x+v)\rVert_{L^{1}(\mathbb{R}^{d})}.$

In particular, for all $t\in[0,T]$,

$\lVert\rho(x,t)\rVert_{BV(\mathbb{R}^{d})}\leq\lVert\rho_{0}(x,t)+\eta(x,t)\rVert_{BV(\mathbb{R}^{d})}.$

###### Proof.

Given an antisymmetric matrix $A$ and $\lambda\geq 0$ we define
$A_{s}:=\exp((\lambda I+A)s)$ for each $s\in[0,1]$. A direct computation
reveals that

$\partial_{s}(A_{s}A_{s}^{\intercal})=A_{s}(\lambda
I+A)A_{s}^{\intercal}+A_{s}(\lambda I+A)^{\intercal}A_{s}^{\intercal}=2\lambda
A_{s}A_{s}^{\intercal}.$

Since $A_{0}A_{0}^{\intercal}=I$, we can integrate directly to see that
$A_{s}A_{s}^{\intercal}=e^{2\lambda s}I$. Given $v\in\mathbb{R}^{d}$, if we
define $T_{s}:\mathbb{R}^{d}\to\mathbb{R}^{d}$ such that $T_{s}(x)=A_{s}x+sv$,
then $T_{s}$ is a smooth curve of affine transformations such that
$\partial_{s}T_{0}(x)=\lambda x+Ax+v$ and
$DT_{s}DT_{s}^{\intercal}=e^{2\lambda s}I$ for all $s\in[0,1]$.

Define $\rho_{s}:=\rho\circ T_{s}$, $w_{s}:=e^{-2\lambda s}w\circ T_{s}$,
$\eta_{s}:=\eta\circ T_{s}$, and $\rho_{0,s}:=\rho_{0}\circ T_{s}$. We can
then compute

$\Delta(w_{s})=DTDT^{\intercal}:(D^{2}w\circ T_{s})+\Delta T_{s}\cdot\nabla
w_{s}=(\Delta w)\circ T_{s},$

where we have used the fact that $DTDT^{\intercal}:(D^{2}w\circ
T_{s})=e^{2\lambda s}tr(e^{-2\lambda s}D^{2}w\circ T_{s})$ and $\Delta
T_{s}=0$. Thus,

$\rho_{s}-\Delta w_{s}=\big{(}\rho-\Delta w\big{)}\circ
T_{s}=(\rho_{0}+\eta)\circ T_{s}=\rho_{0,s}+\eta_{s}.$

It is also clear that $w_{s}(1-\rho_{s})=0$. Proposition 4.3 now implies that
for all $t\in[0,T]$ and every $s\in[0,1]$,

$\lVert\rho(\cdot,t)-\rho_{s}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}\leq\lVert\rho_{0}(\cdot,t)-\rho_{0,s}(\cdot,t)+\eta(\cdot,t)-\eta_{s}(\cdot,t)\rVert_{L^{1}(\mathbb{R}^{d})}.$

Dividing both sides by $s$ and then sending $s\to 0$, the first result
follows. The second result follows from choosing $A=0$ and $\lambda=0$ and
then taking supremum over the vectors $v$ with $|v|\leq 1$. ∎

###  4.2. The full system with bounded initial nutrient

Now we study the full system (P) under the assumption that the initial
nutrient $n_{0}$ is bounded. In what follows, we shall always assume the
initial density $\rho_{0}$ is a _patch_ , namely,

$\rho_{0}\in L^{1}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d}),\mbox{ and
}\rho_{0}(x)\in\\{0,1\\}\mbox{ almost everywhere.}$

While we already know that this must imply that $\rho(x,t)\in\\{0,1\\}$ almost
everywhere for all $t$, the following lemma shows that in this case $\eta$
will be concentrated on the support of $\rho$.

######  Lemma 4.5.

Suppose $\rho_{0}(x)\in\\{0,1\\}$ almost everywhere in $\mathbb{R}^{d}$. Then
for all $t\in[0,T]$, $\rho\in\\{0,1\\}$ and $(1-\rho)\eta=0$ almost everywhere
in $\mathbb{R}^{d}$.

###### Proof.

Since $\eta\geq 0$ and $\partial_{t}\rho\geq 0$ (c.f. Proposition 4.3), a
direct computation shows that

$\frac{d}{dt}\int_{\mathbb{R}^{d}}(1-\rho(x,t))\eta(x,t)\,dx\leq\int_{\mathbb{R}^{d}}(1-\rho(x,t))\rho(x,t)(n_{0}(x)-\eta(x,t))\,dx=0.$

Thus, Gronwall’s inequality and the non-negativity of $(1-\rho)\eta$ implies
that, for all $t\in[0,T]$, $(1-\rho)\eta=0$ almost everywhere in
$\mathbb{R}^{d}$. ∎

The next result gives a more complete characterization of the tumor patch,
i.e., the set $\\{\rho=1\\}$. Formally it says that the tumor patch coincides
with the support of the pressure variable, characterizing the evolution of the
tumor as a free boundary problem driven by the pressure. It will be useful
later in Section 5. When the tumor boundary is regular, this is easy to prove
since the pressure solves the elliptic equation $-\Delta p=n$ in the interior
of the tumor. However such regularity is unknown a priori, so we instead argue
with the more regular variable $w$, the time integral of $p$.

######  Lemma 4.6.

Assume $\rho_{0}(x)\in\\{0,1\\}$ almost everywhere. Suppose
$(\rho,\rho_{0},w,\eta)$ satifies

(4.4) $\rho-\Delta w=\rho_{0}+\eta,\quad
w(1-\rho)=0,\quad\rho\in\\{0,1\\},\quad w\geq
0,\quad\partial_{t}\eta=(n_{0}-\eta)\rho.$

Then for all $t\in(0,T]$,

$\\{x\in\mathbb{R}^{d}:\rho(x,t)=1\\}=\\{x\in\mathbb{R}^{d}:w(x,t)>0\\}\cup\\{x\in\mathbb{R}^{d}:\rho_{0}(x)=1,n_{0}(x)=0\\},$

up to a measure zero set in $\mathbb{R}^{d}$.

###### Proof.

Thanks to the regularity of $w$, $\Delta w$ vanishes almost everywhere on the
set where $w$ vanishes. Thus, combining the elliptic equation and the relation
$w(1-\rho)=0$, it follows that for every $t\in(0,T]$

(4.5) $\rho(x,t)=\chi_{w}(x,t)+(1-\chi_{w}(x,t))(\rho_{0}(x)+\eta(x,t))$

almost everywhere, where $\chi_{w}$ is the characteristic function of the set
$\\{w>0\\}$.

Using Lemma 4.5 we see that $\eta$ solves $\partial_{t}\eta=n_{0}\rho-\eta$.
Therefore, we have the following integral representation for $\eta$:

$\begin{split}\eta(x,t)=&\;n_{0}(x)\int_{0}^{t}\rho(x,s)e^{s-t}\,ds\\\
=&\;n_{0}(x)\int_{0}^{t}\Big{(}\chi_{w}(x,s)+(1-\chi_{w}(x,s))\big{(}\rho_{0}(x)+\eta(x,s)\big{)}\Big{)}e^{s-t}\,ds.\end{split}$

Define

$h(x,t):=\eta(x,t)(1-\chi_{w}(x,t)).$

Since $\chi_{w}$ is increasing in time, it follows that
$\chi_{w}(x,s)(1-\chi_{w}(x,t))=0$ and
$(1-\chi_{w}(x,s))(1-\chi_{w}(x,t))=(1-\chi_{w}(x,s))$ whenever $s\leq t$.
Therefore, multiplying the integral representation of $\eta$ by
$(1-\chi_{w}(x,t))$, we obtain that

$\begin{split}h(x,t)=&\;n_{0}(x)\int_{0}^{t}(1-\chi_{w}(x,s))\big{(}\rho_{0}(x)+\eta(x,s)\big{)}e^{s-t}\,ds\\\
=&\;n_{0}(x)\int_{0}^{t}\big{(}(1-\chi_{w}(x,s))\rho_{0}(x)+h(x,s)\big{)}e^{s-t}\,ds\end{split}$

holds almost everywhere. For almost every $x$ such that
$n_{0}(x)\rho_{0}(x)=0$, Gronwall’s inequality implies that $h(x,t)=0$ for all
$t\in[0,T]$.

By (4.5) and the definition of $h$, we may write
$\\{(x,t):\rho(x,t)=1\\}=\\{(x,t):w(x,t)>0\\}\cup\\{(x,t):\rho_{0}(x)+h(x,t)=1\\}$
up to a measure-zero set. From our work above, when $n_{0}(x)\rho_{0}(x)=0$,
$h(x,t)=0$ almost everywhere. When $n_{0}(x)\rho_{0}(x)\neq 0$ we have
$\rho_{0}(x)=1$ and $n_{0}(x)>0$, which implies that $\rho_{0}(x)+\eta(x,t)>1$
for all $t>0$ and it is non-decreasing. For almost every such $x$, we must
have $w(x,t)>0$ for all $t>0$, as otherwise (4.4) would not hold.

This completes the proof. ∎

Next, we state two results on the stability of the full system (P). They
respectively extend Proposition 4.3 and Proposition 4.4.

######  Theorem 4.7.

Suppose that $(\rho^{0},p^{0},\eta^{0})$ and $(\rho^{1},p^{1},\eta^{1})$ are
solutions to the PDE (P) starting from initial datum
$(\rho_{0}^{0},n_{0}^{0})$ and $(\rho_{0}^{1},n_{0}^{1})$, respectively. Fix a
domain $D\subset\mathbb{R}^{d}$. Suppose $\rho_{0}^{0}$ and $\rho_{0}^{1}$ are
patches, $n^{0}_{0}\in L^{\infty}(D)$, and $p^{0}=p^{1}$ on $\partial
D\times[0,T]$. Denote

$N_{D}(t)=\int_{0}^{t}e^{(\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}-1)\tau}\,d\tau=\begin{cases}\frac{e^{(\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}-1)t}-1}{\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}-1}&\textup{if}\;\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}\neq 1,\\\ t&\textup{otherwise},\end{cases}$

and

$M_{D}(t)=\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}N_{D}(t)+1=\begin{cases}\frac{\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}e^{(\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}-1)t}-1}{\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}-1}&\textup{if}\;\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}\neq 1,\\\ t+1&\textup{otherwise}.\end{cases}$

Then for all $t\in[0,T]$,

$\lVert(\rho^{1}-\rho^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\leq
N_{D}(t)\lVert(n^{1}_{0}-n^{0}_{0})_{+}\rVert_{L^{1}(D)}+M_{D}(t)\lVert(\rho_{0}^{1}-\rho_{0}^{0})_{+}\rVert_{L^{1}(D)}.$

If in addition, $\rho_{0}^{0}\leq\rho_{0}^{1}$ and $n_{0}^{0}\leq n_{0}^{1}$
for almost every $x\in D$, then for all $t\in[0,T]$ we have
$\rho^{0}(x,t)\leq\rho^{1}(x,t)$ for almost every $x\in D$.

######  Remark 4.8.

Compared with the proof of Proposition 3.13, this theorem provides an improved
estimate by taking advantage of the conditions $b=0$ and $\rho^{i}$ $(i=0,1)$
are patch solutions.

###### Proof.

Thanks to Proposition 4.3, for every $t\in[0,T]$ we have the contraction
inequality

$\lVert(\rho^{1}-\rho^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\leq\lVert(\rho^{1}_{0}-\rho_{0}^{0})_{+}\rVert_{L^{1}(D)}+\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}.$

On the other hand, by direct computation and Lemma 4.5,

$\begin{split}&\;\frac{d}{dt}\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}+\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\\\
\leq&\;\lVert(n^{1}_{0}\rho^{1}-n^{0}_{0}\rho^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\\\
\leq&\;\lVert(n_{0}^{1}-n^{0}_{0})_{+}\rVert_{L^{1}(D)}+\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}\lVert(\rho^{1}-\rho^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}.\end{split}$

Plugging in the contraction inequality, we see that

$\begin{split}\frac{d}{dt}\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\leq&\;\big{(}\lVert
n^{0}_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1\big{)}\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\\\
&\;+\lVert(n_{0}^{1}-n^{0}_{0})_{+}\rVert_{L^{1}(D)}+\lVert
n^{0}_{0}\rVert_{L^{\infty}(D)}\lVert(\rho^{1}_{0}-\rho_{0}^{0})_{+}\rVert_{L^{1}(D)}.\end{split}$

Then Gronwall’s inequality implies that

$\lVert(\eta^{1}-\eta^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\leq
N_{D}(t)\Big{(}\lVert(n^{1}_{0}-n^{0}_{0})_{+}\rVert_{L^{1}(D)}+\lVert
n^{0}_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert(\rho^{1}_{0}-\rho_{0}^{0})_{+}\rVert_{L^{1}(D)}\Big{)},$

and thus

$\lVert(\rho^{1}-\rho^{0})_{+}\rVert_{L^{1}(D\times\\{t\\})}\leq
N_{D}(t)\lVert(n_{0}^{1}-n^{0}_{0})_{+}\rVert_{L^{1}(D)}+M_{D}(t)\lVert(\rho^{1}_{0}-\rho_{0}^{0})_{+}\rVert_{L^{1}(D)}.$

The final comparison property is automatic from the above estimates. ∎

######  Proposition 4.9.

Let $(\rho,p,\eta)$ be a solution to (P) starting from initial data $\rho_{0}$
and $n_{0}$. Given an antisymmetric matrix $A\in\mathbb{R}^{d\times d}$ and
some $g\in L^{1}(\mathbb{R}^{d})$ define

$\lVert g\rVert_{BV_{A}(\mathbb{R}^{d})}=\sup_{\varphi\in
C^{\infty}_{c}(\mathbb{R}^{d}),\lVert\varphi\rVert_{L^{\infty}(\mathbb{R}^{d})}\leq
1}\int_{\mathbb{R}^{d}}g(x)Ax\cdot\nabla\varphi\,dx.$

If $\rho_{0}\in BV(\mathbb{R}^{d})$ is a patch and $n_{0}\in
L^{\infty}(\mathbb{R}^{d})\cap BV(\mathbb{R}^{d})$, then for all $t\in[0,T]$,

$\lVert\rho\rVert_{BV(\mathbb{R}^{d}\times\\{t\\})}\leq N(t)\lVert
n_{0}\rVert_{BV(\mathbb{R}^{d})}+M(t)\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})},$

and

$\lVert\rho\rVert_{BV_{A}(\mathbb{R}^{d}\times\\{t\\})}\leq N(t)\lVert
n_{0}\rVert_{BV_{A}(\mathbb{R}^{d})}+M(t)\lVert\rho_{0}\rVert_{BV_{A}(\mathbb{R}^{d})}.$

Here $N(t):=N_{\mathbb{R}^{d}}(t)$ and $M(t):=M_{\mathbb{R}^{d}}(t)$ are
defined as in Theorem 4.7.

###### Proof.

Proposition 4.4 immediately gives us the bound

$\lVert\rho(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\leq\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert\eta(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}.$

Using (4.3), it is clear that $\rho,\eta\in
L^{\infty}([0,T];BV(\mathbb{R}^{d}))$. To obtain a better bound, we
differentiate the $\eta$-equation

$\partial_{t}\nabla\eta=\rho\nabla n_{0}+n_{0}\nabla\rho-\nabla\eta,$

which allows us to estimate

$\begin{split}\frac{d}{dt}\big{(}e^{t}\lVert\eta(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\big{)}\leq&\;e^{t}\left(\lVert
n_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\rho(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\right)\\\
\leq&\;e^{t}\left(\lVert n_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\eta(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\right).\end{split}$

Rearranging this yields

$\frac{d}{dt}\big{(}e^{(1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})})t}\lVert\eta(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\big{)}\leq
e^{(1-\lVert n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})})t}\Big{(}\lVert
n_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})}\Big{)}.$

Thus,

$\lVert\eta(\cdot,t)\rVert_{BV(\mathbb{R}^{d})}\leq N(t)\Big{(}\lVert
n_{0}\rVert_{BV(\mathbb{R}^{d})}+\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\lVert\rho_{0}\rVert_{BV(\mathbb{R}^{d})}\Big{)}.$

The $BV$-estimate now follows. The estimate for the $BV_{A}$-norm can be
proved analogously. ∎

By the comparison principle Theorem 4.7, one can see that any solution
starting with an initial nutrient such that $\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}<1$ must have bounded mass for all
time. Indeed, we may take in Theorem 4.7 that $\rho_{0}^{0}=0$,
$\rho_{0}^{1}=\rho_{0}$, and $n_{0}^{0}=n_{0}^{1}=n_{0}$. In such case, the
solution approaches a stationary state as $t\to\infty$. The following theorem
characterizes the stationary state and provides a convergence rate.

######  Theorem 4.10.

Suppose $\lVert n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}<1$ and
$\rho_{0}\in\\{0,1\\}$ and that $(\rho,p,\eta)$ is a solution to the PDE (P)
starting from the initial data $(\rho_{0},n_{0})$. Let $\rho_{\infty}$ solve
the elliptic equation

(4.6) $(1-n_{0}(x))\rho_{\infty}-\Delta w_{\infty}=\rho_{0},\quad
w_{\infty}(1-\rho_{\infty})=0,\quad\rho_{\infty}\in[0,1],\quad w_{\infty}\geq
0.$

Then we have

$\displaystyle\lVert\rho(\cdot,t)-\rho_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq$
$\displaystyle\;\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t},$
$\displaystyle\lVert\eta(\cdot,t)-n_{0}\rho_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq$
$\displaystyle\;\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t}.$

###### Proof.

Recall that $\eta=\rho\eta$, since $\rho$ is a patch solution. A direct
computation shows that

$\int_{\mathbb{R}^{d}}\partial_{t}\rho\,dx=\int_{\mathbb{R}^{d}}n_{0}\rho-\eta\,dx,\quad\int_{\mathbb{R}^{d}}\partial_{tt}\rho\,dx=\int_{\mathbb{R}^{d}}n_{0}\partial_{t}\rho-(n_{0}\rho-\eta)\,dx.$

Therefore,

$\int_{\mathbb{R}^{d}}\partial_{tt}\rho\,dx+\big{(}1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}\big{)}\int_{\mathbb{R}^{d}}\partial_{t}\rho\,dx\leq
0.$

Hence,

$\int_{\mathbb{R}^{d}}\partial_{t}\rho(x,t)\,dx\leq e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t}\int_{\mathbb{R}^{d}}\rho_{0}n_{0}\,dx.$

A similar calculation gives

$\int_{\mathbb{R}^{d}}\partial_{t}\eta(x,t)\,dx\leq e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t}\int_{\mathbb{R}^{d}}\rho_{0}n_{0}\,dx.$

Since $\partial_{t}\rho,\partial_{t}\eta\geq 0$ almost everywhere, the above
bounds imply that there exists $\rho_{\infty},\eta_{\infty}\in
L^{1}(\mathbb{R}^{d})$ such that

$\displaystyle\lVert\rho(\cdot,t)-\rho_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t},$
$\displaystyle\lVert\eta(\cdot,t)-\eta_{\infty}\rVert_{L^{1}(\mathbb{R}^{d})}\leq\frac{\lVert
n_{0}\rho_{0}\rVert_{L^{1}(\mathbb{R}^{d})}}{1-\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}}e^{(\lVert
n_{0}\rVert_{L^{\infty}(\mathbb{R}^{d})}-1)t}.$

Finally, it is clear that
$0=\lim_{t\to\infty}\partial_{t}\eta(\cdot,t)=\rho_{\infty}n_{0}-\eta_{\infty}.$
Hence, $\eta_{\infty}=n_{0}\rho_{\infty}$. Since all of the variables are
increasing with respect to time, we can pass to the limit in the elliptic
equation (4.2) to obtain (4.6). ∎

##  5\. Regularity of the Tumor Patch

In this section we analyze the free boundary regularity of the tumor region
$\\{\rho=1\\}$. For simplicity we only consider patch solutions $\rho$ of the
tumor growth model (P).

When $\rho_{0}$ equals the characteristic function $\chi_{\Omega_{0}}$ with
$\Omega_{0}$ being a compact set, Proposition 4.3 and Lemma 4.5 imply that for
any $t>0$, $\rho(\cdot,t)=\chi_{\Omega_{t}}$ for some $\Omega_{t}$ that
increases in time. Lemma 4.6 gives a characterization of $\Omega_{t}$. We will
focus on the boundary regularity of $\Omega_{t}$ under suitable geometric
assumptions.

###  5.1. Reflection geometry

First, we recall from [FK14] that sets having _reflection geometries_ , i.e.,
those satisfying ordering properties when reflected with respect to a family
of hyperplanes, have locally Lipschitz boundaries. For both simplicity and
relevance to the case $n_{0}(x)\geq 1$, we focus on isotropic reflection
geometry defined as _$r$ -reflection property_ below. For more general types
of reflection geometries that ensure Lipschitz regularity of the set boundary,
see [KKP21].

For a hyperplane $H$ in $\mathbb{R}^{d}$ with unit normal vector $\nu_{H}$,
the reflection with respect to $H$ is given by

$\phi_{H}(x):=x-2(x-y,\nu_{H})\nu_{H}\hbox{ for some }y\in H.$

######  Definition 5.1 (Definition 3.10, [FK14]).

For $r>0$, we say a bounded open set $\Omega\subset\mathbb{R}^{d}$ satisfies
$r$-reflection property if, for any hyperplane $H$ such that
$H^{-}:=\\{x:(x-y)\cdot\nu_{H}<0,\,y\in H\\}$ contains $B_{r}(0)$,

$\phi_{H}(\Omega\cap H^{+})\subset\Omega\cap H^{-}.$

Here $H^{+}:=\phi_{H}(H^{-})$.

######  Lemma 5.2 (Lemma 3.24, [FK14]).

Let $C(x,\theta):=\\{y:\langle x,y\rangle\geq\cos\theta|x||y|\\}$ denote the
cone with direction $x$ and angle $\theta\in[0,\pi/2]$. Suppose $\Omega$
satisfies $r$-reflection property, and that $\Omega$ contains the closure of
$B_{r}(0)$. Then for all $x\in\partial\Omega$, there is an exterior cone
$C(x,\phi_{x})$ to $\Omega$ at $x$, such that

$x+C(x,\phi_{x})\subset\Omega^{c},\hbox{ where }\cos\phi_{x}=\frac{r}{|x|}.$

######  Definition 5.3.

Let $\mathcal{P}$ be a family of hyperplanes in $\mathbb{R}^{d}$ that do not
go through the origin. $u$ is said to have reflection properties with respect
to $\mathcal{P}$, if for any $H\in\mathcal{P}$, $u_{H}\leq u$ in $H^{+}$. Here
$u_{H}(x):=u\circ\phi_{H}(x)$ is the reflection of $u$ with respect to $H$,
and $H^{+}$ is the half-space generated by $H$ that contains the origin.

Lemma 5.2 implies that, if $u$ has reflection properties with respect to
$\mathcal{P}$ and there is enough of these hyperplanes in $\mathcal{P}$, then
the super level sets of $u$ should have locally Lipschitz boundary.

###  5.2.  Lipschitz regularity for the tumor patch

Using the invariance of the system with respect to reflections, we will show
that, under suitable assumptions on $n_{0}$, if $r$-reflection property holds
for $\Omega_{0}$ initially, it remains to be true for $\Omega_{t}$ for all
$t>0$ (see Proposition 5.4). This allows us to prove that $\partial\Omega_{t}$
enjoys Lipschitz regularity in various typical cases. For example, if
$n_{0}(x)\geq 1$ on $\mathbb{R}^{d}$, we will show that for any initial tumor
region $\Omega_{0}$, $\Omega_{t}$ has locally Lipschitz boundary for all
sufficiently large times: see Corollary 5.5. We will build on the Lipschitz
regularity to later show $C^{1,\alpha}$-regularity of the free boundary, using
the elliptic equation that $w$ solves.

We first prove the so-called reflection comparison result.

######  Proposition 5.4 (Reflection comparison).

For a given hyperplane $H$ in $\mathbb{R}^{d}$, define
$\rho_{H}:=\rho\circ\phi_{H}$. Let $H^{+}$ be one of the half-spaces generated
by $H$. If $(\rho_{0})_{H}\leq\rho_{0}$ and $(n_{0})_{H}\leq n_{0}$ a.e. in
$H^{+}$, then $(\rho)_{H}\leq\rho$, and $(\eta)_{H}\leq\eta$ a.e. in
$H^{+}\times[0,\infty)$.

###### Proof.

Let $w_{H}:=w\circ\phi_{H}$ and there holds trivially $w|_{H}=w_{H}|_{H}$.
Thanks to the symmetries of the Laplacian, we have

$\rho_{H}-\Delta w_{H}=(\rho_{0})_{H}+\eta_{H},\quad
w_{H}(1-\rho_{H})=0,\quad\rho_{H}\in[0,1],\quad w_{H}\geq 0$

almost everywhere in $H^{+}$. Thus, the stability estimate in Proposition 4.3
implies that for every $t\in[0,T]$,

$\lVert(\rho_{H}-\rho)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}\leq\lVert(\eta_{H}-\eta)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}.$

Since $\partial_{t}\eta=n_{0}\rho-\eta$, it follows that
$\partial_{t}\eta_{H}=(n_{0})_{H}\rho_{H}-\eta_{H}$. Thus,

$\begin{split}&\;\frac{d}{dt}\lVert(\eta_{H}-\eta)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}+\lVert(\eta_{H}-\eta)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}\\\
\leq&\;\lVert((n_{0})_{H}\rho_{H}-n_{0}\rho)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}\\\
\leq&\;\lVert(n_{0})_{H}\rVert_{L^{\infty}(H^{+})}\lVert(\rho_{H}-\rho)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}\\\
\leq&\;\lVert(n_{0})_{H}\rVert_{L^{\infty}(H^{+})}\lVert(\eta_{H}-\eta)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}.\end{split}$

Now it follows from Gronwall’s inequality that, for all $t\in[0,T]$,
$\lVert(\eta_{H}-\eta)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}=0$, and thus
$\lVert(\rho_{H}-\rho)_{+}\rVert_{L^{1}(H^{+}\times\\{t\\})}=0$. ∎

The following is a consequence of the reflection comparison and Lemma 5.2.

######  Corollary 5.5.

Assume that the super-level sets of $n_{0}$ satisfy the $r$-reflection for
some $r>0$. Then the followings hold:

  1. (a)

Suppose $n_{0}(x)<1$ on $\mathbb{R}^{d}$. If $\Omega_{0}$ satisfies
$r$-reflection, then so does $\Omega_{t}$ for all $t>0$. In this case, suppose
$\overline{B_{R}(0)}\subset\\{w_{\infty}>0\\}$ for some $R>r$, where
$w_{\infty}$ is defined in Theorem 4.10. Then there is $T>0$ such that
$\Omega_{t}$ contains $B_{R}(0)$ whenever $t>T$. Consequently, for $t>T$,
$\Omega_{t}$ has uniformly Lipschitz boundary, with Lipschitz constant less
than $O(\frac{r}{R-r})$.

  2. (b)

Suppose $n_{0}(x)\geq 1$ on $\mathbb{R}^{d}$ and that $\Omega_{0}$ is a
bounded open set contained in $B_{r}(0)$. Then for any $R>0$,
$T(R):=\inf\\{t:B_{R}(0)\subset\Omega_{t}\\}$ is finite. Consequenctly, for
any $R>r$, when $t\geq T(R)$, the set $\Omega_{t}$ has Lipschitz boundary with
respect to radial direction with Lipschitz constant less than
$O(\frac{r}{R-r})$.

######  Remark 5.6.

Note that the $r$-reflection condition does not restrict the shape of either
$\Omega_{0}$ or super-level sets of $n_{0}$ inside $B_{r}(0)$. Hence, we can
start with any initial data $\Omega_{0}\subset B_{r}(0)$ in both cases of the
above corollary, where the evolution of the set may go through topological
singularities such as merging of the free boundaries. The above results state
that, given that the initial nutrient $n_{0}$ is “well-prepared” outside
$B_{r}(0)$, once $\partial\Omega_{t}$ moves outside $B_{r}(0)$, there will be
no further topological changes in the evolution, and $\partial\Omega_{t}$
remains being Lipschitz.

######  Remark 5.7.

One could relax the assumption on $n_{0}$ so that for level sets lying between
$B_{R}(0)$ and $B_{2R}(0)$, the corresponding super-level set only satisfies
$R$-reflection. That would allow us to study the case, e.g., where level sets
of $n_{0}$ are ellipses. That may admit possibly non-radial asymptotic shapes
for $\Omega_{t}$. We leave such discussion to interested readers.

###### Proof.

For (a), $\Omega_{t}$ satisfies $r$-reflection due to Proposition 5.4.
Moreover, both $\rho$ and $\eta$ monotone increases in time and converge to
$\rho_{\infty}$ and $n_{0}\rho_{\infty}$ in $L^{1}(\mathbb{R}^{d})$, due to
Theorem 4.10. In the patch case we have
$\rho_{\infty}=\chi_{\Omega_{\infty}}$, where $\Omega_{\infty}$ is bounded.
This implies the above convergence also holds in $L^{p}(\mathbb{R}^{d})$ for
any $p\in[1,+\infty)$. It then follows from (4.2) and (4.6) that $w$ then
uniformly converges to $w_{\infty}$ as $t\to\infty$. Hence if we know that
$\overline{B_{R}(0)}\subset\\{w_{\infty}>0\\}$, from the fact that
$w(\cdot,t)$ monotonically increases to converge to $w_{\infty}$, we conclude
that $\Omega_{t}$ contains $\overline{B_{r}(0)}$ for sufficiently large $t$.
Then we can conclude (a) by Lemma 5.2.

For (b), that $T(R)$ being finite follows from comparison with the case
$n_{0}(x)\equiv 1$; see Theorem 4.7 and Remark 6.3 below. Also note that for
any hyperplane $H$ such that $H^{-}$ contains $B_{r}(0)$, $\Omega_{0}\cap
H^{+}$ is an empty set, and thus its reflected image $\phi_{H}(\Omega_{0}\cap
H^{+})$ is trivially contained by $\Omega_{0}\cap H^{-}$. Thus $\Omega_{t}$
satisfies $r$-reflection for all $t\geq 0$, and once $t\geq T(R)$ with $R>r$,
we can apply Lemma 5.2. ∎

Let us now assume that, for a domain
$D\subset\mathbb{R}^{n}\setminus\Omega_{0}$, $\partial\Omega_{t}\cap D$ is
non-empty but has Lebesgue measure zero. Then Lemma 4.6 yields that
$\rho(\cdot,t)=\chi_{\\{w(\cdot,t)>0\\}}$ in $D$, and thus (4.2) can be
written as

$\Delta w=\rho-\eta=\chi_{\\{w>0\\}}-\eta\hbox{ in }D.$

Moreover, by its definition $\eta=n_{0}-n$ is supported only in
$\chi_{\\{w>0\\}}$. Hence, (4.1)-(4.2) yields the following $(w,\eta)$-system
in $D$:

(W)
$\left\\{\begin{array}[]{lll}\eta_{t}=(n_{0}-\eta)\rho=n_{0}\chi_{\\{w>0\\}}-\eta;\\\
\\\ \chi_{\\{w>0\\}}-\Delta w=\eta\chi_{\\{w>0\\}}.\end{array}\right.$

In what follows, we will explore the second equation in (W) to prove further
regularity of the free boundary, in the setting given in Corollary 5.5.

###  5.3. $C^{1,\alpha}$-regularity of the free boundary

According to (W), away from the support of $\rho_{0}$, $w$ solves

$(1-\eta)\chi_{\\{w>0\\}}-\Delta w=0.$

As long as $\eta<1$ near $\partial\\{w>0\\}$ and is $C^{\alpha}$, this problem
falls into the category of standard obstacle problem, whose singular points
feature a blow-up profile of a quadratic polymonial with sub-quadratic error
term. This is impossible if the set $\\{w>0\\}$ is known to have locally
Lipschitz free boundary in $D$. Hence, if we know a priori that $\Omega_{t}$
has Lipschitz boundary, then the free boundary $\partial\\{w>0\\}$ consists of
only regular points, i.e., it is $C^{1,\alpha}$; see [Bla01, Theorem 7.2].
This would be the conclusion we will obtain at the end of this section, in
Corollary 5.13, with a class of initial data discussed in Corollary 5.5.

To complete this argument, we need to prove that $\eta$ is indeed Hölder
continuous in space. Note that $\eta(x,t)$ starts evolving in time once the
set $\Omega_{t}$ reaches $x$. Regularity of $\eta$ thus is directly related to
the dynamics of the set $\Omega_{t}$, or equivalently that of
$\\{w(\cdot,t)>0\\}$. We will show a version of non-degeneracy for the
pressure variable near the boundary (see Proposition 5.12), which further
implies non-degeneracy of the propagation speed of the boundary, ensuring that
the tumor patch reaches nearby points with small time difference.

To study the dynamics of the tumor patch, we will need a direct comparison
principle for $(\rho,p)$, with barriers with a fixed $n$, as follows.

######  Lemma 5.8.

Let $(\rho,p,n)$ solve the original system (P), and suppose
$(\bar{\rho},\bar{p})$ weakly solves

$\partial_{t}\bar{\rho}-\Delta\bar{p}\leq n\bar{\rho}\;\hbox{ in
}D\times[t_{0},t_{1}],\quad$

with $\bar{p}\in L^{2}([0,T];H^{1}(\mathbb{R}^{d}))$ such that
$\bar{\rho}=\operatorname{sgn}_{+}(\bar{p})$. If $\bar{\rho}\leq\rho$ at
$t=t_{0}$, $\bar{p}\leq p$ on $\partial D\times[t_{0},t_{1}]$, and
$\partial_{t}\bar{\rho}\geq 0$, then $\bar{\rho}\leq\rho$ and $\bar{p}\leq p$
almost everywhere on $D\times[t_{0},t_{1}]$.

######  Remark 5.9.

We expect the result to still hold if one drops the requirement that
$\partial_{t}\bar{\rho}\geq 0$ and
$\bar{\rho}=\operatorname{sgn}_{+}(\bar{p})$. Nonetheless, these assumptions
make the proof substantially easier and the above statement is sufficiently
strong for our purposes.

######  Remark 5.10.

Let us note that this result is closely related to the comparison and
uniqueness statements Lemma 3.11 and Proposition 3.13. Nonetheless, neither
statement directly applies to prove the above result. Lemma 3.11 assumes that
one already has an ordering for the source terms, while Proposition 3.13
proves uniqueness through stability rather than comparison. On the other hand,
it is very likely that one could obtain a strengthened version of this result
by appropriately tweaking the argument in Lemma 3.11.

###### Proof.

Let $\bar{w}(x,t)=\int_{t_{0}}^{t}\bar{p}(x,s)\,ds$. Since $\bar{\rho}$ is
increasing and $\bar{\rho}=\operatorname{sgn}_{+}(\bar{p})$, it follows that
$\bar{\rho}=\operatorname{sgn}_{+}(\bar{w})$ for any $t>t_{0}$. Hence,
integrating in time, we see that $\bar{w}$ weakly solve

$\operatorname{sgn}_{+}(\bar{w})(x,t)-\Delta\bar{w}(x,t)\leq\bar{\rho}(x,t_{0})+\int_{t_{0}}^{t}n\operatorname{sgn}_{+}(\bar{w})(x,s)\,ds.$

If we set $w(x,t)=\int_{t_{0}}^{t}p(x,s)\,ds$ with abuse of notations, we
recall that the original system solves

$\rho(x,t)-\Delta w(x,t)=\rho(x,t_{0})+\int_{t_{0}}^{t}n\rho(x,s)\,ds,\quad
w(1-\rho)=0,\quad\rho\in[0,1],\quad w\geq 0$

almost everywhere in $D\times[t_{0},t_{1}]$. Then we argue as in the proof of
Proposition 4.3. Fix $\delta>0$ and let $f_{\delta}:\mathbb{R}\to\mathbb{R}$
be a smooth increasing function such that $f_{\delta}(a)=0$ if $a\leq 0$ and
$f_{\delta}(a)=1$ if $a\geq\delta$. Fix some time $t>t_{0}$. Taking the
difference of the above two formulas and integrating against
$f_{\delta}(\bar{w}-w)$ along $D\times[t_{0},t]$, we see that

$\begin{split}&\;\int_{D\times[t_{0},t]}(\operatorname{sgn}_{+}(\bar{w})-\rho)f_{\delta}(\bar{w}-w)+f^{\prime}_{\delta}(\bar{w}-w)|\nabla(\bar{w}-w)|^{2}\,dx\,ds\\\
\leq&\;\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}\int_{t_{0}}^{t}\lVert(\operatorname{sgn}_{+}(\bar{w})-\rho)_{+}\rVert_{L^{1}(D\times[t_{0},s])}\,ds,\end{split}$

where we have used the fact that $f_{\delta}\in[0,1]$ almost everywhere, that
$f_{\delta}(\bar{w}-w)$ vanishes almost everywhere on $\partial
D\times[t_{0},t_{1}]$, and that $\bar{\rho}\leq\rho$ at $t=t_{0}$. Since
$\operatorname{sgn}_{+}(\bar{w})$ and $\rho$ only take the values $0$ or $1$
almost everywhere, we have

$(\operatorname{sgn}_{+}(\bar{w})-\rho)_{+}=(\operatorname{sgn}_{+}(\bar{w})-\rho)\operatorname{sgn}_{+}(\bar{w}-w).$

Hence, we can send $\delta\to 0$ in the above inequality to obtain

$\lVert(\operatorname{sgn}_{+}(\bar{w})-\rho)_{+}\rVert_{L^{1}(D\times[t_{0},t])}\leq\|n_{0}\|_{L^{\infty}(\mathbb{R}^{d})}\int_{t_{0}}^{t}\lVert(\operatorname{sgn}_{+}(\bar{w})-\rho)_{+}\rVert_{L^{1}(D\times[t_{0},s])}\,ds.$

Now Gronwall’s inequality implies that
$\bar{\rho}=\operatorname{sgn}_{+}(\bar{w})\leq\rho$ almost everywhere in
$D\times[t_{0},t_{1}]$.

It remains to prove $\bar{p}\leq p$ almost everywhere in
$D\times[t_{0},t_{1}]$. Choose some $\psi\in L^{2}([t_{0},t_{1}];H^{1}(D))$
such that $\psi\geq 0$ almost everywhere in $D\times[t_{0},t_{1}]$, $\psi=0$
at $t=t_{1}$, $\psi(1-\bar{\rho})=0$, and $\psi|_{\partial
D\times[t_{0},t_{1}]}=0$. If we fix $\epsilon>0$ and let
$\psi^{\epsilon}(x,t)=\epsilon^{-1}\int_{t}^{\min(t_{1},t+\epsilon)}\psi(x,s)\,ds$,
then $\psi^{\epsilon}$ is a valid test function for the weak equation
$\partial_{t}\bar{\rho}-\Delta\bar{p}\leq\bar{\rho}n$ since
$\psi^{\epsilon}=0$ at $t=t_{1}$. Hence,

$\int_{D\times[t_{0},t_{1}]}-\bar{\rho}\partial_{t}\psi^{\epsilon}+\nabla\psi^{\epsilon}\cdot\nabla\bar{p}-\psi^{\epsilon}\bar{\rho}n\,dx\,dt\leq\int_{D}\bar{\rho}(x,t_{0})\psi^{\epsilon}(x,t_{0})\,dx.$

Note that for almost every $t\in[t_{0},t_{1}]$,

$\bar{\rho}\partial_{t}\psi^{\epsilon}=\bar{\rho}(x,t)\frac{\psi(x,\min\\{t_{1},t+\epsilon\\})-\psi(x,t)}{\epsilon}\leq\frac{\psi(x,\min\\{t_{1},t+\epsilon\\})-\psi(x,t)}{\epsilon}=\partial_{t}\psi^{\epsilon},$

where we use $\psi(1-\bar{\rho})=0$ and the non-negativity of
$\psi^{\epsilon}$ to justify the inequality. Therefore,

$\int_{D\times[t_{0},t_{1}]}\nabla\psi^{\epsilon}\cdot\nabla\bar{p}-\psi^{\epsilon}\bar{\rho}n\,dx\,dt\leq\int_{D}(\bar{\rho}(x,t_{0})-1)\psi^{\epsilon}(x,t_{0})\,dx\leq
0.$

Sending $\epsilon\to 0$ and once again using $\psi(1-\bar{\rho})=0$, we see
that

$\int_{D\times[t_{0},t_{1}]}\nabla\psi\cdot\nabla\bar{p}-\psi n\,dx\,dt\leq
0.$

Since $\bar{\rho}\leq\rho$, we also have $\psi(1-\rho)=0$. Therefore the
complementarity condition, Proposition 3.7, implies that

$\int_{D\times[t_{0},t_{1}]}\nabla\psi\cdot\nabla p-\psi n\,dx\,dt=0.$

Combining the above two formulas yields

$\int_{D\times[t_{0},t_{1}]}\nabla\psi\cdot\nabla(\bar{p}-p)\,dx\,dt\leq 0.$

Choose $\psi=\omega(t)(\bar{p}-p)_{+}$, where $\omega(t)$ is a smooth non-
negative function such that $\omega>0$ on $[t_{0},t_{1})$ and
$\omega(t_{1})=0$. It is clear that this choice satisfies our assumptions on
$\psi$. Hence, we obtain that

$\int_{D\times[t_{0},t_{1}]}\omega(t)|\nabla(\bar{p}-p)_{+}|^{2}\,dx\,dt\leq
0.$

Since $\bar{p}\leq p$ almost everywhere on $\partial D\times[t_{0},t_{1}]$, we
obtain $\bar{p}\leq p$ almost everywhere on $D\times[t_{0},t_{1}]$. ∎

Based on above comparison principle, we will build a radial barrier
$(\bar{\rho},\bar{p})$ to compare with $(\rho,p)$, to show that the pressure
support spreads with a uniform rate. To construct the barrier it is useful to
recall Dahlberg’s lemma:

######  Lemma 5.11 (Dahlberg’s lemma, [WW79]).

Let $u_{1},u_{2}$ be two non-negative harmonic functions in
$D\subset\mathbb{R}^{d}$ of the form

$D=\big{\\{}(x^{\prime},x_{n})\in\mathbb{R}^{n-1}\times\mathbb{R}:|x^{\prime}|<2,|x_{n}|<2M,x_{n}>f(x^{\prime})\big{\\}}$

with a Lipschitz function $f$ with Lipschitz constant less than $M$ and
$f(0)=0$. Assume further that $u_{1}=u_{2}=0$ along the graph of $f$. Then
there exist constants $C_{1},C_{2}$ only depending on $M$, such that

$0<C_{1}\leq\dfrac{u_{1}(x^{\prime},x_{n})}{u_{2}(x^{\prime},x_{n})}\cdot\dfrac{u_{2}(0,M)}{u_{1}(0,M)}\leq
C_{2}$

in the smaller domain

$D_{1/2}=\big{\\{}|x^{\prime}|<1,|x_{n}|<M,x_{n}>f(x^{\prime})\big{\\}}.$

Let us now define $C_{\theta}:=C(-e_{n},\theta)=\\{x_{n}<f(x^{\prime})\\}$
with $f(x^{\prime}):=-\cot\theta|x^{\prime}|$, and let $h$ solve

$-\Delta h=\delta_{\\{x=-3e_{n}\\}}\hbox{ in }C_{\theta},\quad h=0\hbox{ on
}\\{x_{n}=f(x^{\prime})\\}.$

Then $h$ has a polynomial growth from the origin, namely,

(5.1) $\begin{split}&\mbox{there exists some }k=k(\theta)>0\mbox{ that
decreases in }\theta,\mbox{ such that }\\\
&h(x)\simeq(f(x^{\prime})-x_{n})_{+}^{k}\hbox{ in }C_{\theta}\cap
B_{2}(0).\end{split}$

Let $\theta_{d}$ be the angle such that $k(\theta_{d})=2$, where $h$ has
quadratic growth near the origin. For instance, $\theta_{2}=\pi/4$, and
$\theta_{d}$ increases as $d$ increases.

######  Proposition 5.12 (Nondegeracy of the pressure).

Suppose $x_{0}\in\partial\\{p(\cdot,t_{0})>0\\}$, and suppose that
$\\{p(\cdot,t_{0})>0\\}$ contains $(x_{0}+C_{\theta})\cap B_{2}(x_{0})$ with
$\theta>\theta_{d}\geq\pi/4$. Let
$\bar{n}:=\min\\{n(x,t_{0}):|x-(x_{0}-2e_{n})|<1\\}$. Then there exists a
dimensional constant $C>0$ such that the following holds: for any $0<r<1/C$,
we have

$B_{r}(x_{0})\subset\left\\{p\left(\cdot,t_{0}+\frac{Cr^{\alpha}}{\bar{n}}\right)>0\right\\},$

where $\alpha=2-k(\theta)>0$, with $k(\theta)$ given in (5.1).

###### Proof.

Throughout the proof, we use $C$ to denote various dimensional constants. We
may assume $x_{0}=0$ by shifting the coordinates properly. By our assumption
$\theta>\pi/4$, we have $B_{1}(-2e_{n})\subset\\{x_{n}<f(x^{\prime})\\}$,
where $f(x^{\prime}):=-\cot\theta|x^{\prime}|$ as above. Since $-\Delta
p(\cdot,t_{0})\geq\bar{n}$ in $B_{1}(-2e_{n})$, it follows that
$p(x,t_{0})\geq C\bar{n}$ in $B_{1/2}(-2e_{n})$.

Next consider a harmonic function $q$ in the domain
$\\{x_{n}<f(x^{\prime})\\}\cap(B_{2}(0)\backslash B_{1/2}(-2e_{n}))$, with
boundary data

$q=0\hbox{ on }\\{x_{n}=f(x^{\prime})\\}\cup\partial B_{2}(0)\hbox{ and
}q=C\bar{n}\hbox{ on }\partial B_{1/2}(-2e_{n}).$

Since $p(\cdot,t_{0})$ is superharmonic, it follows that $p(\cdot,t_{0})\geq
q$ in $B_{1/2}(0)\cap\\{x_{n}<f(x^{\prime})\\}$.

To obtain explicit lower bound for $q$ near $\\{x_{n}=f(x^{\prime})\\}$, let
us compare $q$ with $h$ given in (5.1). By Dahlberg’s lemma, $q\geq C\bar{n}h$
in $\\{x_{n}<f(x^{\prime})\\}\cap B_{1}(x_{0})$. In particular, it follows
that

(5.2) $p(x,t_{0})\geq C\bar{n}(f(x^{\prime})-x_{n})_{+}^{2-\alpha}\hbox{ in
}\\{x_{n}<f(x^{\prime})\\}\cap B_{1}(x_{0}).$

Since $\\{p(\cdot,t)>0\\}$ increases in time, by the same logic and the
property of the $n$-equation $n_{t}=-n\rho$, we have

(5.3) $p(x,t)\geq
C\bar{n}e^{-(t-t_{0})}(f(x^{\prime})-x_{n})_{+}^{2-\alpha}\hbox{ in
}B_{1}(x_{0})\times[t_{0},\infty).$

Now we will use Lemma 5.8 to estimate the time it takes for the positive set
of $p$ to fully cover $B_{r}(x_{0})$. Let $x_{r}:=x_{0}-re_{n}$ and
$D:=\\{x:|x-x_{r}|>\frac{r}{4}\\}$. Let $\bar{p}$ be a barrier solving

$\left\\{\begin{array}[]{lll}-\Delta\bar{p}=0&\hbox{ in
}&\frac{r}{4}<|x-x_{r}|<r(t),\\\ \\\
\bar{p}=C\bar{n}e^{-(t-t_{0})}r^{2-\alpha}&\hbox{ on
}&\\{|x-x_{r}|=\frac{r}{4}\\},\\\ \\\ \bar{p}=0&\hbox{ on
}&\\{|x-x_{r}|=r(t)\\}.\end{array}\right.$

If we set $r(t)$ so that $r^{\prime}(t)\leq|D\bar{p}|$ on
$\\{|x-x_{r}|=r(t)\\}$, then $\bar{p}$ will solve

$\bar{\rho}_{t}-\Delta\bar{p}\leq 0\hbox{ in }D\times[t_{0},\infty),$

where $\bar{\rho}=\chi_{\\{\bar{p}>0\\}}$. If in addition
$r(t_{0})=\frac{r}{2}$ so that
$\\{\bar{p}(\cdot,t_{0})>0\\}\subset\\{p(\cdot,t_{0})>0\\}$, we will apply
Lemma 5.8 and (5.3) to conclude that $\bar{p}\leq p$ for $t>t_{0}$, which will
then yield a lower bound for the time it takes for the positive set of $p$ to
include $B_{r}(x_{0})$.

Now we choose a specific $r(t)$. Observe that $|D\bar{p}|\geq
C\bar{n}e^{-(t-t_{0})}r^{1-\alpha}$ on $\\{|x-x_{r}|=r(t)\\}$ as long as
$r(t)\in[\frac{1}{2}r,2r]$. Hence we can set

$r(t):=\frac{r}{4}+C\bar{n}\big{(}1-e^{-(t-t_{0})}\big{)}r^{1-\alpha}.$

Now we conclude, since $r(t)=2r$ when $\bar{n}(1-e^{-(t-t_{0})})\sim
r^{\alpha}$, i.e., when $t-t_{0}\sim\frac{r^{\alpha}}{\bar{n}}$. ∎

######  Corollary 5.13.

Suppose that, given a point $x_{0}\in\partial\\{w(\cdot,t_{1})>0\\}$ with
$t_{1}>0$, the set $\partial\\{w(\cdot,t)>0\\}$ is a Lipschitz graph in
$B_{2}(x_{0})$ for all $0<t\leq t_{1}$ with respect to a fixed direction.
Further suppose that the Lipschitz constants of the graphs are all smaller
than a dimensional constant. Then the followings hold:

  1. (a)

$\eta(\cdot,t_{1})\in C^{\alpha}(B_{1}(x_{0}))$ ;

  2. (b)

$\partial\\{w(\cdot,t_{1})>0\\}$ is $C^{1,\alpha}$ in $B_{1}(x_{0})$.

Here the $C^{\alpha}$ and $C^{1,\alpha}$ norms only depend on
$\bar{n}:=\min_{B_{2}(x_{0})}n(\cdot,t_{1})$.

###### Proof.

Define $T_{x}:=\inf\\{t:w(x,t)>0\\}$ for all $x\in\mathbb{R}^{d}$. Then for
$x,y\in B_{1}(x_{0})\cap\\{w(\cdot,t_{1})>0\\}$, Proposition 5.12 implies
that, with $r$ being smaller than a universal constant, $T_{y}\leq
T_{x}+\frac{Cr^{\alpha}}{\bar{n}}$ as long as $|x-y|<r$. Hence we have

$|T_{x}-T_{y}|\leq\frac{C|x-y|^{\alpha}}{\bar{n}}\hbox{ for all }x,y\in
B_{1}(x_{0}).$

Note that $\eta(x,t)=n_{0}(x)(1-\exp^{-(t-T_{x})_{+}})$, and $n_{0}$ is Hölder
continuous by our assumption. Thus we conclude (a) by deriving that

$|\eta(x,t_{1})-\eta(y,t_{1})|\leq\frac{C|x-y|^{\alpha}}{\bar{n}}\hbox{ for
all }x,y\in B_{1}(x_{0}).$

(b) then follows from [Bla01, Theorem 7.2]. ∎

######  Remark 5.14.

Once we obtain $C^{1,\alpha}$-regularity of the free boundary
$\partial\\{w(\cdot,t)>0\\}$, the interior cone angle $\theta$ given in
Proposition 5.12 can be chosen as close to $\frac{\pi}{2}$ as desired, in a
smaller scale. We can thus improve the regularity of $\eta$ to
$C^{1-\epsilon}$ for any $\epsilon>0$, which in turn improves the free
boundary regularity to $C^{1,1-\epsilon}$ for any $\epsilon>0$.

Combining the above results with Corollary 5.5, we arrive at the following
conclusion.

######  Corollary 5.15.

Suppose $n_{0}$ is $C^{1}$, and its super-level sets satisfy $r$-reflection
for some $r>0$. Let $\Omega_{0}$ be an open bounded set in $\mathbb{R}^{d}$
contained in $B_{r}(0)$. Let $T(R)$ be defined as in Corollary 5.5. Then the
followings hold for any $0<\alpha<1$ and any $R\geq Cr$ with $C>1$ being a
universal constant:

  1. (a)

If $n_{0}(x)\geq 1$ on $\mathbb{R}^{d}$, then
$\partial\\{\rho(\cdot,t)=1\\}=\partial\\{w(\cdot,t)>0\\}$ is uniformly
$C^{1,\alpha}$ in a unit neighborhood for any finite time range within
$[T(R),\infty)$.

  2. (b)

If $n_{0}(x)<1$ on $\mathbb{R}^{d}$, the same holds in for any finite time
range within $[T(R),\infty)$ if $R$ additionally satisfies
$B_{R}(0)\subset\\{w_{\infty}>0\\}$. Here $w_{\infty}$ is defined in Theorem
4.10.

Both Corollary 5.13 and Corollary 5.15 only apply to finite time ranges. It is
natural to ask what can be said about uniform regularity of the free boundary
up to $t=+\infty$. This is a nontrivial question due to the possible decay of
$\bar{n}$ (defined in Corollary 5.13) as $t$ tends to infinity. Note that,
roughly speaking, at a boundary point $x_{0}$, $\bar{n}\geq n_{0}e^{-t}$ if
the free boundary does not move significantly. When $n_{0}<1$, this bound is
close to optimal since the tumor patch converges to a bounded set. When
$n_{0}\geq 1$, however, it is possible to improve this bound by comparison
with radial barriers. We will discuss this in Section 6.2, under the
additional assumption that $n_{0}$ is constant.

##  6\. The Constant $n_{0}$ Case

In this section, we further look into the case when $n_{0}$ is a positive
constant, while still inheriting the assumptions that $b=D=0$ and $\rho_{0}$
is a patch with compact support. What is special and surprising in this case
is that the dynamics of the system (P) can be fully characterized by simpler
parameter-free and nutrient-free model problems, which produce the so-called
_master dynamics_. As an application of that, we will address the question of
uniform regularity of the free boundary $\partial\Omega_{t}$.

Recall that for patch solutions, we have $\eta\rho=\eta$. Then $\eta$ solves
(c.f. (4.1))

$\partial_{t}\eta+\eta=n_{0}\rho,\quad\eta(x,0)=0.$

This gives

(6.1) $\eta(x,t)=n_{0}\int_{0}^{t}e^{-(t-\tau)}\rho(x,\tau)\,d\tau,$

and therefore,

(6.2) $n(x,t)=n_{0}-n_{0}\int_{0}^{t}e^{-(t-\tau)}\rho(x,\tau)\,d\tau.$

This holds even for non-constant $n_{0}$. As a result, in the following, we
will focus on the $\rho$-evolution in (P). Let us also recall the elliptic
formulation of the $\rho$-equation in $\mathbb{R}^{d}$ that is derived from
(P)

(6.3) $\rho-\Delta w=\rho_{0}+\eta,\quad w(1-\rho)=0,\quad\rho\in[0,1],\quad
w\geq 0.$

Here $\rho_{0}$ is a patch with compact support.

###  6.1. Two master dynamics

We first prove a growth law of total mass of $\rho$ and a generalization of
it. The latter will be the key of proving the master dynamics.

######  Lemma 6.1.

Define (c.f. Theorem 4.7)

(6.4)
$m(t)=1+n_{0}\int_{0}^{t}e^{(n_{0}-1)\tau}\,d\tau=\begin{cases}\frac{n_{0}e^{(n_{0}-1)t}-1}{n_{0}-1},&\mbox{if
}n_{0}\neq 1,\\\ t+1,&\mbox{otherwise}.\end{cases}$

  1. (a)

For any $t\geq 0$,

$\int_{\mathbb{R}^{d}}\rho(x,t)\,dx=m(t)\int_{\mathbb{R}^{d}}\rho_{0}(x)\,dx.$

  2. (b)

For any arbitrary smooth function $g=g(x)$ in $\mathbb{R}^{d}$, we have for
any $t\geq 0$,

(6.5)
$\begin{split}\int_{\mathbb{R}^{d}}\rho(x,t)g(x)\,dx=&\;m(t)\int_{\mathbb{R}^{d}}\rho_{0}(x)g(x)\,dx\\\
&\;+\int_{\mathbb{R}^{d}}w(x,t)\Delta
g(x)\,dx+\int_{0}^{t}m^{\prime}(t-\tau)\left[\int_{\mathbb{R}^{d}}w(x,\tau)\Delta
g(x)\,dx\right]d\tau.\end{split}$

In particular, if $g$ is harmonic in an open neighborhood of the support of
$\rho(x,t)$, then

(6.6)
$\int_{\mathbb{R}^{d}}\rho(x,t)g(x)\,dx=m(t)\int_{\mathbb{R}^{d}}\rho_{0}(x)g(x)\,dx.$

Combined with (a), this implies that when $\rho$ is a patch solution, the
average of any harmonic $g$ on $\\{\rho(\cdot,t)=1\\}$ is time-invariant.

###### Proof.

Take an arbitrary smooth function $g=g(x)$ in $\mathbb{R}^{d}$. We integrate
(6.3) in $\mathbb{R}^{d}$ against $g(x)$ and use (6.1) to find that

$\int_{\mathbb{R}^{d}}\rho(x,t)g(x)\,dx=\int_{\mathbb{R}^{d}}\rho_{0}g+w(x,t)\Delta
g(x)\,dx+n_{0}\int_{0}^{t}e^{-(t-\tau)}\int_{\mathbb{R}^{d}}\rho(x,\tau)g(x)\,dx\,d\tau.$

Solving this integral equation, we obtain

$\begin{split}\int_{\mathbb{R}^{d}}\rho(x,t)g(x)\,dx=&\;\int_{\mathbb{R}^{d}}\rho_{0}g+w(x,t)\Delta
g(x)\,dx\\\
&\;+n_{0}\int_{0}^{t}e^{(n_{0}-1)(t-\tau)}\int_{\mathbb{R}^{d}}\rho_{0}g+w(x,\tau)\Delta
g(x)\,dx\,d\tau\\\
=&\;\left(1+n_{0}\int_{0}^{t}e^{(n_{0}-1)(t-\tau)}\,d\tau\right)\int_{\mathbb{R}^{d}}\rho_{0}g\,dx\\\
&\;+\int_{\mathbb{R}^{d}}w(x,t)\Delta
g(x)\,dx+n_{0}\int_{0}^{t}e^{(n_{0}-1)(t-\tau)}\left[\int_{\mathbb{R}^{d}}w(x,\tau)\Delta
g(x)\,dx\right]d\tau.\end{split}$

Taking $g(x)=1$ yields (6.4). That in turn implies (6.5). ∎

######  Remark 6.2.

Taking $g(x)=x_{i}$ $(i=1,\cdots,d)$ in (6.6), we find the center of mass of
$\rho(x,t)$ is time-invariant.

######  Remark 6.3.

Suppose $\rho_{0}=\chi_{B_{r_{0}}(0)}$. Then thanks to (6.4), Remark 6.2, and
the uniqueness (c.f. Proposition 3.13 and Theorem 4.7),

$\rho(x,t)=\chi_{B_{r(t)}(0)},\mbox{ where }r(t)=m(t)^{\frac{1}{d}}r_{0}.$

Here $m(t)$ is given in (6.4). In this paper, such radial solutions are
repeatedly used as barriers for comparison.

######  Lemma 6.4.

Let $\Gamma$ be the fundamental solution of $-\Delta$ in $\mathbb{R}^{d}$,
i.e., $-\Delta\Gamma=\delta_{x=0}$. For any $x\in\\{w(\cdot,t)=0\\}$,

(6.7)
$\big{(}\Gamma*\rho(\cdot,t)\big{)}(x)=m(t)\big{(}\Gamma*\rho_{0}\big{)}(x).$

###### Proof.

We introduce a smooth mollifier $\varphi\in C_{0}^{\infty}(\mathbb{R}^{d})$
such that $\varphi\geq 0$, $\varphi$ is radially symmetric, and
$\int_{\mathbb{R}^{d}}\varphi=1$. With $\epsilon>0$, define
$\varphi_{\epsilon}(x):=\epsilon^{-d}\varphi(\frac{x}{\epsilon})$. Let
$\Gamma_{\epsilon}:=\Gamma*\varphi_{\epsilon}$, which is clearly smooth in
$\mathbb{R}^{d}$ and which satisfies
$-\Delta\Gamma_{\epsilon}=\varphi_{\epsilon}$. Applying Lemma 6.1 with
$g(\cdot)=\Gamma_{\epsilon}(x-\cdot)$, we find that, for any
$x\in\mathbb{R}^{d}$,

$\begin{split}\big{(}\Gamma_{\epsilon}*\rho(\cdot,t)\big{)}(x)=&\;\int_{\mathbb{R}^{d}}\rho(y,t)\Gamma_{\epsilon}(x-y)\,dy\\\
=&\;m(t)\int_{\mathbb{R}^{d}}\rho_{0}(y)\Gamma_{\epsilon}(x-y)\,dy\\\
&\;+\int_{\mathbb{R}^{d}}w(y,t)\Delta_{y}\Gamma_{\epsilon}(x-y)\,dy+\int_{0}^{t}m^{\prime}(t-\tau)\left[\int_{\mathbb{R}^{d}}w(y,\tau)\Delta_{y}\Gamma_{\epsilon}(x-y)\,dy\right]d\tau\\\
=&\;m(t)\big{(}\Gamma_{\epsilon}*\rho_{0}\big{)}(x)-\varphi_{\epsilon}*\left[w(\cdot,t)+\int_{0}^{t}m^{\prime}(t-\tau)w(\cdot,\tau)\,d\tau\right].\end{split}$

Then we take $\epsilon\to 0$. Since $w(x,t)$ is non-decreasing in time (c.f.
Proposition 4.1 and the fact $\eta$ is non-decreasing in time), the term in
the bracket in the last line can be dominated by $m(t)w(\cdot,t)$. Hence, for
any $x\in\mathbb{R}^{d}$ satisfying

(6.8) $\lim_{\epsilon\to 0}\big{(}\varphi_{\epsilon}*w(\cdot,t)\big{)}(x)=0,$

we can show (6.7) by using the spatial continuity of $\Gamma*\rho$ and
$\Gamma*\rho_{0}$. Since $w(\cdot,t)$ is continuous, the condition (6.8) holds
in the set $\\{w(\cdot,t)=0\\}$. This completes the proof. ∎

######  Lemma 6.5.

Let $m(t)$ be given by (6.4). Consider the elliptic equation in
$\mathbb{R}^{d}$

(6.9)
$\tilde{\rho}-\Delta\tilde{w}=m(t)\rho_{0},\quad\tilde{w}(1-\tilde{\rho})=0,\quad\tilde{\rho}\in[0,1],\quad\tilde{w}\geq
0.$

Then $\tilde{\rho}=\rho$ for all $t\geq 0$ and almost everywhere in space.

###### Proof.

By (6.1) and Lemma 6.4, For all $x\in\\{w(\cdot,t)=0\\}$,

$\big{(}\Gamma*\eta(\cdot,t)\big{)}(x)=n_{0}\int_{0}^{t}e^{-(t-\tau)}m(\tau)\,d\tau\cdot(\Gamma*\rho_{0})(x)=n_{0}\left(m(t)-e^{(n_{0}-1)t}\right)(\Gamma*\rho_{0})(x).$

Here we calculated by integration by parts that

$\begin{split}\int_{0}^{t}e^{-(t-\tau)}m(\tau)\,d\tau=&\;e^{-t}\left(e^{t}m(t)-1-\int_{0}^{t}e^{\tau}m^{\prime}(\tau)\,d\tau\right)\\\
=&\;m(t)-e^{-t}-e^{-t}\int_{0}^{t}e^{\tau}\cdot
n_{0}e^{(n_{0}-1)\tau}\,d\tau\\\ =&\;m(t)-e^{(n_{0}-1)t}.\end{split}$

In the view of right-hand sides of (6.3) and (6.9), let

$\Phi:=\Gamma*\left(\rho_{0}+\eta(\cdot,t)-m(t)\rho_{0}\right).$

One can check that, for all $x\in\\{w(\cdot,t)=0\\}$,

$\Phi(x)=\left[(n_{0}-1)m(t)-n_{0}e^{(n_{0}-1)t}+1\right](\Gamma*\rho_{0})(x)=0.$

Now we claim that

$\tilde{\rho}=\rho,\quad\tilde{w}=w-\Phi,$

gives a solution of (6.9). One only has to verify that
$\tilde{w}(1-\tilde{\rho})=0$. From what has been proved, whenever $w=0$, we
have $\Phi=0$ and thus $\tilde{w}=0$. On the other hand, $\tilde{\rho}=1$
whenever $\rho=1$. Therefore, $\tilde{w}(1-\tilde{\rho})=0$ holds whenever
$w(1-\rho)=0$. This justifies the claim.

Since the solution of (6.9) is unique (c.f. Proposition 4.1), we conclude that
$\tilde{\rho}=\rho$. ∎

This leads to the following equivalent characterization of
$\\{\rho(\cdot,t)\\}_{t}$.

######  Proposition 6.6 (Master dynamics I).

Let $\rho_{*}=\rho_{*}(x,t)$ and $p_{*}=p_{*}(x,t)$ be a weak solution (in the
sense of Definition 3.9) of

(6.10) $\partial_{t}\rho_{*}-\nabla\cdot(\rho_{*}\nabla
p_{*})=\rho_{0},\quad\rho_{*}\leq 1,\quad p_{*}\in
P_{\infty}(\rho_{*}),\quad\rho_{*}|_{t=0}=\rho_{0}.$

Then for any given $n_{0}>0$, $\\{\rho(x,t)\\}_{t}$ defined by (6.1) and (6.3)
(or equivalently, by (P)) satisfies

$\rho(x,t)=\rho_{*}(x,m(t)-1)\mbox{ for all }t\geq 0,$

where $m(t)$ is defined by (6.4). $\\{\rho_{*}(x,t)\\}_{t\geq 0}$ is thus
called the first master dynamics.

In particular, when $n_{0}=1$, $\rho(x,t)=\rho_{*}(x,t)$.

We stress that this is a highly non-trivial property of the model under the
given assumptions. It cannot be obtained by a nonlinear change of the time
variable.

Following the spirit of Lemmas 6.1, 6.4, and 6.5, we may derive a second
equivalent characterizations of $\\{\rho(x,t)\\}_{t}$. It will be particularly
helpful for understanding long-time behavior of $\rho$ when $n_{0}\geq 1$
because it comes with a suitable spatial re-scaling.

######  Proposition 6.7 (Master dynamics II).

Let $\rho_{\dagger}=\rho_{\dagger}(x,t)$ and $p_{\dagger}=p_{\dagger}(x,t)$ be
a weak solution of

(6.11)
$\partial_{t}\rho_{\dagger}-\nabla\cdot\Big{[}\rho_{\dagger}\nabla\big{(}p_{\dagger}+V(x)\big{)}\Big{]}=0,\quad\rho_{\dagger}\leq
1,\quad p_{\dagger}\in
P_{\infty}(\rho_{\dagger}),\quad\rho_{\dagger}|_{t=0}=\rho_{0},$

where $V(x)=\frac{|x|^{2}}{2d}$. Then for any given $n_{0}>0$,
$\\{\rho(x,t)\\}_{t}$ defined by (6.1) and (6.3) (or equivalently, by (P))
satisfies

$\rho(x,t)=\rho_{\dagger}\left(m(t)^{-\frac{1}{d}}x,\,\ln m(t)\right)\mbox{
for all }t\geq 0,$

where $m(t)$ is defined by (6.4). Thus, $\\{\rho_{\dagger}(x,t)\\}_{t\geq 0}$
is called the second master dynamics.

###### Proof.

Following a similar argument as in Lemmas 6.1, 6.4, and 6.5, we can prove
that, if one defines

$\tilde{\tilde{\rho}}-\Delta\tilde{\tilde{w}}=\rho_{0}+\int_{0}^{t}\frac{m^{\prime}(\tau)}{m(\tau)}\tilde{\tilde{\rho}}(\cdot,\tau)\,d\tau,\quad\tilde{\tilde{w}}\big{(}1-\tilde{\tilde{\rho}}\big{)}=0,$

then $\tilde{\tilde{\rho}}=\rho$ for all $t\geq 0$. Indeed, it suffices to
show that, on $\\{\tilde{\tilde{w}}(\cdot,t)=0\\}$,

$\big{(}\Gamma*\tilde{\tilde{\rho}}(\cdot,t)\big{)}(x)=m(t)\big{(}\Gamma*\rho_{0}\big{)}(x),$

and thus

$\Gamma*\left(\rho_{0}+\int_{0}^{t}\frac{m^{\prime}(\tau)}{m(\tau)}\tilde{\tilde{\rho}}(\cdot,\tau)\,d\tau-m(t)\rho_{0}\right)\equiv
0$

on the same set. We skip the details.

Hence, $\rho$ satisfies

(6.12)
$\partial_{t}\rho-\nabla\cdot\big{(}\rho\nabla\tilde{\tilde{p}}\big{)}=\frac{m^{\prime}(t)}{m(t)}\rho,\quad\rho\leq
1,\quad\tilde{\tilde{p}}\in P_{\infty}(\rho),\quad\rho|_{t=0}=\rho_{0}.$

It is then straightforward to verify that, if $\rho_{\dagger}$ and
$p_{\dagger}$ solve (6.11), then

$\displaystyle\rho(x,t):=\rho_{\dagger}\left(m(t)^{-\frac{1}{d}}x,\,\ln
m(t)\right),$
$\displaystyle\tilde{\tilde{p}}(x,t):=\frac{m^{\prime}(t)}{m(t)}\cdot
m(t)^{\frac{2}{d}}p_{\dagger}\left(m(t)^{-\frac{1}{d}}x,\,\ln m(t)\right)$

satisfy (6.12). ∎

Using Proposition 6.7, one can readily characterize the long-time behavior of
the $\rho$-patch when $n_{0}>0$ is constant. Indeed, when $n_{0}\geq 1$, the
long-time dynamics of $\rho$ in (6.1) and (6.3) corresponds to infinite-time
asymptotics of $\rho_{\dagger}$ in (6.11), which has been well-studied the
literature, see e.g. [AKY14, Theorem 5.6]. While for $n_{0}\in(0,1)$, we knew
from Theorem 4.10 that $\rho$ converges to a compactly supported
$\rho_{\infty}$ as $t\to+\infty$. On the other hand, under the rescaling of
Proposition 6.7, such dynamics of $\rho$ actually corresponds to an excerpt of
the master $\rho_{\dagger}$-dynamics up to a finite time. Therefore, we may
characterize the long-time behavior of $\rho$ in the model (6.1) and (6.3)
with any value of $n_{0}$ in the following unified way.

######  Proposition 6.8.

Suppose $n_{0}>0$ is constant in $\mathbb{R}^{d}$. Let $m(t)$ be defined in
(6.4), and denote $\beta(t):=m(t)^{\frac{1}{d}}$. Assume $\Omega_{0}$ to be a
bounded open set, such that $B_{r_{1}}(0)\subset\Omega_{0}\subset
B_{r_{2}}(0)$ for some $r_{1},r_{2}>0$. Let $r_{\infty}>0$ be defined such
that $|B_{r_{\infty}}(0)|=|\Omega_{0}|$.

Let $\rho(x,t)$ solve (6.1) and (6.3) with $\rho_{0}=\chi_{\Omega_{0}}$. Then
there exists a constant $C>0$ only depending on $r_{1}$ and $r_{2}$, but not
on $n_{0}$, such that

$\beta(t)W_{2}\Big{(}\rho\big{(}\beta(t)x,t\big{)},\chi_{B_{r_{\infty}}(0)}(x)\Big{)}\leq
C$

for all $t\geq 0$. Here $W_{2}$ denotes the 2-Wasserstein distance.

###### Proof.

Let us recall that, by virtue of Lemma 6.1, $\rho(\beta(t)x,t)$ has the same
total mass as $\chi_{B_{r_{\infty}}(0)}(x)$.

First we assume $n_{0}\geq 1$.

By comparison with the radial barriers, $\rho(\beta(t)x,t)$ and
$\chi_{B_{r_{\infty}}(0)}(x)$ are both supported in $B_{r_{2}}(0)$ for all
time. Let $T>0$ satisfy $\beta(T)=\frac{2r_{2}}{r_{1}}$. Then it is obvious
that the above inequality holds on $[0,T]$, with $C>0$ only depending on
$r_{1}$ and $r_{2}$.

Next we consider the case $t\geq T$. Assume $\rho(x,t)=\chi_{\Omega_{t}}(x)$.
By Corollary 5.5(b) and comparison with the radial barriers (see Remark 6.3),
for all $t\geq T$, $\Omega_{t}$ contains $B_{2r_{2}}(0)$ and thus it has
Lipschitz boundary. By Proposition 6.7,

$\rho_{\dagger}(x,\ln
m(T))=\rho\big{(}\beta(T)x,T\big{)}=\chi_{\beta(T)^{-1}\Omega_{T}}(x).$

Note that the rescaled set $\beta(T)^{-1}\Omega_{T}$ has Lipschitz boundary.
We then solve (6.11) starting from $t=\ln m(T)$ with “initial data”
$\rho_{\dagger}(x,\ln m(T))$. By [AKY14, Theorem 5.6 and its proof], for all
$t^{\prime}\geq 0$,

$W_{2}\Big{(}\rho_{\dagger}\big{(}x,\ln
m(T)+t^{\prime}\big{)},\chi_{B_{r_{\infty}}(0)}(x)\Big{)}\leq
e^{-\frac{t^{\prime}}{d}}W_{2}\Big{(}\rho_{\dagger}\big{(}x,\ln
m(T)\big{)},\chi_{B_{r_{\infty}}(0)}(x)\Big{)}.$

Then the desired result follows from suitable change of variables.

By Proposition 6.7,

(6.13) $\rho_{\dagger}(x,\ln m(t))=\rho(\beta(t)x,t)\mbox{ holds for all
}n_{0}>0.$

Hence, the above argument essentially proves that, without assuming
$\Omega_{0}$ has Lipschitz boundary (c.f. [AKY14]), for all $\tau\geq 0$,

$e^{\frac{\tau}{d}}W_{2}\big{(}\rho_{\dagger}(x,\tau),\chi_{B_{r_{\infty}}(0)}(x)\big{)}\leq
C.$

Hence, using (6.13) again, the case $n_{0}\in(0,1)$ is proved immediately. ∎

###  6.2. Uniform free boundary regularity

In this section, we are going to prove uniform free boundary regularity up to
$t=+\infty$ under the assumption that $n_{0}$ is constant in $\mathbb{R}^{d}$.

Thanks to the master dynamics, the case $n_{0}<1$ is trivial. This is because,
up to a re-scaling in time (see Proposition 6.6), its $\rho$-evolution in a
time range of the form $(t,+\infty)$ corresponds to the $\rho$-evolution in a
finite time range with a different $n_{0}\geq 1$. The latter has already been
characterized in Corollary 5.15(a).

Also by the master dynamics, it suffices to study regularity of
$\partial\Omega_{t}$ of the solution corresponding to _one_ arbitrary
$n_{0}\geq 1$. When $n_{0}>1$, one can apply comparison principle and radial
barriers constructed in Remark 6.3 to show that $\Omega_{t}$ expands
exponentially fast. It is thus reasonable to expect uniform regularity of the
free boundary after a suitable re-scaling. Indeed, the key lies in the uniform
Lipschitz estimate for $\partial\\{w(\cdot,t)>0\\}$ in Corollary 5.5(b).

Let us point out that, while the master dynamics in rescaled variable (6.11)
corresponds to a Hele-Shaw flow, the presence of the drift prevents us from
directly applying existing regularity results (e.g. [CJK07]) to our problem.

######  Theorem 6.9.

Fix $n_{0}>1$. Let $\Omega_{0}$, $r_{1}$, $r_{2}$, and $\beta(t)$ be given as
in Proposition 6.8. Let $\rho_{0}=\chi_{\Omega_{0}}$. Then there is
$\alpha\in(0,1)$ and $T>0$ depending on $r_{1}$, $r_{2}$, $d$, and $n_{0}$,
such that the followings hold for all $t\geq T$.

  1. (a)

The rescaled set $\tilde{\Omega}_{t}:=\beta(t)^{-1}\Omega_{t}$ has uniformly
$C^{1,\alpha}$-boundary;

  2. (b)

The rescaled nutrient variable $\tilde{n}(x,t):=n(\beta(t)x,t)$ is uniformly
bounded in $C^{\alpha}(\\{|x|\geq 2\beta^{-1}(t)r_{2}^{2}/r_{1}\\})$.

###### Proof.

From comparison with radial barriers (see Theorem 4.7 and Remark 6.3), we find
that

(6.14) $B_{\beta(t)r_{1}}(0)\subset\Omega_{t}=\\{\rho(\cdot,t)=1\\}\subset
B_{\beta(t)r_{2}}(0)$

up to measure-zero set. Also, for some $C>0$ depending on $n_{0}$ and $r_{1}$,

$p(x,t)\geq C\beta^{2}(t)\hbox{ if }|x|\leq\frac{r_{1}\beta(t)}{2}.$

Hence, the re-scaled pressure variable
$\tilde{p}(x,t):=\beta^{-2}(t)p(\beta(t)x,t)$ satisfies

$\tilde{p}(x,t)\geq C\hbox{ if }|x|\leq\frac{r_{1}}{2}.$

On the other hand, let $T(\cdot)$ be introduced in Corollary 5.5. From
Corollary 5.5, for $t\geq T(2r_{2})$, we know that
$\partial\\{p(\cdot,t)>0\\}$ is a Lipschitz graph with respect to the radial
direction, with the Lipschitz constant less than
$O((\beta(t)r_{1})^{-1}r_{2})$. Therefore,
$\tilde{\Omega}_{t}=\\{\tilde{p}(\cdot,t)>0\\}\subset B_{r_{2}}(0)$, with
$\partial\tilde{\Omega}_{t}$ being uniformly Lipschitz with respect to the
radial direction. When $t$ is suitably large, depending on $r_{1}$, $r_{2}$,
and $n_{0}$, we have the Lipschitz constant to be small enough for applying a
similar argument as in Proposition 5.12.

Since $\tilde{p}$ is superharmonic in its positive set, arguing with
Dahlberg’s Lemma as in the proof of Proposition 5.12, we conclude that there
is a constant $\alpha\in(0,1)$ that is independent of the time such that, for
given $\tilde{x}_{*}\in\partial\tilde{\Omega}_{t}$ and for $t\geq T(2r_{2})$
we have

$\tilde{p}(x,t)\geq Cd(x,\tilde{\Omega}_{t}^{c})^{2-\alpha}\hbox{ in
}B_{1}(\tilde{x}_{*}).$

In the original coordinate, this corresponds to

$p(x,t)\geq C\beta(t)^{\alpha}d(x,\Omega_{t}^{c})^{2-\alpha}$

for any $x\in B_{1}(\beta(t)\tilde{x}_{*})$ and $t\geq T(2r_{2})$. Thus the
barrier argument as in the proof of Proposition 5.12, with $\bar{n}$ replaced
by $\beta(t)^{\alpha}$, yields that, for any $x_{*}\in\partial\Omega_{t}$ with
$t\geq T(2r_{2})$,

(6.15) $B_{r}(x_{*})\subset\Omega_{t+C(\frac{r}{\beta(t)})^{\alpha}}.$

for sufficiently small $r$. This further implies

(6.16) $|T_{x}-T_{y}|\leq
C\left(\frac{|x-y|}{\min\\{\beta(T_{x}),\beta(T_{y})\\}}\right)^{\alpha}\hbox{
for any }x,y\mbox{ with }T_{x},T_{y}\geq T(2r_{2}).$

To justify this, we assume $T_{x}<T_{y}$ without loss of generality. We first
consider the case where $|x-y|$ is large. By (6.14), if $T_{*}>0$ satisfies
that

$\beta(T_{x}+T_{*})r_{1}\geq|y|,$

then $|T_{x}-T_{y}|\leq T_{*}$. Also by (6.14),

$|y|\leq|x-y|+|x|\leq|x-y|+\beta(T_{x})r_{2}.$

Hence, we let $T_{*}$ satisfy

$\beta(T_{x}+T_{*})r_{1}=|x-y|+\beta(T_{x})r_{2},$

which implies (c.f. (6.4))

$e^{(n_{0}-1)T_{*}}\leq\frac{m(T_{x}+T_{*})}{m(T_{x})}=\left(\frac{|x-y|}{\beta(T_{x})r_{1}}+\frac{r_{2}}{r_{1}}\right)^{d},$

Therefore,

$|T_{x}-T_{y}|\leq\frac{d}{n_{0}-1}\ln\left(\frac{|x-y|}{\beta(T_{x})r_{1}}+\frac{r_{2}}{r_{1}}\right),$

which implies (6.16) whenever $|x-y|/\beta(T_{x})$ is sufficiently large.
Otherwise, if $|x-y|/\beta(T_{x})\leq C$ where $C$ depends on $r_{1}$,
$r_{2}$, $d$, $n_{0}$, and $\alpha$, we may apply (6.15) to obtain (6.16).

Under the assumption $T_{x}<T_{y}$, we have that

$\big{|}\eta(x,t)-\eta(y,t)\big{|}\leq
n_{0}\left(1-e^{-|T_{x}-T_{y}|}\right)e^{-(t-T_{y})_{+}}\leq
n_{0}\min\big{\\{}1,|T_{x}-T_{y}|\big{\\}}e^{-(t-T_{y})_{+}}.$

Hence, thanks to (6.16), for any $x,y$ such that $T_{x},T_{y}\geq T(2r_{2})$,

$\frac{|\eta(x,t)-\eta(y,t)|}{|x-y|^{\alpha}}\leq
C\min\big{\\{}|x-y|^{-\alpha},\beta(T_{x})^{-\alpha}\big{\\}}e^{-(t-T_{y})_{+}}.$

Since we defined $\tilde{\eta}(x,t):=\eta(\beta(t)x,t)$,

(6.17)
$\frac{|\tilde{\eta}(\beta(t)^{-1}x,t)-\eta(\beta(t)^{-1}y,t)|}{\beta(t)^{-\alpha}|x-y|^{\alpha}}\leq
C\beta(t)^{\alpha}\max\big{\\{}|x-y|,\beta(T_{x})\big{\\}}^{-\alpha}e^{-(t-T_{y})_{+}}.$

We claim that

$\beta(T_{y})\leq C\big{(}|x-y|+\beta(T_{x})\big{)}\leq
C\max\big{\\{}|x-y|,\beta(T_{x})\big{\\}},$

where $C$ may depend on $r_{1}$, $r_{2}$, $d$, and $n_{0}$. Indeed, by the
estimate for $T_{*}$ derived above

$\frac{\beta(T_{y})}{\beta(T_{x})}\leq
Ce^{\frac{(n_{0}-1)}{d}|T_{y}-T_{x}|}\leq Ce^{\frac{(n_{0}-1)}{d}T_{*}}\leq
C\left(\frac{|x-y|}{\beta(T_{x})}+1\right).$

Hence, when $t\leq T_{y}$, the right-hand side of (6.17) is bounded by a
universal constant that only depends on $r_{1}$, $r_{2}$, $d$, $n_{0}$, and
$\alpha$. We may further assume $\alpha$ to be suitably small so that the
right-hand side of (6.17) is uniformly bounded for $t\geq T_{y}$. Now noticing
that $|x|\geq 2r_{2}^{2}/r_{1}$ guarantees $T_{x}\geq T(2r_{2})$ (c.f.
(6.14)), we can conclude (b), i.e., $\tilde{\eta}(x,t)$ has uniform-in-time
Hölder regularity for $t\geq T(2r_{2})$ for $|x|\geq
2\beta^{-1}(t)r_{2}^{2}/r_{1}$.

Lastly, observe that $\tilde{w}(x,t):=w(\beta(t)x,t)$ solves the obstacle
problem

$\Delta\tilde{w}=\tilde{f}\chi_{\\{\tilde{w}(\cdot,t)>0\\}},\mbox{ where
}\tilde{f}=1-\tilde{\rho}_{0}-\tilde{\eta}.$

Then using the regularity of $\tilde{\eta}$, and Theorem 7.2 of [Bla01], we
can conclude (a). ∎

Before ending this section, let us briefly discuss the uniform boundary
regularity issue in the case of non-constant $n_{0}$.

If $\|n_{0}\|_{L^{\infty}}\geq 1$, long-time asymptotics of the $\rho$-patches
can be rather complicated, as it does not rule out that $n_{0}$ could be less
than $1$ in some areas. It is not even clear whether the total mass of the
tumor would diverge. Suitable conditions need to be imposed on $n_{0}$ in
order to make the question of uniform regularity more meaningful.

For $\|n_{0}\|_{L^{\infty}}<1$, Theorem 4.10 states that $w(\cdot,t)$ monotone
increases to converge to $w_{\infty}$, which features bounded support. Thus,
the pressure as well as $n$ vanishes in $\Omega_{t}$ as time tends to
infinity, and the regularizing effect of the pressure variable vanishes over
time. On the other hand, under suitable assumptions, $\\{w_{\infty}>0\\}$
features smooth free boundary. To see this, recall that $w_{\infty}$ solves
the obstacle problem (c.f. (4.6))

$\Delta w_{\infty}=(1-\rho_{0}-n_{0})\chi_{\\{w_{\infty}>0\\}}.$

When $\rho_{0}=\chi_{\Omega_{0}}$ with $\Omega_{0}$ having smooth (say
$C^{1,1}$) boundary, the set $\\{w_{\infty}>0\\}$ lies strictly outside of the
support of $\rho_{0}$, and thus near its free boundary $w_{\infty}$ solves
$\Delta w_{\infty}=(1-n_{0})\chi_{\\{w_{\infty}>0\\}}$. It follows that, in
the setting of Corollary 5.15(b), $\partial\\{w_{\infty}>0\\}$ is $C^{\infty}$
provided that $n_{0}$ is smooth. Nevertheless, it remains open whether one can
use this asymptotic regularity of the free boundary to show uniform regularity
of $\partial\\{w(\cdot,t)>0\\}$ in time.

##  Appendix A The Proof of Lemma 3.11

The proof closely follows the argument in [PQV14, Sections 3 and 5] (also see
[GKM22, Proposition 5.1]), with some extra efforts for handling unboundedness
of the spatial domain.

###### Proof.