# Hölder regularity results for parabolic nonlocal double phase problems

Jacques Giacomoni 1 111e-mail<EMAIL_ADDRESS>Deepak Kumar
2222e-mail<EMAIL_ADDRESS>and K. Sreenadh 3333 e-mail:
<EMAIL_ADDRESS>

${}^{1}\,$Université de Pau et des Pays de l’Adour, LMAP (UMR E2S-UPPA CNRS
5142)
Bat. IPRA, Avenue de l’Université F-64013 Pau, France
${}^{2}\,$Research Institute of Mathematics, Seoul National University,
Seoul 08826, South Korea
${}^{3}\,$Department of Mathematics, Indian Institute of Technology Delhi,
Hauz Khaz, New Delhi-110016, India

###### Abstract

In this article, we obtain higher Hölder regularity results for weak solutions
to nonlocal problems driven by the fractional double phase operator

$\displaystyle\mathcal{L}u(x):=$ $\displaystyle 2\;{\rm
P.V.}\int_{\mathbb{R}^{N}}\frac{|u(x)-u(y)|^{p-2}(u(x)-u(y))}{|x-y|^{N+ps_{1}}}dy$
$\displaystyle+2\;{\rm
P.V.}\int_{\mathbb{R}^{N}}a(x,y)\frac{|u(x)-u(y)|^{q-2}(u(x)-u(y))}{|x-y|^{N+qs_{2}}}dy,$

where $1<p\leq q<\infty$, $0<s_{2},s_{1}<1$ and the modulating coefficient
$a(\cdot,\cdot)$ is a non-negative bounded function. Specifically, we prove
higher space-time Hölder continuity result for weak solutions of time
depending nonlocal double phase problems for a particular subclass of the
modulating coefficients. Using suitable approximation arguments, we further
establish higher (global) Hölder continuity results for weak solutions to the
stationary problems involving the operator $\mathcal{L}$ with modulating
coefficients that are locally continuous.

Key words: Nonlocal double phase operators, parabolic equations, higher
regularity results, Hölder continuity.

2010 Mathematics Subject Classification: 35B65, 35J60, 35K55, 35R11.

##  1 Introduction

In this article, our aim is to obtain space-time Hölder continuity results for
the local weak solutions to the following nonlocal parabolic problem:

$\partial_{t}u+\mathcal{L}u=f\;\;\text{ in }\Omega\times I,$
(${\mathcal{P}}_{t}$)

where $\Omega\subset\mathbb{R}^{N}$ is a bounded domain,
$I:=(t_{0},t_{1}]\subset\mathbb{R}$, $f\in
L^{\infty}_{\mathrm{loc}}(\Omega\times I)$ and the nonlocal operator
$\mathcal{L}$ is defined as

$\displaystyle\mathcal{L}u(x):=$ $\displaystyle 2\;{\rm
P.V.}\int_{\mathbb{R}^{N}}\frac{|u(x)-u(y)|^{p-2}(u(x)-u(y))}{|x-y|^{N+ps_{1}}}dy$
$\displaystyle+2\;{\rm
P.V.}\int_{\mathbb{R}^{N}}a(x,y)\frac{|u(x)-u(y)|^{q-2}(u(x)-u(y))}{|x-y|^{N+qs_{2}}}dy,$

with $1<p\leq q<\infty$, $0<s_{2},s_{1}<1$ and $a(\cdot,\cdot)\in
L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ is a non-negative symmetric
function.

The operator $\mathcal{L}$ is known as the fractional double phase operator
(as introduced in [17]), which is a nonlocal analogue of the double phase
operator $\mathcal{D}u:={\rm div}(|\nabla u|^{p-2}\nabla u+a(x)|\nabla
u|^{q-2}\nabla u)$ with $1<p<q$ and $0\leq a(\cdot)$ bounded. The zero set of
the modulating coefficient $a(\cdot)$ dictates the ellipticity of operator
$\mathcal{D}$, that is, when $a(x)>0$, the problem has the $q$-Laplacian
behaviour while on the zero set it has the $p$-Laplacian behaviour. The
corresponding energy functional falls into the category of the so-called
functionals with nonstandard growth conditions of $(p,q)$-type, according to
Marcellini’s terminology [33]. These kinds of functionals involve integrals of
the form

$\displaystyle\mathcal{E}(u)=\int_{\Omega}F(x,\nabla u(x))~{}dx,$

where the energy density $F$, satisfies the unbalanced growth conditions

$c_{1}|\xi|^{p}\leq|F(x,\xi)|\leq c_{2}(|\xi|^{q}+1),\quad 1\leq p\leq q.$

The physical significance of these models lies in the field of nonlinear
elasticity, specifically in homogenization theory. A particular form of the
functional $\mathcal{E}$ is the double phase energy functional given by

$u\mapsto\int_{\Omega}(|\nabla u|^{p}+a(x)|\nabla u|^{q})dx,\quad 0\leq
a(x)\leq M,\ 1<p<q.$

This functional was first introduced by Zhikov in [42], to model the
Lavrentiev phenomenon on strongly anisotropic materials.

The regularity for weak solutions to problems involving the operator
$\mathcal{D}$ or the minimizers of corresponding functionals has been
extensively studied in recent years, see [3, 7, 13, 16, 25, 33] and the survey
article [35]. Particularly, it is known that for $q\leq p+\beta$ (where
$a(\cdot)\in C^{0,\beta}$), the gradient of the bounded weak solutions (or the
minimizers) are locally Hölder continuous. Concerning the regularity results
for time dependent problems involving nonstandard growth operators, we refer
to [15, 26]. For the existence results of problems involving these kinds of
operators, we mention [2, 4, 22, 34, 40].

In the stationary case and when $a\equiv 0$, the problem (${\mathcal{P}}_{t}$)
is driven by the fractional $p$-Laplacian ($(-\Delta)_{p}^{s}$), which is a
nonlinear analogue of the well known fractional Laplacian. The study of
problems involving fractional Sobolev spaces and corresponding nonlocal
operators has gained ample attention due to their wide applications in the
real world problems, such as game theory, image processing, anomalous
diffusion, finance, conservation laws, phase transition, obstacle problems and
material science, see [19]. The regularity results for weak solutions to
problems involving these kinds of nonlocal operators (or minimizers of
corresponding energy functional) are well known, see for instance [5, 10, 14,
18, 28, 36] and their references. Specifically, the weak solution is known to
be $C^{0,s}$-regular up to the boundary, which is optimal.

As far as the nonlocal problems involving the fractional double phase operator
is concerned, the corresponding regularity theory is a relatively new topic
and has not been explored much barring a few works in the stationary case. In
this direction, we mention [17], where De Filippis and Palatucci obtained an
interior Hölder continuity result for viscosity solutions to the problem
$\mathcal{L}u=f\in L^{\infty}$. For the case of weak solutions, in [21], Fang
and Zhang proved the local Hölder continuity result (with unspecified Hölder
exponent) of bounded weak solutions by using the De Giorgi-Nash-Moser type
method, for $qs_{2}\leq ps_{1}$ and a non-negative bounded modulating
coefficient. Here, additionally authors proved that bounded weak solutions are
viscosity solutions. Subsequently, by obtaining a suitable Logarithmic lemma
and applying the De Giorgi-Nash-Moser iteration technique, Byun et al. in [8]
obtained a similar local Hölder regularity result for the case $ps_{1}<qs_{2}$
when the modulating coefficient is assumed to be in the space
$C^{0,\beta}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ with
$qs_{2}-ps_{1}<\beta<1$. We also quote [38] for self improving properties of
the bounded weak solutions to problems involving the operator $\mathcal{L}$,
for the case $qs_{2}\leq ps_{1}$. For the constant modulating coefficient
case, we mention the regularity results for weak solutions obtained in [23,
24, 27]. Particularly, in [23], the almost optimal Hölder continuity result,
up to the boundary, is proved in the super-quadratic case, while in [24], the
global Hölder continuity result is established for all $p,q>1$. See also the
recent work of Chaker et al. [11] for the regularity results of nonlocal
problems with nonstandard growth operators. However, our problem does not fall
into the category of [11] due to assumptions on the nonlinear growth and
fractional differentiability orders.

The regularity theory for time dependent nonlocal problems has been developing
significantly since the last decade. For instance, [9, 30, 12] deals with the
Hölder regularity of weak solutions by following different methods (weak
Harnack inequality or De Giorgi’s type approach) for the linear or non-
degenerate nonlocal operator case. However, the development concerning the
nonlinear nonlocal operators case is recent. In [41], Stömqvist obtained the
local boundedness result for parabolic problems for a large class of nonlocal
problems, whose prototype includes the fractional $p$-Laplacian for $p>2$. For
similar problems and involving a non-zero reaction term, as a parabolic
counterpart for [18], Ding et al. in [20] proved the local boundedness
property for all $p\geq 2N/(N+2s)$. In continuation to [5], Brasco et al. [6]
obtained the higher interior space-time Hölder continuity result with specific
exponents for the local weak solutions to fractional $p$-Laplacian parabolic
problem in the super-quadratic case. Recently, in [1, 32], authors have
obtained the local Hölder regularity results for weak solution to fractional
$p$-Laplacian type problems, for all $p>1$, by means of suitable De Giorgi
type iteration method and intrinsic scaling arguments. Concerning the
nonstandard growth parabolic problems, in [37], authors have obtained a local
boundedness result for variational solutions to a class of nonlocal parabolic
double phase equations.

After these works the regularity results for parabolic problems driven by the
fractional double phase operators were completely open. Addressing this issue,
we establish space-time Hölder continuity result with specific exponents for
the local weak solutions to problem (${\mathcal{P}}_{t}$) when the modulating
coefficient is locally translation invariant in $\Omega\times\Omega$. To
achieve this aim, we first obtain an iterative scheme involving the discrete
differential of the solution with respect to the space variable using Moser’s
iteration technique much in the spirit of [6]. Here, we highlight that
handling the terms involving the exponent $q$ is more tactical and tricky as
compared to the constant modulating coefficient case (see [23]). Subsequently,
we obtain higher time regularity for weak solutions by means of the higher
space regularity estimates. In order to motivate our main regularity results,
we sketch the proof of the existence of weak solutions to the initial value
parabolic problem in the Appendix. The novelty in this regard lies in
introducing a new reflexive Banach space and proving fundamental properties
such coercivity and continuity for the operator involved. Indeed, because of
the varying ellipticity of the operator $\mathcal{L}$, it is not expected that
the solution should lie on the fractional space $W^{s_{1},p}\cap W^{s_{2},q}$.
Additionally, due to the absence of a suitable weighted fractional Poincaré
type inequality with respect to the exponent $q$, we take into account
appropriate fractional Musielak-Sobolev spaces and corresponding time
dependent Banach spaces to define the notion of (local) weak solutions.

Regarding the stationary case, the higher Hölder regularity for solutions to
fractional double phase problems was not explored before. To this aim, we
establish Hölder continuity results with explicit Hölder exponent for weak
solutions to problems involving a more general source term and a modulating
coefficient that is continuous along the diagonal in $\Omega\times\Omega$. For
this, we perform approximation arguments based on the methods presented in [5,
10]. However, the lack of scaling property of the leading operators makes it
difficult to pass the higher Hölder continuity results of problems with
locally translation invariant modulating coefficients and homogeneous right
hand side to our case. To overcome this, we fix the radius of the ball where
Hölder continuity is obtained appropriately depending apriori on the data of
the problem. We emphasize that such results to the best of our knowledge were
not known even for non-homogeneous operators where the modulating coefficients
are bounded away from zero. Finally, we obtain a higher global Hölder
continuity result for the weak solutions to the stationary problem when the
right hand side is bounded, which is new to the current literature on the
subject.

Turning to the layout of the rest of the paper, in Section 2, we mention the
functional framework and our main results. Section 3 contains some technical
results which will be needed in the subsequent sections. In Section 4, we
establish higher space-time Hölder continuity result for locally translation
invariant modulating coefficients. In Section 5, we obtain higher Hölder
continuity results in the stationary case for more general modulating
coefficients and source terms. Finally, in the appendix, we sketch the proof
of the existence of weak solutions to the parabolic problem.

##  2 Functions space framework and main results

###  2.1 Notations

For $1<q<\infty$, we set $[\xi]^{q-1}:=|\xi|^{q-2}\xi$, for all
$\xi\in\mathbb{R}^{N}$ and

$(-\Delta)_{q,a}^{s_{2}}u(x):=2{\rm
P.V.}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x)-u(y)]^{q-1}}{|x-y|^{N+qs_{2}}}dy.$

For $x_{0}\in\mathbb{R}^{N}$ and $v\in L^{1}(B_{r}(x_{0}))$, we set

$(v)_{x_{0},r}:=\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}v(x)dx=\frac{1}{|B_{r}(x_{0})|}\int_{B_{r}(x_{0})}v(x)dx$

and when the center is clear from the context, we will denote it as $(v)_{r}$.
We denote the parabolic cylinders by

$\displaystyle Q_{R,r}(x_{0},t_{0}):=B_{R}(x_{0})\times(t_{0}-r,t_{0}],$

where $B_{R}(x_{0})\subset\mathbb{R}^{N}$ and $r>0$. For simplicity in
writing, we will often use the notation $(\ell,s)\in\\{(p,s_{1}),(q,s_{2})\\}$
together with $a_{\ell}$, where $a_{p}\equiv 1$ and $a_{q}=a$ in
$\mathbb{R}^{N}\times\mathbb{R}^{N}$, unless otherwise specified. We
abbreviate the quantity
$\|a\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}$ by $\|a\|_{\infty}$.
In several places, we will use

$\displaystyle d\mu_{\ell,s}(x,y):=\frac{dxdy}{|x-y|^{N+\ell
s}},\quad\mbox{for }(\ell,s)\in\\{(p,s_{1}),(q,s_{2})\\},$

with $d\mu_{1}=d\mu_{p,s_{1}}$ and $d\mu_{2}=d\mu_{q,s_{2}}$. For
$1<p<\infty$, we denote its Hölder conjugate by $p^{\prime}$, i.e.,
$1/p+1/p^{\prime}=1$. Further, for $s\in(0,1)$, we denote the Critical Sobolev
exponent by $p^{*}_{s}:=Np/(N-ps)$ if $N>ps$, otherwise it is an arbitrarily
large number. The constants $c$ and $C$ appearing in the proofs may vary line
to line. For a Banach space $(X,\|\cdot\|)$, we denote its topological dual by
$X^{*}$.

###  2.2 Functions spaces and definitions

In this subsection, we mention the function spaces used in the paper and
define the notion of local weak solution to problem (${\mathcal{P}}_{t}$). For
$E\subset\mathbb{R}^{N}$, $x,y\in E$, $\xi\in\mathbb{R}$ and
$u:E\to\mathbb{R}$, set

$\displaystyle
H_{a}(x,y,\xi):=\frac{|\xi|^{p}}{|x-y|^{ps_{1}}}+a(x,y)\frac{|\xi|^{q}}{|x-y|^{qs_{2}}}\quad\mbox{and
}$ $\displaystyle\mathcal{H}_{a}(u,E\times E):=\iint_{E\times
E}H_{a}(x,y,|u(x)-u(y)|)\frac{dxdy}{|x-y|^{N}}.$

We define the following space:

$\displaystyle\mathcal{W}_{a}(\mathbb{R}^{N}):=\big{\\{}w\in
L^{p}(\mathbb{R}^{N})\ :\
\mathcal{H}_{a}(w,\mathbb{R}^{N}\times\mathbb{R}^{N})<\infty\big{\\}}$

endowed with the norm

$\displaystyle\|u\|_{\mathcal{W}_{a}(\mathbb{R}^{N})}:=\|u\|_{L^{p}(\mathbb{R}^{N})}+[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})},$

where the Luxembourg seminorm is given by

$\displaystyle[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}:=\inf\big{\\{}\lambda>0\
:\mathcal{H}_{a}\Big{(}\frac{u}{\lambda},\mathbb{R}^{N}\times\mathbb{R}^{N}\Big{)}\leq
1\big{\\}}.$ (2.1)

It is a routine exercise to verify that $\mathcal{W}_{a}(\mathbb{R}^{N})$ is a
uniformly convex Banach space. Moreover, by the properties of modular spaces,
it is easy to prove that

  * (i)

$\min\\{[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}^{p},[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}^{q}\\}\leq\mathcal{H}_{a}(u,\mathbb{R}^{N}\times\mathbb{R}^{N})\leq\max\\{[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}^{p},[u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}^{q}\\}$,

  * (ii)

$[u_{n}-u]_{\mathcal{W}_{a}(\mathbb{R}^{N})}\to 0$ if and only if
$\mathcal{H}_{a}(u_{n}-u,\mathbb{R}^{N}\times\mathbb{R}^{N})\to 0$, as
$n\to\infty$.

For any proper open set $\Omega\subset\mathbb{R}^{N}$, we define the space
$\mathcal{W}_{a}(\Omega)$ as below:

$\displaystyle\mathcal{W}_{a}(\Omega):=\big{\\{}u\in L^{p}(\Omega)\ :\
\int_{\Omega}W_{a}(x)|u(x)|^{q}dx+\mathcal{H}_{a}(u,\Omega\times\Omega)<\infty\big{\\}},$

equipped with the norm

$\displaystyle\|u\|_{\mathcal{W}_{a}(\Omega)}:=\|u\|_{L^{p}(\Omega)}+\|u\|_{L^{q}(W_{a},\>\Omega)}+[u]_{\mathcal{W}_{a}(\Omega)},$

where

$\displaystyle\|u\|_{L^{q}(W_{a},\>\Omega)}^{q}:=\int_{\Omega}W_{a}(x)|u(x)|^{q}dx\quad\mbox{with
}W_{a}(x):=\int_{\mathbb{R}^{N}\setminus\Omega}\frac{a(x,y)}{|x-y|^{N+qs_{2}}}dy,$
(2.2)

and $[u]_{\mathcal{W}_{a}(\Omega)}$ is defined as in (2.1) by replacing
$\mathbb{R}^{N}\times\mathbb{R}^{N}$ with $\Omega\times\Omega$. Note that the
set of functions in $W^{s_{1},p}(\Omega)\cap W^{s_{2},q}(\Omega)$ with compact
support are obviously contained in $\mathcal{W}_{a}(\Omega)$. Moreover, it is
also evident that $\mathcal{W}_{a}(\Omega)$ is continuously embedded into
$W^{s_{1},p}(\Omega)$. Here in after, we will suppress the subscript $a$ from
the definitions of $\mathcal{W}_{a}$ and $W_{a}$, if it is clear from the
context. We call a function $u\in\mathcal{W}_{\rm loc}(\Omega)$ if
$u\in\mathcal{W}(\Omega^{\prime})$, for all $\Omega^{\prime}\Subset\Omega$.

Concerning the time dependent spaces, for $I\subset\mathbb{R}$ and a Banach
space $X$ with a modular $\varrho$, we define

$\displaystyle\mathcal{Y}(I;X):=\bigg{\\{}u(t)\in X\mbox{ for a.e. }t\in
I\mbox{ and }\int_{I}\varrho(u(t))dt<\infty\bigg{\\}}.$

For the space $\mathcal{Y}(I;\mathcal{W}(\Omega))$, we take the modular
$\varrho$ as

$\varrho(u(t)):=\mathcal{H}(u(\cdot,t),\Omega\times\Omega)+\int_{\Omega}\big{(}|u(x,t)|^{p}+W(x)|u(x,t)|^{q}\big{)}dx$

and equip it with the Luxembourg norm on $I$ similar to (2.1).

######  Definition 2.1.

Let $u:\mathbb{R}^{N}\to\mathbb{R}$ be a measurable function,
$0<m,\sigma<\infty$ and $a(\cdot,\cdot)\in
L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ be a non-negative function.
Then, we define the tail space as below:

$\displaystyle L^{m}_{\sigma,a}(\mathbb{R}^{N})=\bigg{\\{}u\in L^{m}_{\rm
loc}(\mathbb{R}^{N}):\sup_{x\in\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{|u(y)|^{m}}{(1+|y|)^{N+\sigma}}dy<\infty\bigg{\\}}$

and for $a(\cdot,\cdot)\equiv 1$, we simply denote it by
$L^{m}_{\sigma}(\mathbb{R}^{N})$. For $m>1$, the nonlocal tail centered at
$x_{0}\in\mathbb{R}^{N}$ with radius $R>0$ is defined as

$\displaystyle
T_{\sigma,a}^{m}(u;x_{0},R)=R^{\sigma}\sup_{x\in\mathbb{R}^{N}}\int_{\mathbb{R}^{N}\setminus
B_{R}(x_{0})}a(x,y)\frac{|u(y)|^{m-1}}{|x_{0}-y|^{N+\sigma}}dy.$

Furthermore, for $J\subset\mathbb{R}$, we define the time depending nonlocal
tail term by

$\displaystyle T_{\infty,\sigma,a}^{m}(u;x_{0},R,J)=R^{\sigma}\sup_{t\in
J}\sup_{x\in\mathbb{R}^{N}}\int_{\mathbb{R}^{N}\setminus
B_{R}(x_{0})}a(x,y)\frac{|u(y,t)|^{m-1}}{|x_{0}-y|^{N+\sigma}}dy.$

For $1<p,q<\infty$, we will often use the notation

$\displaystyle
T_{\infty,\sigma}^{p,q}(u;x_{0},R,J):=T_{\infty,\sigma}^{p}(u;x_{0},R,J)+T_{\infty,\sigma}^{q}(u;x_{0},R,J).$

######  Definition 2.2 (Local weak solution).

Let $f\in(\mathcal{Y}(I;\mathcal{W}(\Omega)))^{*}$, then we say that $u$ is a
local weak solution to problem (${\mathcal{P}}_{t}$) if for all
$J:=[T_{0},T_{1}]\subset I$,

  * (i)

$u\in\mathcal{Y}(J;\mathcal{W}_{\rm loc}(\Omega))\cap
L^{p-1}(J;L^{p-1}_{ps_{1}}(\mathbb{R}^{N}))\cap
L^{q-1}(J;L^{q-1}_{qs_{2},a}(\mathbb{R}^{N}))\cap C(J;L^{2}_{\rm
loc}(\Omega))$ and

  * (ii)

for all $\phi\in\mathcal{Y}(J;\mathcal{W}(\Omega))\cap C^{1}(J;L^{2}(\Omega))$
with compact spatial support contained in $\Omega$, there holds

$\displaystyle-\int_{J}\int_{\Omega}u(x,t)\partial_{t}\phi(x,t)dxdt+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{[u(x,t)-u(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi(x,t)-\phi(y,t))dxdydt$
$\displaystyle\quad+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x,t)-u(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi(x,t)-\phi(y,t))dxdydt$
$\displaystyle=\int_{\Omega}u(x,T_{0})\phi(x,T_{0})dx-\int_{\Omega}u(x,T_{1})\phi(x,T_{1})dx+\langle
f(\cdot,\cdot),\phi(\cdot,\cdot)\rangle_{\mathcal{Y},\mathcal{Y}^{*}}.$ (2.3)

Local weak sub-solution (resp. super-solution) is defined similarly where
((ii)) holds by replacing the sign ”$=$” with ”$\leq$ (resp. $\geq$)” for all
non-negative test functions.

We say that $u$ is locally $\alpha$-Hölder continuous in space, denoted by
$C_{x,{\rm loc}}^{0,\alpha}(\Omega\times I)$, if for any $K\times
J\Subset\Omega\times I$,

$\displaystyle\sup_{t\in J}\;[u(\cdot,t)]_{C^{0,\alpha}(K)}<\infty$

and analogously the local $\gamma$-Hölder continuity in time is defined.

Finally, we recall the following standard inequality: for $\ell\geq 2$, there
exists a constant $c(\ell)>0$ such that

$\displaystyle
c(\ell)|\xi-\zeta|^{\ell}\leq([\xi]^{\ell-1}-[\zeta]^{\ell-1})(\xi-\zeta)\quad\mbox{for
all }\xi,\zeta\in\mathbb{R}.$ (2.4)

###  2.3 Assumptions and main results

We establish a higher Hölder continuity result for locally translation
invariant modulating coefficients. Precisely, we assume the following
assumption:

A.1

The function $a(\cdot,\cdot)$ is locally translation invariant in
$\Omega\times\Omega$, that is, for any $K_{1}\times
K_{2}\Subset\Omega\times\Omega$, there holds $a(x+h,y+h)=a(x,y)$ for all
$(x,y)\in K_{1}\times K_{2}$ and for all $|h|<\min\\{{\rm
dist}(K_{1},\Omega),{\rm dist}(K_{2},\Omega)\\}$.

Before stating our main result in this regard, we fix the following notations:

$\displaystyle\Theta\equiv\Theta(p,s_{1}):=\min\bigg{\\{}1,\frac{ps_{1}}{p-1}\bigg{\\}}$
(2.5)

and

$\displaystyle\Gamma\equiv\Gamma(p,s_{1}):=\begin{cases}\frac{1}{ps_{1}-(p-2)}&\mbox{if
}s_{1}\geq\frac{p-1}{p},\\\ 1&\mbox{if }s_{1}<\frac{p-1}{p}.\end{cases}$ (2.6)

######  Theorem 2.1.

Suppose that $p\geq 2$ and $qs_{2}\leq ps_{1}$. Let the assumption A.1 holds
and let $u$ be a local weak solution to problem (${\mathcal{P}}_{t}$) such
that $u\in L^{\infty}_{\rm loc}(I;L^{\infty}(\mathbb{R}^{N}))$. Then, for all
$\alpha\in(0,\Theta)$ and all $\gamma\in(0,\Gamma)$,

$\displaystyle u\in C^{0,\alpha}_{x,{\rm loc}}(\Omega\times I)\cap
C^{0,\gamma}_{t,{\rm loc}}(\Omega\times I),$

where $\Theta$ and $\Gamma$ are given by (2.5) and (2.6), respectively.
Particularly, for $Q_{2R,2R^{s_{1}p}}\equiv
Q_{2R,2R^{s_{1}p}}(x_{0},\hat{t}_{0})\Subset\Omega\times I$ and for all
$\alpha\in(0,\Theta)$ and $\gamma\in(0,\Gamma)$, there exists a positive
constant $C$ depending only on $N,s_{1},s_{2},p,q,\alpha,\gamma$ and
$\|a\|_{\infty}$ such that

$\displaystyle|u(x,t)-u(y,\tau)|$
$\displaystyle\leq\frac{C}{R^{\alpha}}\Big{(}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[\hat{t}_{0}-R^{s_{1}p},\hat{t}_{0}])}+1+\|f\|^{\frac{1}{p-1}}_{L^{\infty}(Q_{R,R^{s_{1}p}})}\Big{)}^{\frac{qj_{\infty}}{p}}|x-y|^{\alpha}$
$\displaystyle+\frac{C}{R^{ps_{1}\gamma}}\Big{(}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[\hat{t}_{0}-R^{s_{1}p},\hat{t}_{0}])}+1+\|f\|^{\frac{1}{p-1}}_{L^{\infty}(Q_{R,R^{s_{1}p}})}\Big{)}^{\frac{q(\gamma(p-2)+1)j_{\infty}}{p}}|t-\tau|^{\gamma},$

for all $x,y\in B_{\frac{R}{4}}(x_{0})$ and
$t,\tau\in(\hat{t}_{0}-\frac{R^{ps_{1}}}{4},\hat{t}_{0}]$, where
$j_{\infty}\in\mathbb{N}$ depends only on $N,p,s_{1}$ and $\alpha$.

######  Remark 1.

We remark that Theorem 2.1 is true for local weak solutions that are not
bounded in the whole $\mathbb{R}^{N}$ with respect to the space variable. In
this regard, we need a weaker condition, namely $u\in L^{\infty}_{{\rm
loc}}(\Omega\times I)\cap L^{\infty}_{\rm
loc}(I;L^{p-1}_{ps_{1}}(\mathbb{R}^{N}))\cap L^{\infty}_{\rm
loc}(I;L^{q-1}_{qs_{2}}(\mathbb{R}^{N}))$. However, in this case, the value of
$\Gamma$ (the exponent concerning the time regularity) is reduced and given by
$\Gamma:=\min\\{1,1/(ps_{1})\\}$ (see Theorem 4.5 and Remark 2 for further
details).

Next, we consider the following stationary counterpart of problem
(${\mathcal{P}}_{t}$):

$\mathcal{L}u=f\;\;\text{ in }\Omega,$ (${\mathcal{S}}$)

where $f\in L^{\Bbbk}_{\rm loc}(\Omega)$, for $\Bbbk>\max\\{N/(ps_{1}),1\\}$.
In the result below, we remove the restriction of local translation invariant
nature of $a(\cdot,\cdot)$ in $\Omega\times\Omega$, namely

A.2

The function $a(\cdot,\cdot)$ is continuous along the diagonal in
$\Omega\times\Omega$, that is, for any $R>0$ with $B_{R}(x_{0})\Subset\Omega$,
there exists $\epsilon>0$ such that for all $\delta>0$, there exists
$h_{\delta}>0$ satisfying

$\displaystyle\sup_{\begin{subarray}{c}x,y\in B_{R}(x_{0})\\\
|x-y|\leq\epsilon\end{subarray}}|a(x+h,y+h)-a(x,y)|\leq\delta\quad\mbox{for
all }|h|<{h_{\delta}}.$

######  Theorem 2.2.

Suppose that $2\leq p\leq q<\min\\{p^{*}_{s_{1}},ps_{1}/s_{2}\\}$. Let
$a(\cdot,\cdot)$ satisfies the assumption A.2 and let $u$ be a local weak
solution to problem (${\mathcal{S}}$) such that $u\in
L^{q-1}_{qs_{2}}(\mathbb{R}^{N})$. Then, $u\in C^{0,\alpha}_{\rm loc}(\Omega)$
for all $\alpha\in(0,{\Theta_{1}})$, where

$\displaystyle{\Theta_{1}}:=\min\bigg{\\{}\frac{ps_{1}-N/\Bbbk}{p-1},\frac{qs_{2}}{q-1},1\bigg{\\}}.$
(2.7)

Regarding the global Hölder continuity result, we have the following.

######  Corollary 2.3.

Let $\Omega$ be a bounded domain with $C^{1,1}$ boundary and let the
hypotheses of Theorem 2.2 hold. Additionally, if $s_{1}<q^{\prime}s_{2}$,
assume that $a(\cdot,\cdot)\in C^{0,\beta}(\mathbb{R}^{2N})$ with
$\beta>s_{2}$. Let $u\in\mathcal{W}(\mathbb{R}^{N})$ be a weak solution to
problem (${\mathcal{S}}$) with $f\in L^{\infty}(\Omega)$ and the homogeneous
Dirichlet boundary condition in $\mathbb{R}^{N}\setminus\Omega$. Then, $u\in
C^{0,\alpha}(\mathbb{R}^{N})$, for all $\alpha<s_{1}$.

##  3 Some technical results

###  3.1 The case of parabolic problems

We first define the time regularization of the test functions as below.

######  Definition 3.1 (Regularization of test functions).

Let $\varsigma:\mathbb{R}\to\mathbb{R}$ be a non-negative, even smooth
function with compact support in $(-1/2,1/2)$ satisfying
$\int_{\mathbb{R}}\varsigma(\tau)d\tau=1$. For $g\in L^{1}((a,b))$, we define

$\displaystyle
g^{\varepsilon}(t)=\frac{1}{\varepsilon}\int_{t-\varepsilon/2}^{t+\varepsilon/2}\varsigma\Big{(}\frac{t-\sigma}{\varepsilon}\Big{)}g(\sigma)d\sigma=\frac{1}{\varepsilon}\int_{-\varepsilon/2}^{\varepsilon/2}\varsigma\Big{(}\frac{\sigma}{\varepsilon}\Big{)}g(t-\sigma)d\sigma,$

for $t\in(a,b)$ with $\varepsilon\in(0,\min\\{b-t,t-a\\})$.

######  Lemma 3.1 (Caccioppoli inequality).

Let $u$ be a local sub-solution to problem (${\mathcal{P}}_{t}$) such that
$u\in L^{\infty}_{{\rm loc}}(\Omega\times I)\cap L^{\infty}_{\rm
loc}(I;L^{p-1}_{ps_{1}}(\mathbb{R}^{N}))\cap L^{\infty}_{\rm
loc}(I;L^{q-1}_{qs_{2},a}(\mathbb{R}^{N}))$. Suppose that $B_{r}\equiv
B_{r}(x_{0})\Subset\Omega$, for some $r>0$ and $t_{2}<t_{3}$, $\tau>0$
satisfies $[t_{2}-\tau,t_{3}]\subset I:=(t_{0},t_{1}]$. Then, for all non-
negative functions $\psi\in C_{c}^{\infty}(B_{r})$ and $\eta\in
C^{\infty}(\mathbb{R})$ satisfying $\eta\equiv 0$ if $t\leq t_{2}-\tau$ and
$\eta(t)\equiv 1$ if $t\geq t_{2}$, there exists a constant
$C=C(N,s_{1},s_{2},p,q,\|a\|_{\infty})>0$ such that

$\displaystyle\sup_{t_{2}<t<t_{3}}\int_{B_{r}}w_{+}^{2}(x,t)\psi^{q}(x)dx+\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}\int_{B_{r}}\frac{|w_{+}(x,t)\psi^{\frac{q}{p}}(x)-w_{+}(y,t)\psi^{\frac{q}{p}}(y)|^{p}}{|x-y|^{N+s_{1}p}}\eta^{2}(t)dxdydt$
$\displaystyle+\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}\int_{B_{r}}a(x,y)\frac{|w_{+}(x,t)\psi(x)-w_{+}(y,t)\psi(y)|^{q}}{|x-y|^{N+s_{2}q}}\eta^{2}(t)dxdydt$
$\displaystyle\leq
C\sum_{(\ell,s)}\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{|\psi(x)-\psi(y)|^{\ell}}{|x-y|^{N+s\ell}}(w_{+}(x,t)+w_{+}(y,t))^{\ell}\eta^{2}(t)dxdydt$
$\displaystyle\
+C\sum_{(\ell,s)}\sup_{\begin{subarray}{c}t_{2}-\tau<t<t_{3}\\\ x\in{\rm
supp}(\psi)\end{subarray}}\Bigg{(}\int_{\mathbb{R}^{N}\setminus
B_{r}}a_{\ell}(x,y)\frac{w_{+}(y,t)^{\ell-1}}{|x-y|^{N+s\ell}}dy\Bigg{)}\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdt$
$\displaystyle\
+C\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}|f(x,t)|w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdt+C\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}w_{+}^{2}(x,t)\psi^{q}(x)\eta(t)|\partial_{t}\eta(t)|dxdt,$

where $w_{+}=\max\\{u-k,0\\}$, for $k\in\mathbb{R}$.

Proof.  First we assume $t_{3}<t_{1}$ and take

$\displaystyle
0<\epsilon<\frac{\epsilon_{0}}{2}:=\frac{1}{4}\min\\{t_{2}-\tau-
t_{0},t_{1}-t_{3},t_{3}-t_{2}+\tau\\}.$

Set
$\hat{t}_{2}:=t_{2}-\tau-\epsilon_{0},\quad\hat{t}_{3}:=t_{3}+\epsilon_{0}$.
Let $\phi\in\mathcal{Y}((\hat{t}_{2},\hat{t}_{3});\mathcal{W}(B_{r}))\cap
C^{1}((\hat{t}_{2},\hat{t}_{3});L^{2}(B_{r}))$ be such that its spatial
support is compactly contained in $B_{r}$. Let $\phi^{\epsilon}$ be the
regularization of $\phi$, as defined in Definition 3.1. Thus, from ((ii)) (for
sub-solution), we have

$\displaystyle-\int_{J}\int_{B_{r}}u(x,t)\partial_{t}\phi^{\epsilon}(x,t)dxdt+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{[u(x,t)-u(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt$
$\displaystyle\quad+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x,t)-u(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt$
$\displaystyle\leq\int_{B_{r}}u(x,T_{2})\phi^{\epsilon}(x,T_{2})dx-\int_{B_{r}}u(x,T_{3})\phi^{\epsilon}(x,T_{3})dx+\int_{J}\int_{B_{r}}f(x,t)\phi^{\epsilon}(x,t)dxdt,$
(3.1)

where $J=[T_{2},T_{3}]$ with $T_{2}=t_{2}-\tau$ and $T_{3}\in(t_{2},t_{3}]$ to
be determined later (note that $T_{2}-\epsilon$ and
$T_{3}+\epsilon\in(\hat{t}_{2},\hat{t}_{3})$ for all
$\epsilon<\epsilon_{0}/2$). Then, by the properties of convolution, Fubini’s
Theorem and integration by parts, it follows that

$\displaystyle-\int_{J}\int_{B_{r}}u(x,t)\partial_{t}\phi^{\epsilon}(x,t)dxdt$
$\displaystyle=\int_{T_{2}+\frac{\epsilon}{2}}^{T_{3}-\frac{\epsilon}{2}}\int_{B_{r}}\partial_{t}u^{\epsilon}(x,t)\phi(x,t)dxdt+\Upsilon(\epsilon)$
$\displaystyle\
-\int_{B_{r}}u^{\epsilon}(x,T_{3}-\frac{\epsilon}{2})\phi(x,T_{3}-\frac{\epsilon}{2})dx+\int_{B_{r}}u^{\epsilon}(x,T_{2}+\frac{\epsilon}{2})\phi(x,T_{2}+\frac{\epsilon}{2})dx,$

where

$\displaystyle\Upsilon(\epsilon):=$
$\displaystyle-\int_{B_{r}}\int_{T_{2}-\frac{\epsilon}{2}}^{T_{2}+\frac{\epsilon}{2}}\Bigg{(}\frac{1}{\epsilon}\int_{T_{2}}^{\sigma+\frac{\epsilon}{2}}u(x,t)\varsigma(\frac{\sigma-t}{\epsilon})dt\Bigg{)}\partial_{\sigma}\phi(x,\sigma)d\sigma
dx$
$\displaystyle-\int_{B_{r}}\int_{T_{3}-\frac{\epsilon}{2}}^{T_{3}+\frac{\epsilon}{2}}\Bigg{(}\frac{1}{\epsilon}\int_{\sigma-\frac{\epsilon}{2}}^{T_{3}}u(x,t)\varsigma(\frac{\sigma-t}{\epsilon})dt\Bigg{)}\partial_{\sigma}\phi(x,\sigma)d\sigma
dx.$

Using this in (3.1), we get

$\displaystyle\int_{T_{2}+\frac{\epsilon}{2}}^{T_{3}-\frac{\epsilon}{2}}\int_{B_{r}}\partial_{t}u^{\epsilon}(x,t)\phi(x,t)dxdt+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{[u(x,t)-u(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt$
$\displaystyle+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x,t)-u(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt+\Upsilon(\epsilon)$
$\displaystyle\leq\int_{B_{r}}u(x,T_{2})\phi^{\epsilon}(x,T_{2})dx-\int_{B_{r}}u(x,T_{3})\phi^{\epsilon}(x,T_{3})dx+\int_{J}\int_{B_{r}}f(x,t)\phi^{\epsilon}(x,t)dxdt$
$\displaystyle\
+\int_{B_{r}}u^{\epsilon}(x,T_{3}-\frac{\epsilon}{2})\phi(x,T_{3}-\frac{\epsilon}{2})dx-\int_{B_{r}}u^{\epsilon}(x,T_{2}+\frac{\epsilon}{2})\phi(x,T_{2}+\frac{\epsilon}{2})dx.$

We take $\phi(x,t)=v^{\epsilon}(x,t)\psi^{q}(x)\eta^{2}(t)$ in the above
equation, where $v(x,t):=(u^{\epsilon}-k)_{+}(x,t)$. Thus, for
$(\ell,s)\in\\{(p,s_{1}),(q,s_{2})\\}$,

$\displaystyle\underbrace{\int_{T_{2}+\frac{\epsilon}{2}}^{T_{3}-\frac{\epsilon}{2}}\int_{B_{r}}\partial_{t}(u^{\epsilon}(x,t)-k)_{+}(v^{\epsilon}\psi^{q}\eta^{2})(x,t)dxdt}_{\mathfrak{I}_{1}^{\epsilon}}+\Upsilon(\epsilon)$
$\displaystyle+\sum_{(\ell,s)}\underbrace{\int_{J}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{[w(x,t)-w(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}((v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,t)-(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(y,t))dxdydt}_{\mathfrak{I}_{2}^{\epsilon}(\ell)}$
$\displaystyle+2\sum_{(\ell,s)}\underbrace{\int_{J}\int_{\mathbb{R}^{N}\setminus
B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{[w(x,t)-w(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,t)dxdydt}_{\mathfrak{I}_{3}^{\epsilon}(\ell)}$
$\displaystyle\leq\underbrace{\int_{J}\int_{B_{r}}f(x,t)(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,t)dxdt}_{\mathfrak{I}_{4}^{\epsilon}}$
$\displaystyle\
+\underbrace{\int_{B_{r}}u^{\epsilon}(x,T_{3}-\frac{\epsilon}{2})(v^{\epsilon}\psi^{q}\eta^{2})(x,T_{3}-\frac{\epsilon}{2})dx-\int_{B_{r}}u(x,T_{3})(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,T_{3})dx}_{\mathfrak{I}_{5}^{\epsilon}}$
$\displaystyle\
+\underbrace{\int_{B_{r}}u(x,T_{2})(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,T_{2})dx-\int_{B_{r}}u^{\epsilon}(x,T_{2}+\frac{\epsilon}{2})(v^{\epsilon}\psi^{q}\eta^{2})(x,T_{2}+\frac{\epsilon}{2})dx}_{\mathfrak{I}_{6}^{\epsilon}}.$

An integration by parts formula yields

$\displaystyle\mathfrak{I}_{1}^{\epsilon}=$
$\displaystyle\frac{1}{2}\int_{B_{r}}(u^{\epsilon}(x,T_{3}-\frac{\epsilon}{2})-k)_{+}^{2}\psi^{q}(x)\eta^{2}(T_{3}-\frac{\epsilon}{2})dx$
$\displaystyle-\frac{1}{2}\int_{B_{r}}(u^{\epsilon}(x,T_{2}+\frac{\epsilon}{2})-k)_{+}^{2}\psi^{q}(x)\eta^{2}(T_{2}+\frac{\epsilon}{2})dx$
$\displaystyle-\int_{T_{2}+\frac{\epsilon}{2}}^{T_{3}-\frac{\epsilon}{2}}\int_{B_{r}}(u^{\epsilon}(x,t)-k)_{+}^{2}\psi^{q}(x)\eta(t)\partial_{t}\eta(t)dxdt.$

Noting $u\in C([\hat{t}_{2},\hat{t}_{3}];L^{2}(B_{r}))$, we pass through the
limit as $\epsilon\to 0$ in the above equation to get

$\displaystyle\mathfrak{I}_{1}^{\epsilon}\to\mathfrak{I}_{1}:=$
$\displaystyle\frac{1}{2}\int_{B_{r}}w_{+}(x,T_{3})^{2}\psi^{q}(x)\eta^{2}(T_{3})dx-\frac{1}{2}\int_{B_{r}}w_{+}(x,T_{2})^{2}\psi^{q}(x)\eta^{2}(T_{2})dx$
$\displaystyle-\int_{J}\int_{B_{r}}w_{+}(x,t)^{2}\psi^{q}(x)\eta(t)\partial_{t}\eta(t)dxdt.$

To estimate $\mathfrak{I}_{2}^{\epsilon}(\ell)$, using some elementary
properties of convolution, we can prove that the sequence
$\\{(v^{\epsilon}\eta^{2})^{\epsilon}\psi^{q}\\}_{\epsilon}$ is bounded in
$\mathcal{Y}(J;\mathcal{W}(B_{r}))$. Indeed,

$\displaystyle\int_{J}[(v^{\epsilon}\eta^{2})^{\epsilon}\psi^{q}]^{\ell}_{W^{s,\ell}_{a_{\ell}}(B_{r})}dt$
$\displaystyle\leq
C\|\psi\|_{L^{\infty}}^{q\ell}\int_{J}[(v^{\epsilon}\eta^{2})^{\epsilon}]^{\ell}_{W^{s,\ell}_{a_{\ell}}(B_{r})}dt$
$\displaystyle\quad+C\|\nabla\psi\|^{\ell}_{L^{\infty}}\|\psi\|^{(q-1)\ell}_{L^{\infty}}\int_{J}\|(v^{\epsilon}\eta^{2})^{\epsilon}\|^{\ell}_{L^{\ell}(a_{\ell},B_{r})}dt$

and then we use the uniform $L^{\infty}$ bound on $v^{\epsilon}$. Therefore,
by the uniform convexity of the space $\mathcal{Y}(J;\mathcal{W}(B_{r}))$, up
to a subsequence, we get

$\displaystyle(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}\rightharpoonup(u-k)_{+}\psi^{q}\eta^{2}\quad\mbox{weakly
in }\mathcal{Y}(J;\mathcal{W}(B_{r})),$

which implies that

$\displaystyle
a_{\ell}(x,y)^{1/\ell}\frac{(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(x,t)-(v^{\epsilon}\psi^{q}\eta^{2})^{\epsilon}(y,t)}{|x-y|^{N/\ell+s}}$

converges weakly in $L^{\ell}(J;L^{\ell}(B_{r}\times B_{r}))$. Moreover, on
account of $u\in\mathcal{Y}(J;\mathcal{W}(B_{r}))$, we have

$\displaystyle
a_{\ell}^{1-1/\ell}(x,y)[w(x,t)-w(y,t)]^{\ell-1}|x-y|^{1-1/\ell}\in
L^{\frac{\ell}{\ell-1}}(J;L^{\frac{\ell}{\ell-1}}(B_{r}\times B_{r})).$

Thus, we get that
$\mathfrak{I}_{2}^{\epsilon}(\ell)\to\mathfrak{I}_{2}(\ell)$, as $\epsilon\to
0$, where

$\displaystyle\mathfrak{I}_{2}(\ell)=\int_{J}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{[w(x,t)-w(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}(w_{+}(x,t)\psi^{q}(x)-w_{+}(y,t)\psi^{q}(y))\eta^{2}(t)dxdydt.$

For $\mathfrak{I}_{3}^{\epsilon}$, we proceed similarly by taking into account
boundedness of the sequence
$\\{(v^{\epsilon}\eta^{2})^{\epsilon}\psi^{q}\\}_{\epsilon}$ in the space
$L^{p-1}(J;L^{p-1}_{ps_{1}}(\mathbb{R}^{N}))\cap
L^{q-1}(J;L^{q-1}_{qs_{2},a}(\mathbb{R}^{N}))$ instead of
$\mathcal{Y}(J;\mathcal{W}(B_{r}))$. Consequently,

$\displaystyle\mathfrak{I}_{3}^{\epsilon}(\ell)\to\mathfrak{I}_{3}(\ell):=\int_{J}\int_{\mathbb{R}^{N}\setminus
B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{[w(x,t)-w(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdydt.$

By using similar reasoning, we can prove that $\Upsilon(\epsilon)\to 0$, as
$\epsilon\to 0$. Further, on account of properties of convolution, we get
$\mathfrak{I}_{5}^{\epsilon}\to 0$ and $\mathfrak{I}_{6}^{\epsilon}\to 0$ as
$\epsilon\to 0$. Additionally,

$\displaystyle\mathfrak{I}_{4}^{\epsilon}\to\mathfrak{I}_{4}:=\int_{J}\int_{B_{r}}f(x,t)w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdt.$

Therefore,

$\displaystyle\mathfrak{I}_{1}+\mathfrak{I}_{2}+\mathfrak{I}_{3}\leq\mathfrak{I}_{4}.$
(3.2)

Step 1: Estimate of $\mathfrak{I}_{1}$. Noting the assumption that
$\eta(T_{2})=0$ and $\eta\equiv 1$ in $[t_{2},t_{3}]$, we get

$\displaystyle\mathfrak{I}_{1}=\frac{1}{2}\int_{B_{r}}w_{+}(x,T_{3})^{2}\psi^{q}(x)dx-\int_{T_{2}}^{T_{3}}\int_{B_{r}}w_{+}(x,t)^{2}\psi^{q}(x)\eta(t)\partial_{t}\eta(t)dxdt.$
(3.3)

Step 2: Estimate of $\mathfrak{I}_{2}(\ell)$. Setting
$\psi_{\ell}=\psi^{q/\ell}$, we proceed exactly as in [23, Lemma 3.2] to get

$\displaystyle\mathfrak{I}_{2}(\ell)\geq$
$\displaystyle\frac{1}{4}\int_{T_{2}}^{T_{3}}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{|w_{+}(x,t)-w_{+}(y,t)|^{\ell}}{|x-y|^{N+s\ell}}\eta^{2}(t)(\psi_{\ell}(x)+\psi_{\ell}(y))^{\ell}dxdydt$
$\displaystyle-C\int_{T_{2}}^{T_{3}}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{|\psi_{\ell}(x)-\psi_{\ell}(y)|^{\ell}}{|x-y|^{N+s\ell}}(w_{+}(x,t)+w_{+}(y,t))^{\ell}\eta^{2}(t)dxdydt.$
(3.4)

Step 3: Estimate of $\mathfrak{I}_{3}(\ell)$. Observe that

$\displaystyle[w(x,t)-w(y,t)]^{\ell-1}w_{+}(x,t)\geq-
w_{+}(y,t)^{\ell-1}w_{+}(x,t).$

Thus, from the definition of $T_{2}$, we deduce that

$\displaystyle\mathfrak{I}_{3}(\ell)\geq-C\sup_{\begin{subarray}{c}t_{2}-\tau<t<T_{3}\\\
x\in supp(\psi)\end{subarray}}\Bigg{(}\int_{\mathbb{R}^{N}\setminus
B_{r}}a_{\ell}(x,y)\frac{w_{+}(y,t)^{\ell-1}dy}{|x-y|^{N+s\ell}}\Bigg{)}\int_{t_{2}-\tau}^{T_{3}}\int_{B_{r}}w_{+}(x,t)\psi^{p}(x)\eta^{2}(t)dxdt.$
(3.5)

Now, noting

$\displaystyle|w_{+}(x,t)\psi_{\ell}(x)-w_{+}(y,t)\psi_{\ell}(y)|^{\ell}$
$\displaystyle\leq
2^{\ell-1}|w_{+}(x,t)-w_{+}(y,t)|^{\ell}(\psi_{\ell}(x)^{\ell}+\psi_{\ell}(y)^{\ell})$
$\displaystyle\
+2^{\ell-1}(w_{+}(x,t)^{\ell}+w_{+}(y,t)^{\ell})|\psi_{\ell}(x)-\psi_{\ell}(y)|^{\ell}$

and combining (3.3), (3.1) and (3.5) in (3.2), we get

$\displaystyle\int_{t_{2}-\tau}^{T_{3}}\int_{B_{r}}\int_{B_{r}}\frac{|w_{+}(x,t)\psi_{p}(x)-w_{+}(y,t)\psi_{p}(y)|^{p}}{|x-y|^{N+s_{1}p}}\eta^{2}(t)dxdydt+\int_{B_{r}}w_{+}^{2}(x,T_{3})\psi^{p}(x)dx$
$\displaystyle+\int_{t_{2}-\tau}^{T_{3}}\int_{B_{r}}\int_{B_{r}}a(x,y)\frac{|w_{+}(x,t)\psi(x)-w_{+}(y,t)\psi(y)|^{q}}{|x-y|^{N+s_{2}q}}\eta^{2}(t)dxdydt$
$\displaystyle\leq
C\sum_{(\ell,s)}\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{|\psi_{\ell}(x)-\psi_{\ell}(y)|^{\ell}}{|x-y|^{N+s\ell}}(w_{+}(x,t)+w_{+}(y,t))^{\ell}dxdydt$
$\displaystyle\
+C\sum_{(\ell,s)}\sup_{\begin{subarray}{c}t_{2}-\tau<t<t_{3}\\\ x\in
supp(\psi)\end{subarray}}\Bigg{(}\int_{\mathbb{R}^{N}\setminus
B_{r}}a_{\ell}(x,y)\frac{w_{+}(y,t)^{\ell-1}dy}{|x-y|^{N+s\ell}}\Bigg{)}\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdt$
$\displaystyle\
+C\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}|f(x,t)|w_{+}(x,t)\psi^{q}(x)\eta^{2}(t)dxdt+C\int_{t_{2}-\tau}^{t_{3}}\int_{B_{r}}w_{+}^{2}(x,t)\psi^{q}(x)\eta(t)|\partial_{t}\eta(t)|dxdt.$
(3.6)

To get the required result of the Lemma, we take $T_{3}\in(t_{2},t_{3}]$ with
the property

$\displaystyle\int_{B_{r}}w_{+}^{2}(x,T_{3})\psi^{p}(x)dx=\sup_{t_{2}<t<t_{3}}\int_{B_{r}}w_{+}^{2}(x,t)\psi^{p}(x)dx$

and once again we take $T_{3}=t_{3}$, then we sum the resulting inequalities.
This proves the Lemma for $t_{3}<t_{1}$. Note that the constant $C$ appearing
in (3.1) does not depend on $t_{3}$. Therefore, by following an approximation
argument much in the spirit of [6, Proposition 4.1, Step 7], we get a similar
result for the case $t_{3}=t_{1}$. This completes the proof of the lemma.
$\square$

######  Corollary 3.2.

Suppose that $qs_{2}\leq ps_{1}$ and $R\in(0,1)$. Let $u$ be a local weak
solution to

$\displaystyle\partial_{t}u+\mathcal{L}u=f\quad\mbox{in
}B_{2}\times(-2R^{ps_{1}},0]$

with $f\in L^{\infty}(B_{1}\times(-R^{ps_{1}},0])$ such that $u\in
L^{\infty}(\mathbb{R}^{N}\times(-R^{ps_{1}},0])$. Then,

$\displaystyle\Big{(}R^{-N}\int_{-\frac{7}{8}R^{s_{1}p}}^{0}[u]_{W^{s_{1},p}(B_{R})}^{p}dt\Big{)}^{\frac{1}{q}}\leq
C\big{(}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-R^{s_{1}p},0])}+1+\|f\|_{L^{\infty}(B_{1}\times[-R^{s_{1}p},0])}\big{)},$

where $C=C(N,s_{1},s_{2},p,\|a\|_{\infty})>0$ is a constant.

Proof.  Set

$\displaystyle
K_{0}:=\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-R^{s_{1}p},0])}+1+\|f\|_{L^{\infty}(B_{1}\times[-R^{s_{1}p},0])}\quad\mbox{and
}\tilde{u}=u+K_{0}.$

Then, $\tilde{u}$ is still a local weak solution and satisfies $\tilde{u}\geq
1$ in $\mathbb{R}^{N}\times[-R^{ps_{1}},0]$. Now, applying Lemma 3.1 for the
choice $t_{2}=-7R^{s_{1}p}/8$, $\tau=R^{s_{1}p}/8$, $t_{3}=0$ together with
$k=0$ and the test functions $\psi\in C^{\infty}_{c}(B_{R})$, $\eta\in
C_{c}^{\infty}(\mathbb{R})$ satisfying

$\displaystyle\psi\equiv 1\mbox{ in }B_{3R/4},\ 0\leq\psi\leq 1,\
|\nabla\psi|\leq CR^{-1};$ $\displaystyle\eta\equiv 1\mbox{ in
}[-7R^{ps_{1}}/8,0],\ 0\leq\eta\leq 1,\ |\eta^{\prime}|\leq CR^{-ps_{1}},$
(3.7)

we obtain

$\displaystyle\int_{-\frac{7}{8}R^{s_{1}p}}^{0}[\tilde{u}]_{W^{s_{1},p}(B_{3R/4})}^{p}dt$
$\displaystyle\leq
C\sum_{(\ell,s)}\int_{-R^{s_{1}p}}^{0}\int_{B_{R}}\int_{B_{R}}a_{\ell}(x,y)\frac{|\psi(x)-\psi(y)|^{\ell}}{|x-y|^{N+s\ell}}(\tilde{u}(x,t)+\tilde{u}(y,t))^{\ell}dxdydt$
$\displaystyle\ +C\sum_{(\ell,s)}\sup_{\begin{subarray}{c}-R^{s_{1}p}<t<0\\\
x\in supp(\psi)\end{subarray}}\Bigg{(}\int_{\mathbb{R}^{N}\setminus
B_{R}}a_{\ell}(x,y)\frac{u(y,t)_{+}^{\ell-1}dy}{|x-y|^{N+s\ell}}\Bigg{)}\int_{-R^{s_{1}p}}^{0}\int_{B_{R}}\tilde{u}(x,t)\psi^{\ell}(x)dxdt$
$\displaystyle\
+C\int_{-R^{s_{1}p}}^{0}\int_{B_{R}}|f(x,t)|\tilde{u}(x,t)\psi^{q}(x)dxdt+C\int_{-R^{s_{1}p}}^{0}\int_{B_{R}}\tilde{u}^{2}(x,t)\psi^{q}(x)|\partial_{t}\eta(t)|dxdt.$

Using (3.1), it is easy to deduce that

$\displaystyle\int_{-\frac{7}{8}R^{s_{1}p}}^{0}[\tilde{u}]_{W^{s_{1},p}(B_{3R/4})}^{p}dt$
$\displaystyle\leq
C\|a\|_{\infty}\sum_{(\ell,s)}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-R^{s_{1}p},0])}^{\ell}R^{N}R^{ps_{1}-\ell
s}$
$\displaystyle\quad+CR^{N+ps_{1}}\|f\|_{L^{\infty}(B_{1}\times[-R^{s_{1}p},0])}+CR^{N}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-R^{s_{1}p},0])}^{2}$
$\displaystyle\leq CK_{0}^{q}R^{N},$

where we have used $qs_{2}\leq ps_{1}$. This proves the corollary. $\square$

###  3.2 The case of stationary problems

We mention some preliminary results for the stationary problem. We start with
the following Caccioppoli type inequality without any non-negativity
assumption on the solution. The proof runs along the similar lines of [23,
Lemma 3.2] (see also [8, Lemma 4.2]).

######  Lemma 3.3.

Let $u$ be a local weak solution to problem (${\mathcal{S}}$) such that either
$u\in L^{\infty}_{{\rm loc}}(\Omega)$ or $q\leq p^{*}_{s_{1}}$. Then, for
$B_{r}\equiv B_{r}(x_{0})\Subset\Omega$ with $r>0$ and for all $\psi\in
C_{c}^{\infty}(B_{r})$ non-negative, we have

$\displaystyle\int_{B_{r}}\int_{B_{r}}\frac{|u(x)\psi^{\frac{q}{p}}(x)-u(y)\psi^{\frac{q}{p}}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy$
$\displaystyle\leq
C\int_{B_{r}}|f(x)||u(x)|\psi^{q}(x)dx+C\sum_{(\ell,s)}\int_{B_{r}}\int_{B_{r}}a_{\ell}(x,y)\frac{|\psi(x)-\psi(y)|^{\ell}}{|x-y|^{N+s\ell}}(|u(x)|+|u(y)|)^{\ell}dxdy$
$\displaystyle\ +C\sum_{(\ell,s)}\int_{B_{r}}\int_{\mathbb{R}^{N}\setminus
B_{r}}a_{\ell}(x,y)\frac{|u(y)|^{\ell-1}|u(x)|+|u(y)||u(x)|^{\ell-1}}{|x-y|^{N+s\ell}}\psi^{q}(x)dydx.$

Subsequently, following the proof of [24, Proposition 3.1], we have the
following local boundedness result (see also [8, Theorem 1.1]).

######  Lemma 3.4.

Suppose that $q\leq p^{*}_{s_{1}}$. Let $u$ be a local weak solution to
problem (${\mathcal{S}}$) in $\Omega$. Then, $u\in L^{\infty}_{{\rm
loc}}(\Omega)$.
Additionally, suppose that $q<p^{*}_{s_{1}}$. Then, for $B_{r}\equiv
B_{r}(x_{0})\Subset\Omega$ with $r\in(0,1)$, there exists a constant
$C=C(N,s_{1},p,s_{2},q,\Bbbk,\|a\|_{\infty})>0$, as a non-decreasing function
of $\|a\|_{\infty}$, such that

$\displaystyle\sup_{B_{3r/4}}|u|\leq
C\Big{[}\Big{(}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}}|u(x)|^{\vartheta}dx\Big{)}^{\frac{q}{p\vartheta}}+T^{p}_{ps_{1}}(u;x_{0},\frac{3r}{4})^{\frac{1}{p-1}}+T^{q}_{qs_{2},a}(u;x_{0},\frac{3r}{4})^{\frac{1}{q-1}}+\|f\|_{L^{\Bbbk}(B_{r})}^{\frac{1}{p-1}}+1\Big{]},$

where $\vartheta:=\max\\{q,p\sigma\\}$ with
$\sigma=\frac{p\Bbbk-(p^{*}_{s_{1}})^{\prime}}{p(\Bbbk-(p^{*}_{s_{1}})^{\prime})}\Big{(}<\frac{p^{*}_{s_{1}}}{p}\Big{)}$
if $\Bbbk<\infty$, while $\sigma=1$ if $\Bbbk=\infty$.

Now, we mention the notion of a weak solution to the exterior data Dirichlet
problems and related energy spaces. Let $E$ and $E^{\prime}$ are bounded
domains such that $E\Subset E^{\prime}\subset\mathbb{R}^{N}$. For given $g\in
L^{p-1}_{ps_{1}}(\mathbb{R}^{N})\cap L^{q-1}_{qs_{2},a}(\mathbb{R}^{N})$, we
define

$\displaystyle X_{g,a}(E,E^{\prime}):=\\{v\in\mathcal{W}_{a}(E^{\prime})\cap
L^{p-1}_{ps_{1}}(\mathbb{R}^{N})\cap
L^{q-1}_{qs_{2},a}(\mathbb{R}^{N}):v=g\quad\mbox{a.e. in
}\mathbb{R}^{N}\setminus E\\},$

equipped with the norm of $\mathcal{W}_{a}(E^{\prime})$. For $f\in
L^{\Bbbk}(E^{\prime})$ with $\Bbbk>\max\\{N/(ps_{1}),1\\}$ and
$g\in\mathcal{W}_{a}(E^{\prime})\cap L^{p-1}_{s_{1}p}(\mathbb{R}^{N})\cap
L^{q-1}_{s_{2}q,a}(\mathbb{R}^{N})$, consider the problem:

$\displaystyle\mathcal{L}v=f\quad\mbox{in }E;\quad v=g\quad\mbox{in
}\mathbb{R}^{N}\setminus E.$ (3.8)

######  Definition 3.2.

A function $v\in X_{g,a}(E,E^{\prime})$ is said to be a weak solution of
problem (3.8) if, for all $\phi\in X_{0,a}(E,E^{\prime})$,

$\displaystyle\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\Bigg{(}\frac{[v(x)-v(y)]^{p-1}}{|x-y|^{N+s_{1}p}}+a(x,y)\frac{[v(x)-v(y)]^{q-1}}{|x-y|^{N+s_{2}q}}\Bigg{)}(\phi(x)-\phi(y))dxdy=\int_{E}f(x)\phi(x)dx.$

The existence of a weak solution $v\in X_{g,a}(E,E^{\prime})$ to problem (3.8)
is guaranteed by an analogous result of Theorem A.1. Again, in the sequel, we
will drop the term $a$ from the definition of the space $X_{g,a}$, if the
context is clear. Concerning the boundedness properties of a weak solution $v$
to problem (3.8), we see that since it also satisfies the same Caccioppoli
inequality of Lemma 3.3, the local boundedness result of Lemma 3.4 holds for
$v$, too. For the boundary estimates, following the idea of the proof of [31,
Theorem 5] and Lemma 3.4 (note that $v\in\mathcal{W}(E^{\prime})$), we have

######  Lemma 3.5.

Suppose that $q\leq p^{*}_{s_{1}}$. Let $v\in X_{g}(E,E^{\prime})$ be a weak
solution to (3.8). Additionally, assume that for $x_{0}\in\partial E$ and
$0<r<\min\\{1,{\rm dist}(x_{0},E^{\prime})\\}$, $f\in L^{\Bbbk}(B_{r}(x_{0}))$
and $g\in L^{\infty}(B_{r}(x_{0}))$. Then, $v\in L^{\infty}(B_{r/4}(x_{0}))$.
Furthermore, for the case $q<p^{*}_{s_{1}}$, there exists a constant
$C=C(N,s_{1},s_{2},p,q,\Bbbk,\|a\|_{\infty})>0$ such that

$\displaystyle\|v\|_{L^{\infty}(B_{r/4}(x_{0}))}$ $\displaystyle\leq
C\Big{[}\Big{(}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}|u(x)|^{\vartheta}dx\Big{)}^{\frac{q}{p\vartheta}}+T^{p}_{ps_{1}}(u;x_{0},\frac{r}{4})^{\frac{1}{p-1}}+T^{q}_{qs_{2},a}(u;x_{0},\frac{r}{4})^{\frac{1}{q-1}}+1+\|f\|_{L^{\Bbbk}(B_{r}(x_{0}))}^{\frac{1}{p-1}}\Big{]},$

where $\vartheta$ is given by Lemma 3.4.

Next, we mention the following Hölder continuity result with unspecified
Hölder exponent for local weak solutions to problem (${\mathcal{S}}$).

######  Proposition 3.6.

Suppose that $p\geq 2$ and $qs_{2}\leq ps_{1}$. Let $u$ be a locally bounded
local weak solution to problem (${\mathcal{S}}$) such that $u\in
L^{q-1}_{qs_{2}}(\mathbb{R}^{N})$. Then, there exists $\alpha_{1}>0$ such that
$u\in C^{0,\alpha_{1}}_{{\rm loc}}(\Omega)$. More precisely, if $B_{2R}\equiv
B_{2R}(x_{0})\Subset\Omega$, for some $R\in(0,1)$, then for all $r\in(0,R/2)$,

$\displaystyle[u]_{C^{0,\alpha_{1}}(B_{r/2})}\leq\frac{C}{R^{\alpha_{1}}}\bigg{(}\|u\|_{L^{\infty}(B_{r})}+T^{p}_{ps_{1}}(u;x_{0},\frac{R}{2})^{\frac{1}{p-1}}+T^{q}_{qs_{2}}(u;x_{0},\frac{R}{2})^{\frac{1}{q-1}}+1+\|f\|_{L^{\Bbbk}(B_{r})}^{\frac{1}{p-1}}\bigg{)},$

where $C$ is a positive constant depending only on $N,s_{1},s_{2},p,q,\Bbbk$
and $\|a\|_{\infty}$.

Proof.  Let $v\in X_{u}(B_{3R/2},B_{7R/4})\cap
L^{q-1}_{qs_{2}}(\mathbb{R}^{N})$ be the unique solution to the problem:

$\left\\{\begin{array}[]{rllll}\mathcal{L}v&=0\quad\mbox{in }B_{3R/2},\\\
v&=u\quad\mbox{in }\mathbb{R}^{N}\setminus B_{3R/2}.\end{array}\right.$

Then, by a careful inspection of the proof [21, Lemma 3.5], specifically the
definition of $\omega(r_{0})$, there (see also the proof of [8, Lemma 5.3]),
for all $\rho<3R/4$, we have

$\displaystyle\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{\rho}}|v(x)-(v)_{x_{0},\rho}|^{p}dx\leq
C\Big{(}\frac{\rho}{R}\Big{)}^{p\alpha}\bigg{(}\|v\|_{L^{\infty}(B_{2\rho})}^{p}+\sum_{(\ell,s)}T^{\ell}_{\ell
s,a_{\ell}}(v;x_{0},\rho)^{\frac{p}{\ell-1}}+1\bigg{)},$

where we have further used the relation

$\displaystyle\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{\rho}}|v(x)-(v)_{x_{0},\rho}|^{p}dx\leq
C\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{\rho}}({\rm osc}_{B_{\rho}}v)^{p}dx.$

Moreover, taking $u-v$ as a test function and applying Hölder’s inequality, we
can prove that

$\displaystyle[u-v]_{W^{s_{1},p}(\mathbb{R}^{N})}^{p}\leq
CR^{N\big{(}\frac{p}{\Bbbk(p-1)}-\frac{N-ps_{1}}{N(p-1)}\big{)}}\|f\|_{L^{\Bbbk}(B_{R/2})}^{\frac{1}{p-1}}.$

Consequently, using the local boundedness result of Lemma 3.4, the tail terms
can be estimated as

$\displaystyle T^{\ell}_{\ell s,a_{\ell}}(v;x_{0},R/2)^{\frac{p}{\ell-1}}\leq
CT^{\ell}_{\ell
s,a_{\ell}}(u;x_{0},3R/2)^{\frac{p}{\ell-1}}+C\Big{(}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{3R/2}}|u|^{\ell}dx+\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{3R/2}}|u-v|^{\ell}dx\Big{)}^{p/\ell}.$

Finally, using the relation

$\displaystyle\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}}|u(x)-(u)_{x_{0},r}|^{p}dx\leq
C\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}}|u-v|^{p}dx+C\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}}|v(x)-(v)_{x_{0},r}|^{p}dx,$

we get the required result of the proposition (see [5, Theorem 3.6] for the
further details). $\square$

##  4 Higher interior regularity

In this section, we obtain higher space-time Hölder continuity result and
prove Theorem 2.1.

###  4.1 Spatial regularity: locally translation invariant case

We first fix some notations which will be used in this section. For a
measurable function $u:\mathbb{R}^{N}\to\mathbb{R}$ and $h\in\mathbb{R}^{N}$,
we define

$\displaystyle u_{h}(x)=u(x+h),\ \ \delta_{h}u(x)=u_{h}(x)-u(x),\ \
\delta^{2}_{h}u(x)=\delta_{h}(\delta_{h}u(x))=u_{2h}(x)+u(x)-2u_{h}(x).$

Now we present a regularization lemma involving the discrete differentials.

######  Lemma 4.1.

Let $u$ be a local weak solution to problem (${\mathcal{P}}_{t}$) in
$B_{2}\times(-2,0]$. Suppose that $u\in L^{\infty}(E\times[-1,0])$ for all
$E\Subset B_{2}$. Let $\eta\in C^{0,1}_{c}(B_{2})$ be a non-negative function
and let $\tau\in C^{\infty}(\mathbb{R})$ be such that $0\leq\tau\leq 1$,
$\tau\equiv 0$ in $(-\infty,T_{0}]$ and $\tau\equiv 1$ in $[T_{1},\infty)$,
for some $-1<T_{0}<T_{1}<0$. Then, for any locally Lipschitz function
$\Psi:\mathbb{R}\to\mathbb{R}$ and $h\in\mathbb{R}^{N}$ such that $0<|h|<{\rm
dist}({\rm supp}\eta,\partial B_{2})/4$, we have

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\big{(}[u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\qquad\times(\Psi(\delta_{h}u(x,t))\eta^{p}(x)-\Psi(\delta_{h}u(y,t))\eta^{p}(y))\tau(t)d\mu_{1}dt$
$\displaystyle+\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\big{(}a_{h}(x,y)[u_{h}(x,t)-u_{h}(y,t)]^{q-1}-a(x,y)[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\qquad\times(\Psi(\delta_{h}u(x,t))\eta^{p}(x)-\Psi(\delta_{h}u(y,t))\eta^{p}(y))\tau(t)d\mu_{2}dt$
$\displaystyle+\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u(x,T_{1}))\eta^{p}(x)dxdt$
$\displaystyle=\int_{T_{0}}^{T_{1}}\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u(x,t))\eta^{p}(x)\tau^{\prime}(t)dxdt+\int_{T_{0}}^{T_{1}}\int_{B_{2}}\big{(}f_{h}-f\big{)}(x,t)\Psi(u_{h}(x,t))dxdt,$
(4.1)

where $\widetilde{\Psi}(z)=\int_{0}^{z}\Psi(\sigma)d\sigma$,
$a_{h}(x,y)=a(x+h,y+h)$ and $f_{h}(x,t)=f(x+h,t)$.

Proof.  We only sketch the proof as it is similar to the one provided in [6,
Lemma 3.3]. Set $J:=[T_{0},T_{1}]\subset(-1,0)$ and let
$\phi\in\mathcal{Y}(J;\mathcal{W}(B_{2}))\cap C^{1}(J;L^{2}(B_{2}))$ be such
that its spatial support is compactly contained in $B_{2}$. Let
$\phi^{\epsilon}$ be the regularization of $\phi$, as defined in Definition
3.1, for

$\displaystyle
0<\epsilon<\frac{\epsilon_{0}}{2}:=\frac{1}{4}\min\\{T_{0}+1,-T_{1},T_{1}-T_{0}\\}.$

Then, proceeding similar to the proof of [6, Lemma 3.3], we have

$\displaystyle\int_{T_{0}+\frac{\epsilon}{2}}^{T_{1}-\frac{\epsilon}{2}}\int_{B_{2}}\partial_{t}u^{\epsilon}(x,t)\phi(x,t)dxdt+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{[u(x,t)-u(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt$
$\displaystyle+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x,t)-u(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi^{\epsilon}(x,t)-\phi^{\epsilon}(y,t))dxdydt+\Upsilon_{u}(\epsilon)$
$\displaystyle=\int_{J}\int_{B_{2}}f(x,t)\phi^{\epsilon}(x,t)dxdt+\int_{B_{2}}u(x,T_{0})\phi^{\epsilon}(x,T_{0})dx-\int_{B_{2}}u^{\epsilon}(x,T_{0}+\frac{\epsilon}{2})\phi(x,T_{0}+\frac{\epsilon}{2})dx$
$\displaystyle\
+\int_{B_{2}}u^{\epsilon}(x,T_{1}-\frac{\epsilon}{2})\phi(x,T_{1}-\frac{\epsilon}{2})dx-\int_{B_{2}}u(x,T_{1})\phi^{\epsilon}(x,T_{1})dx,$
(4.2)

where

$\displaystyle\Upsilon_{u}(\epsilon):=$
$\displaystyle-\int_{B_{2}}\Bigg{(}\frac{1}{\epsilon}\int_{T_{0}}^{T_{0}+\epsilon}u(x,t)\varsigma\Big{(}\frac{T_{0}-t}{\epsilon}+\frac{1}{2}\Big{)}dt\Bigg{)}\phi(x,T_{0}+\frac{\epsilon}{2})dx$
$\displaystyle+\int_{B_{2}}\int_{-\frac{1}{2}}^{\frac{1}{2}}\Bigg{(}\int_{-\frac{1}{2}}^{\rho}u(x,\epsilon\rho+T_{0}-\epsilon\sigma)\varsigma^{\prime}(\sigma)d\sigma\Bigg{)}\phi(x,\epsilon\rho+T_{0})d\rho
dx$
$\displaystyle+\int_{B_{2}}\Bigg{(}\frac{1}{\epsilon}\int_{T_{1}-\epsilon}^{T_{1}}u(x,t)\varsigma\Big{(}\frac{T_{1}-t}{\epsilon}-\frac{1}{2}\Big{)}dt\Bigg{)}\phi(x,T_{1}-\frac{\epsilon}{2})dx$
$\displaystyle-\int_{B_{2}}\int_{-\frac{1}{2}}^{\frac{1}{2}}\Bigg{(}\int_{\frac{1}{2}}^{\rho}u(x,\epsilon\rho+T_{1}-\epsilon\sigma)\varsigma^{\prime}(\sigma)d\sigma\Bigg{)}\phi(x,\epsilon\rho+T_{1})d\rho
dx.$ (4.3)

Again testing the weak formulation with $\phi_{-h}(x,t):=\phi(x-h,t)$ and
changing the variables, we see that (4.1) and (4.1) hold for $u_{h}$ with
$a(\cdot,\cdot)$ and $f(x,t)$ replaced by $a_{h}(x,y)$ and $f_{h}(x,t)$,
respectively. On subtracting the resulting equations and taking

$\displaystyle\phi=\Psi(u_{h}^{\epsilon}-u^{\epsilon})\eta^{p}\tau_{\epsilon}=\Psi(\delta_{h}u^{\epsilon})\eta^{p}\tau_{\epsilon},$

where

$\displaystyle\tau_{\epsilon}(t)=\tau\Bigg{(}\frac{T_{1}-T_{0}}{T_{1}-T_{0}-\epsilon}\Big{(}t-T_{1}+\frac{\epsilon}{2}\Big{)}+T_{1}\Bigg{)},$

we obtain

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\big{(}[u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\quad\times\Big{(}\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)-\big{(}\Psi(\delta_{h}u^{\epsilon}(y,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(y)\Big{)}d\mu_{1}dt$
$\displaystyle+\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\big{(}a_{h}(x,y)[u_{h}(x,t)-u_{h}(y,t)]^{q-1}-a(x,y)[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\qquad\times\Big{(}\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)-\big{(}\Psi(\delta_{h}u^{\epsilon}(y,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(y)\Big{)}d\mu_{2}dt$
$\displaystyle+\int_{T_{0}}^{T_{1}}\int_{B_{2}}\partial_{t}(\delta_{h}u^{\epsilon}(x,t))\Psi(\delta_{h}u^{\epsilon}(x,t))\eta^{p}(x)\tau_{\epsilon}(t)dxdt+\big{(}\Upsilon_{u_{h}}(\epsilon)-\Upsilon_{u}(\epsilon)\big{)}$
$\displaystyle=\int_{T_{0}}^{T_{1}}\int_{B_{2}}(f_{h}-f)(x,t)\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)dxdt$
$\displaystyle\
+\int_{B_{2}}\Big{(}\delta_{h}u^{\epsilon}(x,T_{1}-\frac{\epsilon}{2})\Psi(\delta_{h}u^{\epsilon}(x,T_{1}-\epsilon/2))-\delta_{h}u(x,T_{1})\big{(}\Psi(\delta_{h}u^{\epsilon}(x,T_{1}))\big{)}^{\epsilon}\Big{)}\eta^{p}(x)dx.$
(4.4)

Using the integration by parts formula, the definition of $\widetilde{\Psi}$
and convergence arguments similar to the proof of [6, Lemma 3.3], we get, as
$\epsilon\to 0$,

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{B_{2}}\partial_{t}(\delta_{h}u^{\epsilon})\Psi(\delta_{h}u^{\epsilon})\eta^{p}(x)\tau_{\epsilon}(t)dxdt$
$\displaystyle=\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u^{\epsilon}(x,T_{1}))\eta^{p}(x)dx-\int_{T_{0}}^{T_{1}}\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u^{\epsilon})\eta^{p}(x)\tau_{\epsilon}^{\prime}(t)dxdt$
$\displaystyle\to\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u(x,T_{1}))\eta^{p}dxdt-\int_{T_{0}}^{T_{1}}\int_{B_{2}}\widetilde{\Psi}(\delta_{h}u)\eta^{p}\tau^{\prime}dxdt.$

An analogous treatment implies that, as $\epsilon\to 0$,

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{B_{2}}(f_{h}-f)\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)dxdt\to\int_{T_{0}}^{T_{1}}\int_{B_{2}}(f_{h}-f)\Psi(\delta_{h}u(x,t))\tau(t)\eta^{p}(x)dxdt.$

Moreover, similar convergence arguments as in [6, Lemma 3.3] yield that the
second term on the right hand side of (4.1) converges to $0$ along with

$\displaystyle\Upsilon_{u_{h}}(\epsilon)-\Upsilon_{u}(\epsilon)\to
0,\quad\mbox{as }\epsilon\to 0.$

Now, to check the convergence of the first two terms on the left side of
(4.1), we observe that

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\big{(}[u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\quad\times\Big{(}\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)-\big{(}\Psi(\delta_{h}u^{\epsilon}(y,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(y)\Big{)}d\mu_{1}dt$
$\displaystyle=\int_{T_{0}}^{T_{1}}\iint_{B_{2-2h}\times
B_{2-2h}}\big{(}[u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\quad\times\Big{(}\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)-\big{(}\Psi(\delta_{h}u^{\epsilon}(y,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(y)\Big{)}d\mu_{1}dt$
$\displaystyle+2\int_{T_{0}}^{T_{1}}\iint_{B_{2-2h}\times(\mathbb{R}^{N}\times
B_{2-h})}\big{(}[u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1}\big{)}$
$\displaystyle\qquad\qquad\quad\times\big{(}\Psi(\delta_{h}u^{\epsilon}(x,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}(x)d\mu_{1}dt.$

Analogously the second integral on the left side of (4.1) is split into two
parts. As in [6, Appendix B], we obtain that the sequence
$\\{\big{(}\Psi(\delta_{h}u^{\epsilon}(\cdot,t))\tau_{\epsilon}(t)\big{)}^{\epsilon}\eta^{p}\\}$
is bounded in $\mathcal{Y}(J;\mathcal{W}(B_{2-2h}))$ with respect to
$\epsilon$. Then, proceeding similar to the proof of convergence of
$\mathfrak{I}_{2}^{\epsilon}(\ell)$ and $\mathfrak{I}_{3}^{\epsilon}(\ell)$ in
Lemma 3.1, we get the required convergence. $\square$

######  Proposition 4.2.

Suppose that $p\geq 2$, $qs_{2}<ps_{1}+2(1-s_{1})$ and the assumption A.1
holds. Let $u$ be a local weak solution to the following problem:

$\displaystyle\partial_{t}u+(-\Delta)_{p}^{s_{1}}u+\lambda(-\Delta)_{q,a}^{s_{2}}u=f\quad\mbox{in
}B_{2}\times(-2,0],$

where $\lambda>0$ is a parameter and $f\in L^{\infty}(B_{1}\times(-1,0])$,
such that

$\displaystyle\|u\|_{L^{\infty}(B_{1}\times[-1,0])}\leq
1\quad\mbox{and}\quad\sum_{(\ell,s)}T^{\ell}_{\infty,\ell s}(u;0,1,-1,0)\leq
1.$ (4.5)

Additionally, suppose that for some $m\geq p$ and $0<h_{0}<\frac{1}{10}$, we
have

$\displaystyle\int_{T_{0}}^{T_{1}}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}dt<\infty,$

for $4h_{0}<R\leq 1-5h_{0}$ and $-1<T_{0}<T_{1}\leq 0$. Then, for all
$\mu\in(0,T_{1}-T_{0})$, there holds:

$\displaystyle\int_{T_{0}+\mu}^{T_{1}}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{s_{1}}}\bigg{\|}^{m+1}_{L^{m+1}(B_{R-4h_{0}})}dt+\frac{1}{m+3-p}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta_{h}u(\cdot,T_{1})}{|h|^{\frac{(m+2-p)s_{1}}{m+3-p}}}\bigg{\|}^{m+3-p}_{L^{m+3-p}(B_{R-4h_{0}})}$
$\displaystyle\leq
C(1+\lambda+\|f\|_{L^{\infty}(B_{1}\times(-1,0])})\int_{T_{0}}^{T_{1}}\left(\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1\right)dt,$

where $C=C(N,h_{0},p,q,m,s_{1},\|a\|_{\infty})>0$ (which depends inversely on
$h_{0}$).

Proof.  We first assume that $T_{1}<0$ and set $r:=R-4h_{0}$. For $\kappa\geq
1$, we take

$\displaystyle\Psi(z):=[z]^{\kappa}$

in Lemma 4.1, together with

$\displaystyle\eta\equiv 1\quad\mbox{in }B_{r},\ \eta\equiv 0\quad\mbox{in
}B^{c}_{(R+r)/2},\ \ 0\leq\eta\leq 1\quad\mbox{and }|\nabla\eta|\leq
C(R-r)^{-1}=C(4h_{0})^{-1};$ $\displaystyle 0\leq\tau\leq 1,\quad\tau\equiv 0\
\mbox{in }(-\infty,T_{0}],\quad\tau\equiv 1\ \mbox{in
}[T_{0}+\mu,\infty)\quad\mbox{and }|\tau^{\prime}|\leq C\mu^{-1}.$

Dividing (4.1) by $|h|^{1+\nu\kappa}>0$, we have

$\displaystyle\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{([u_{h}(x,t)-u_{h}(y,t)]^{p-1}-[u(x,t)-u(y,t)]^{p-1})}{|h|^{1+\nu\kappa}}$
$\displaystyle\qquad\qquad\qquad\times\big{(}[\delta_{h}u(x,t)]^{\kappa}\eta^{p}(x)-[\delta_{h}u(y,t)]^{\kappa}\eta^{p}(y)\big{)}\tau(t)d\mu_{1}dt$
$\displaystyle+\lambda\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{(a_{h}(x,y)[u_{h}(x,t)-u_{h}(y,t)]^{q-1}-a(x,y)[u(x,t)-u(y,t)]^{q-1})}{|h|^{1+\nu\kappa}}$
$\displaystyle\qquad\qquad\qquad\times\big{(}[\delta_{h}u(x,t)]^{\kappa}\eta^{p}(x)-[\delta_{h}u(y,t)]^{\kappa}\eta^{p}(y)\big{)}\tau(t)d\mu_{2}dt$
$\displaystyle+\frac{1}{\kappa+1}\int_{B_{2}}\frac{|\delta_{h}u(x,T_{1})|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}(x)dx$
$\displaystyle=\frac{1}{\kappa+1}\int_{T_{0}}^{T_{1}}\int_{B_{2}}\frac{|\delta_{h}u(x)|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}(x)\tau^{\prime}(t)dxdt+\int_{T_{0}}^{T_{1}}\int_{B_{2}}\delta_{h}f(x,t)\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}\eta^{p}(x)\tau(t)dxdt,$

where $\delta_{h}f(x,t)=f(x+h,t)-f(x,t)$. We set the following:

$\displaystyle\mathbb{I}_{1,q}(t):=\lambda\int_{B_{R}}\int_{B_{R}}a(x,y)\frac{([u_{h}(x,t)-u_{h}(y,t)]^{q-1}-[u(x,t)-u(y,t)]^{q-1})}{|h|^{1+\nu\kappa}}$
$\displaystyle\qquad\qquad\qquad\times\big{(}[\delta_{h}u(x,t)]^{\kappa}\eta^{p}(x)-[\delta_{h}u(y,t)]^{\kappa}\eta^{p}(y)\big{)}d\mu_{2},$
$\displaystyle\mathbb{I}_{2,q}(t):=2\lambda\int_{B_{\frac{R+r}{2}}}\int_{B_{R}^{c}}\frac{(a(x+h,y+h)[u_{h}(x,t)-u_{h}(y,t)]^{q-1}-a(x,y)[u(x,t)-u(y,t)]^{q-1})}{|h|^{1+\nu\kappa}}$
$\displaystyle\qquad\qquad\qquad\qquad\qquad\times[\delta_{h}u(x,t)]^{\kappa}\eta^{p}(x)d\mu_{2},$
$\displaystyle\mathbb{I}_{3}:=\frac{1}{\kappa+1}\int_{T_{0}}^{T_{1}}\int_{B_{2}}\frac{|\delta_{h}u(x,t)|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}(x)\tau^{\prime}(t)dxdt,$
$\displaystyle\mathbb{I}_{f}:=\int_{T_{0}}^{T_{1}}\int_{B_{2}}|\delta_{h}f(x,t)|\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}\eta^{p}(x)\tau(t)dxdt$

and analogously, $\mathbb{I}_{1,p}(t)$ and $\mathbb{I}_{2,p}(t)$ are defined
by taking $\lambda=1$ and $a\equiv 1$. Therefore, noting $a(x+h,y+h)=a(x,y)$
in $B_{R}\times B_{R}$, we have

$\displaystyle\int_{T_{0}}^{T_{1}}\mathbb{I}_{1,p}(t)\tau(t)dt+\frac{1}{\kappa+1}\int_{B_{2}}\frac{|\delta_{h}u(x,T_{1})|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}dx\leq$
$\displaystyle\int_{T_{0}}^{T_{1}}|\mathbb{I}_{2,q}(t)|\tau(t)dt+\int_{T_{0}}^{T_{1}}|\mathbb{I}_{2,p}(t)|\tau(t)dt$
$\displaystyle-\int_{T_{0}}^{T_{1}}\mathbb{I}_{1,q}(t)\tau(t)dt+\mathbb{I}_{3}+\mathbb{I}_{f}.$
(4.6)

For notational convenience, we suppress the $t$-dependence from $u$ and
$u_{h}$.
Estimate of $\mathbb{I}_{1,p}$. Proceeding similarly to Step 1 of the proof of
[5, Proposition 4.1], we obtain

$\displaystyle\mathbb{I}_{1,p}$ $\displaystyle\geq
c\bigg{[}\frac{[\delta_{h}u]^{\frac{\kappa+p-1}{p}}\eta}{|h|^{\frac{1+\nu\kappa}{p}}}\bigg{]}^{p}_{W^{s_{1},p}(B_{R})}-C\int_{B_{R}}\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\kappa+p-1}+|\delta_{h}u(y)|^{\kappa+p-1}}{|h|^{1+\nu\kappa}}|\eta(x)-\eta(y)|^{p}d\mu_{1}$
$\displaystyle\quad-C\int_{B_{R}}\int_{B_{R}}\left[\frac{|\delta_{h}u(x)|^{\kappa+1}+|\delta_{h}u(y)|^{\kappa+1}}{|h|^{1+\nu\kappa}}\big{(}|u_{h}(x)-u_{h}(y)|^{\frac{p-2}{2}}+|u(x)-u(y)|^{\frac{p-2}{2}}\big{)}^{2}\right.$
$\displaystyle\qquad\qquad\qquad\left.\times\big{|}\eta^{\frac{p}{2}}(x)-\eta^{\frac{p}{2}}(y)\big{|}^{2}\right]d\mu_{1}$
$\displaystyle=:c\bigg{[}\frac{[\delta_{h}u]^{\frac{\kappa+p-1}{p}}\eta}{|h|^{\frac{1+\nu\kappa}{p}}}\bigg{]}^{p}_{W^{s_{1},p}(B_{R})}-\mathbb{I}_{4,p}-\mathbb{I}_{5,p}.$
(4.7)

Moreover, from the proof of [23, Proposition 3.9], we have

$\displaystyle|\mathbb{I}_{4,p}|+|\mathbb{I}_{5,p}|\leq
C\left(\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\right).$ (4.8)

Estimate of $\mathbb{I}_{1,q}$. Following the proof of [23, Equ. (3.16)] (note
that the first term is non-negative there), we obtain

$\displaystyle\mathbb{I}_{1,q}$
$\displaystyle\geq-C\lambda\int_{B_{R}}\int_{B_{R}}a(x,y)\left[\frac{|\delta_{h}u(x)|^{\kappa+1}+|\delta_{h}u(y)|^{\kappa+1}}{|h|^{1+\nu\kappa}}\big{(}|u_{h}(x)-u_{h}(y)|^{\frac{q-2}{2}}+|u(x)-u(y)|^{\frac{q-2}{2}}\big{)}^{2}\right.$
$\displaystyle\left.\qquad\qquad\qquad\times\big{|}\eta^{\frac{p}{2}}(x)-\eta^{\frac{p}{2}}(y)\big{|}^{2}\right]d\mu_{2}$
$\displaystyle=:-C\lambda\mathbb{I}_{4,q}.$ (4.9)

In order to estimate $\mathbb{I}_{4,q}$, we set

$\displaystyle\mathbb{\tilde{I}}_{4,q}:$
$\displaystyle=\lambda\int_{B_{R}}\int_{B_{R}}a(x,y)\frac{|\delta_{h}u(x)|^{\kappa+1}}{|h|^{1+\nu\kappa}}|u(x)-u(y)|^{q-2}\big{|}\eta^{\frac{p}{2}}(x)-\eta^{\frac{p}{2}}(y)\big{|}^{2}d\mu_{2}.$

For the case $q=2$, using the bounds on $\eta$ and $\nabla\eta$, we get

$\displaystyle\mathbb{\tilde{I}}_{4,q}$
$\displaystyle\leq\frac{C(N,s)}{h_{0}^{2}}\lambda\|a\|_{\infty}\|u\|_{L^{\infty}(B_{R+h_{0}})}\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}dx.$

For the case $q>2$, observe that $qs_{2}\leq ps_{1}+\sigma$, for
$\sigma<2(1-s_{1})$. Therefore, employing Young’s inequality and the bounds on
$u$, $a(\cdot,\cdot)$ and $\eta$, we deduce that

$\displaystyle\mathbb{\tilde{I}}_{4,q}$ $\displaystyle\leq CR^{ps_{1}+\sigma-
qs_{2}}\iint_{B_{R}\times B_{R}}\lambda
a(x,y)\frac{|u(x)-u(y)|^{q-2}}{|x-y|^{N+s_{1}p+\sigma}}\frac{|\delta_{h}u(x)|^{\kappa+1}}{|h|^{1+\nu\kappa}}\big{|}\eta^{\frac{p}{2}}(x)-\eta^{\frac{p}{2}}(y)\big{|}^{2}dxdy$
$\displaystyle\leq\frac{C\lambda\|a\|_{\infty}\|u\|_{L^{\infty}(B_{R})}^{q-p+1}}{h_{0}^{2}}\iint_{B_{R}\times
B_{R}}\frac{|u(x)-u(y)|^{p-2}}{|x-y|^{N+s_{1}p+\sigma-2}}\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}dxdy$
$\displaystyle\leq C\left[\iint_{B_{R}\times
B_{R}}\frac{|u(x)-u(y)|^{m}}{|x-y|^{N+ms_{1}\frac{(p-2-\epsilon)}{p-2}}}+\iint_{B_{R}\times
B_{R}}|x-y|^{\frac{m(2-2s_{1}-\epsilon
s_{1}-\sigma)}{m-p+2}-N}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}\right]$ $\displaystyle\leq
C[u]^{m}_{W^{\frac{s_{1}(p-2-\epsilon)}{p-2},m}(B_{R+h_{0}})}+C\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx,$

where we have used the relation $2-2s_{1}-\epsilon s_{1}-\sigma>0$, that is,
$\epsilon\in(0,\frac{2-\sigma}{s_{1}}-2)$ (thanks to $\sigma<2(1-s_{1})$) in
the last inequality. Note that the right hand side of the above expression is
exactly the same as the estimate of $\tilde{I}_{11}$, given by (3.20) in [23].
Thus, similar to [23, Equ. (3.22)], we have

$\displaystyle\mathbb{\tilde{I}}_{4,q}\leq
C\lambda\Bigg{(}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Bigg{)}.$

Consequently, from (4.1), we obtain

$\displaystyle-\mathbb{I}_{1,q}\leq c{\mathbb{I}}_{4,q}\leq
C\lambda\Bigg{(}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Bigg{)}.$ (4.10)

Estimate of $\mathbb{I}_{2,\ell}$. Noting $|u|\leq 1$ on $B_{1}\times[-1,0]$,
we have

$\displaystyle|a_{\ell}[u_{h}(x)-u_{h}(y)]^{\ell-1}[\delta_{h}u(x)]^{\kappa}|\leq
C\|a_{\ell}\|_{\infty}(1+|u_{h}(y)|^{\ell-1})|\delta_{h}u(x)|^{\kappa}.$

Therefore, on account of (4.5), we deduce that

$\displaystyle\int_{B_{\frac{R+r}{2}}}\int_{B_{R}^{c}}\frac{a_{\ell}(x+h,y+h)[u_{h}(x,t)-u_{h}(y,t)]^{\ell-1}}{|h|^{1+\nu\kappa}}[\delta_{h}u(x)]^{\kappa}\eta^{p}(x)d\mu$
$\displaystyle\leq
C\int_{B_{\frac{R+r}{2}}}\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}dx\int_{B_{R}^{c}}\frac{1+|u_{h}(y)|^{\ell-1}}{|y|^{N+s\ell}}dy\leq\frac{C}{(R-r)^{s\ell}}\Bigg{(}1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Bigg{)}.$

Consequently, we have

$\displaystyle|\mathbb{I}_{2,p}|+|\mathbb{I}_{2,q}|\leq
C(1+\lambda)\Bigg{(}1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Bigg{)}.$ (4.11)

Estimate of $\mathbb{I}_{3}$ and $\mathbb{I}_{f}$. Using (4.5) once again, we
get

$\displaystyle\mathbb{I}_{3}\leq\frac{C}{\mu}\int_{T_{0}}^{T_{1}}\Bigg{(}1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Bigg{)}dt$ (4.12)

and noting $f\in L^{\infty}(B_{R}\times[T_{0},T_{1}])$, we infer that

$\displaystyle\mathbb{I}_{f}\leq
2\|f\|_{L^{\infty}(B_{R}\times[T_{0},T_{1}])}\int_{T_{0}}^{T_{1}}\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\kappa}}{|h|^{1+\nu\kappa}}dxdt\leq
C\int_{T_{0}}^{T_{1}}\Big{(}1+\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx\Big{)}dt.$ (4.13)

Combining (4.1), (4.8), (4.10), (4.11), (4.12) and (4.13) with (4.1), we
obtain

$\displaystyle\int_{T_{0}+\mu}^{T_{1}}\bigg{[}\frac{[\delta_{h}u]^{\frac{\kappa+p-1}{p}}\eta}{|h|^{\frac{1+\nu\kappa}{p}}}\bigg{]}^{p}_{W^{s_{1},p}(B_{R})}dt+\frac{1}{\kappa+1}\int_{B_{R}}\frac{|\delta_{h}u(x,T_{1})|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}dx$
$\displaystyle\leq
C_{\lambda,f}\int_{T_{0}}^{T_{1}}\Bigg{(}\int_{B_{R}}\frac{|\delta_{h}u(x)|^{\frac{\kappa
m}{m-p+2}}}{|h|^{\frac{(1+\nu\kappa)m}{m-p+2}}}dx+\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1\Bigg{)}dt,$

where
$C_{\lambda,f}=C(N,p,s_{1},h_{0},q,m,\|a\|_{\infty})(1+\lambda+\|f\|_{L^{\infty}(B_{1}\times(-1,0])})$
is a constant and we have used the fact that $\tau\equiv 1$ on
$[T_{0}+\mu,T_{1}]$. Then, for the choice of $(1+\nu\kappa)/\kappa<1$,
following the calculation to derive [23, Equ. (3.35)], we get

$\displaystyle\int_{T_{0}+\mu}^{T_{1}}\sup_{0<|h|<h_{0}}\int_{B_{r}}\bigg{|}\frac{\delta^{2}_{h}u}{|h|^{\frac{1+\nu\kappa+s_{1}p}{\kappa-1+p}}}\bigg{|}^{\kappa-1+p}dxdt+\frac{1}{\kappa+1}\sup_{0<|h|<h_{0}}\int_{B_{R}}\frac{|\delta_{h}u(x,T_{1})|^{\kappa+1}}{|h|^{1+\nu\kappa}}\eta^{p}dx$
$\displaystyle\leq
C_{\lambda,f}\int_{T_{0}}^{T_{1}}\Bigg{(}\sup_{0<|h|<h_{0}}\int_{B_{R}}\bigg{|}\frac{\delta^{2}_{h}u}{|h|^{\frac{1+\nu\kappa}{\kappa}}}\bigg{|}^{\frac{m\kappa}{m-p+2}}dx+\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1\Bigg{)}dt.$

The above expression is same as [6, Equ. (4.11)], thus, from step 6 of the
same paper, we deduce that

$\displaystyle\int_{T_{0}+\mu}^{T_{1}}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m+1}_{L^{m+1}(B_{R-4h_{0}})}dt+\frac{1}{m+3-p}\sup_{0<|h|<h_{0}}\int_{B_{R-4h_{0}}}\frac{|\delta_{h}u(x,T_{1})|^{m+3-p}}{|h|^{(m+2-p)s_{1}}}dx$
$\displaystyle\leq
C(1+\lambda+\|f\|_{L^{\infty}(B_{1}\times(-1,0])})\int_{T_{0}}^{T_{1}}\Bigg{(}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}+1\Bigg{)}dt.$

For the case $T_{1}=0$, we use an approximation argument as in step 7 of [6,
Proposition 4.1]. $\square$

######  Theorem 4.3.

Suppose that the function $a(\cdot,\cdot)$ satisfies the assumption A.1 and
$qs_{2}\leq ps_{1}$. Let $u$ be a local weak solution to problem
(${\mathcal{P}}_{t}$) such that $u\in L^{\infty}_{{\rm loc}}(\Omega\times
I)\cap L^{\infty}_{\rm loc}(I;L^{p-1}_{ps_{1}}(\mathbb{R}^{N}))\cap
L^{q-1}_{\rm loc}(I;L^{q-1}_{qs_{2}}(\mathbb{R}^{N}))$. Then, for every
$\alpha\in(0,s_{1})$, $u\in C^{0,\alpha}_{x,{\rm loc}}(\Omega\times I)$. More
precisely, for every $\alpha\in(0,s_{1})$ and $R\in(0,1)$ satisfying
$Q_{2R,2R^{s_{1}p}}\equiv Q_{2R,2R^{s_{1}p}}(x_{0},0)\Subset\Omega\times I$,
there exists a constant $C=C(N,s_{1},p,q,\alpha,\|a\|_{\infty})>0$ such that

$\displaystyle\sup_{t\in[-\frac{R^{s_{1}p}}{2},0]}[u(\cdot,t)]_{C^{0,\alpha}(B_{\frac{R}{2}})}\leq\frac{C}{R^{\alpha}}\mathcal{K}_{R}^{(q-p)i_{\infty}+1},$

where $i_{\infty}\in\mathbb{N}$ is such that
$i_{\infty}>\frac{3\alpha+(N-2)s_{1}}{s_{1}-\alpha}$ and

$\displaystyle\mathcal{K}_{R}:=$
$\displaystyle\Big{[}\|u\|_{L^{\infty}(Q_{R,R^{s_{1}p}})}+\Big{(}\frac{1}{R^{N}}\int_{-\frac{7}{8}R^{s_{1}p}}^{0}[u(\cdot,t)]_{W^{s_{1},p}(B_{R})}^{p}dt\Big{)}^{\frac{1}{p}}+1$
$\displaystyle\ +\sum_{(\ell,s)}T^{\ell}_{\infty,\ell
s}(u;x_{0},R,-R^{s_{1}p},0)^{\frac{1}{\ell-1}}+\|f\|^{\frac{1}{p-1}}_{L^{\infty}(Q_{R,R^{s_{1}p}})}\Big{]}.$
(4.14)

Proof.  Let $\sigma\in[-R^{s_{1}p}(1-\mathcal{K}_{R}^{2-p}),0]$ and define

$\displaystyle
u_{R,\sigma}(x,t):=\frac{1}{\mathcal{K}_{R}}u(Rx,\mathcal{K}_{R}^{2-p}R^{s_{1}p}t+\sigma)\quad\mbox{for
}(x,t)\in B_{2}\times(-2,0].$

By the scaling property of the operators involved in problem
(${\mathcal{P}}_{t}$), it is easy to observe that $u_{R,\sigma}$ solves

$\displaystyle\partial_{t}u_{R,\sigma}+(-\Delta)_{p}^{s_{1}}u_{R,\sigma}+\lambda(-\Delta)_{q,\tilde{a}}^{s_{2}}u_{R,\sigma}=\tilde{f}\quad\mbox{in
}B_{2}\times(-2,0],$

where $\lambda=\mathcal{K}_{R}^{q-p}R^{s_{1}p-s_{2}q}$,
$\tilde{a}(x,y)=a(Rx,Ry)$ and
$\tilde{f}(x,t)=K_{R}^{1-p}R^{ps_{1}}f(Rx,\mathcal{K}_{R}^{2-p}R^{s_{1}p}t+\sigma)$.
Moreover, $u_{R,\sigma}$ satisfies

$\displaystyle\|u_{R,\sigma}\|_{L^{\infty}(B_{1}\times[-1,0])}\leq 1,\ \
\sum_{(\ell,s)}T^{\ell}_{\infty,\ell
s}(u_{R,\sigma};0,1,-1,0)^{\frac{1}{\ell-1}}\leq 1\quad\mbox{and
}\int_{-\frac{7}{8}}^{0}[u_{R,\sigma}]_{W^{s_{1},p}(B_{1})}^{p}dt\leq 1.$

We claim that

$\displaystyle\sup_{t\in[-1/2,0]}[u_{R,\sigma}(\cdot,t)]_{C^{0,\alpha}(B_{1/2})}\leq
C(1+\lambda+\|\tilde{f}\|_{L^{\infty}(B_{1}\times(-1,0])})^{i_{\infty}},$
(4.15)

where $i_{\infty}\in\mathbb{N}$ is such that
$i_{\infty}>\frac{3\alpha+(N-2)s_{1}}{s_{1}-\alpha}$ and $C$ is a constant
independent of $R$ and $\sigma$. Indeed, the proof follows by using
Proposition 4.2 and following the iterative scheme presented as in [6, Theorem
4.2] together with the observation that the process runs for $i_{\infty}$
steps only. Scaling back to $u$, from (4.15), we get

$\displaystyle\sup_{\sigma-\frac{1}{2}\mathcal{K}_{R}^{2-p}R^{s_{1}p}\leq
t\leq
0}[u(\cdot,t)]_{C^{0,\alpha}(B_{R/2})}\leq\frac{C}{R^{\alpha}}(1+\lambda+\mathcal{K}_{R}^{1-p}R^{ps_{1}}\|f\|_{L^{\infty}(B_{R}\times(-R^{ps_{1}},0])})^{i_{\infty}}\mathcal{K}_{R},$

where $C=C(N,s_{1},p,q,\|a\|_{\infty},\alpha)>0$ is a constant. Note that in
the above expression $C$ does not depend on $\sigma$ and noting
$\sigma\in[-R^{s_{1}p}(1-\mathcal{K}_{R}^{2-p}),0]$ with
$0<R<1\leq\mathcal{K}_{R}$, we get

$\displaystyle\sup_{-\frac{R^{s_{1}p}}{2}\leq t\leq
0}[u(\cdot,t)]_{C^{0,\alpha}(B_{R/2})}\leq\frac{C}{R^{\alpha}}\mathcal{K}_{R}^{(q-p)i_{\infty}+1}.$

This proves the theorem by recalling the definition of $\mathcal{K}_{R}$.
$\square$

######  Proposition 4.4.

Let the hypotheses of Proposition 4.2 hold. Additionally suppose that for some
$\kappa\geq 2$, $\nu<1$ and $0<h_{0}<\frac{1}{10}$, we have
$(1+\nu\kappa)/\kappa<1$ and

$\displaystyle\int_{T_{0}}^{T_{1}}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{s_{1}}}\bigg{\|}^{m}_{L^{m}(B_{R+4h_{0}})}dt<\infty,$

for $4h_{0}<R\leq 1-5h_{0}$ and $-1<T_{0}<T_{1}\leq 0$. Then, for all
$\mu\in(0,T_{1}-T_{0})$, there holds:

$\displaystyle\int_{T_{0}+\mu}^{T_{1}}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{\frac{1+s_{1}p+\nu\kappa}{\kappa-1+p}}}\bigg{\|}^{\kappa-1+p}_{L^{\kappa-1+p}(B_{R-4h_{0}})}dt+\frac{1}{\kappa+1}\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta_{h}u(\cdot,T_{1})}{|h|^{\frac{1+\nu\kappa}{\kappa+1}}}\bigg{\|}^{\kappa+1}_{L^{\kappa+1}(B_{R-4h_{0}})}$
$\displaystyle\leq
C(1+\lambda+\|f\|_{L^{\infty}(B_{1}\times(-1,0])})\int_{T_{0}}^{T_{1}}\left(\sup_{0<|h|<h_{0}}\bigg{\|}\frac{\delta^{2}_{h}u(\cdot,t)}{|h|^{\frac{1+\nu\kappa}{\kappa}}}\bigg{\|}^{\kappa}_{L^{\kappa}(B_{R+4h_{0}})}+1\right)dt,$

where $C=C(N,h_{0},p,q,\kappa,s_{1},\|a\|_{\infty})>0$ (which depends
inversely on $h_{0}$).

Proof.  For every $\epsilon>0$, by Theorem 4.3, we have $u\in
C^{0,s_{1}-\epsilon}_{x,{\rm loc}}(B_{2}\times(-2,0])$ and

$\displaystyle\sup_{t\in[T_{1},T_{1}]}[u(\cdot,t)]_{C^{0,s_{1}-\epsilon}(B_{R+h_{0}})}\leq
C(N,h_{0},p,s_{1},q,s_{2},\|a\|_{\infty}).$

Consequently, for
$\epsilon\in\Big{(}0,\min\\{\frac{2(1-s_{1})}{p-2},s_{1}\\}\Big{)}$, using the
bounds on $\eta$ and $\nabla\eta$, we obtain

$\displaystyle\int_{B_{R}}a_{\ell}\frac{|u(x)-u(y)|^{\ell-2}|\eta^{\frac{p}{2}}(x)-\eta^{\frac{p}{2}}(y)|^{2}}{|x-y|^{N+s\ell}}dy\leq
C\int_{B_{R}}|x-y|^{-N-s\ell+2+(\ell-2)(s_{1}-\epsilon)}dy\leq C,$

where we have used $2+(\ell-2)(s_{1}-\epsilon)-s\ell>0$. Therefore, with the
notation of Proposition 4.2, we have

$\displaystyle|\mathbb{I}_{5,p}(t)+\mathbb{I}_{4,q}(t)|\leq
C(1+\lambda)\int_{B_{R}}\frac{|\delta_{h}u(x,t)|^{\kappa+1}}{|h|^{1+\nu\kappa}}dx\leq
C\int_{B_{R}}\frac{|\delta_{h}u(x,t)|^{\kappa}}{|h|^{1+\nu\kappa}}dx,$

where we used the local boundedness of $u$. The other terms can be estimated
similarly. Then, proceeding as in [6, Proposition 5.1], we can complete the
proof of the proposition. $\square$

Using Proposition 4.4 and following the iteration argument as in the proof of
[6, Theorem 5.2], we can prove the following higher space regularity result.

######  Theorem 4.5.

Let the hypothesis of Theorem 4.3 be true. Then, for every
$\alpha\in(0,\Theta)$, $u\in C^{0,\alpha}_{x,{\rm loc}}(\Omega\times I)$,
where $\Theta$ is given by (2.5). More precisely, for every
$\alpha\in(0,\Theta)$ and $R\in(0,1)$ satisfying
$Q_{2R,2R^{s_{1}p}}\Subset\Omega\times I$, there exists a constant
$C=C(N,s_{1},s_{2},p,q,\alpha,\|a\|_{\infty})>0$ such that

$\displaystyle\sup_{t\in[-\frac{R^{s_{1}p}}{2},0]}[u(\cdot,t)]_{C^{0,\alpha}(B_{\frac{R}{2}})}\leq\frac{C}{R^{\alpha}}\mathcal{K}_{R}^{j_{\infty}},$
(4.16)

for some $j_{\infty}\in\mathbb{N}$ depending only on $N,p,s_{1}$ and $\alpha$,
where the constant $\mathcal{K}_{R}$ is given by (4.3).

###  4.2 Regularity in time

To move forward, for $u\in L^{1}(Q_{r,\theta}(x_{0},t_{0}))$, we set

$\displaystyle(u)_{(x_{0},t_{0}),r,\theta}=\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}u(x,t)dxdt,$

and when the center is clear, we denote it by $(u)_{r,\theta}$.

######  Proposition 4.6.

Let $u$ be a local weak solution to problem (${\mathcal{P}}_{t}$) in
$B_{2}\times(-2,0]$ with $f\in L^{\infty}(B_{1}\times(-1,0])$ such that

$\displaystyle\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-1,0])}\leq
1\quad\mbox{and }\sup_{t\in[-1/2,0]}[u(\cdot,t)]_{C^{0,\alpha}(B_{1/2})}\leq
K_{\alpha}\quad\mbox{for }\alpha\in(s_{1},\Theta),$ (4.17)

where $\Theta$ is defined in (2.5). Then, there exists a positive constant $C$
depending only on $N,s_{1},p,q,K_{\alpha},\alpha,\|a\|_{\infty}$ and
$\|f\|_{L^{\infty}(B_{1}\times(-1,0])}$ such that

$\displaystyle|u(x,t)-u(x,\tau)|\leq C|t-\tau|^{\gamma},\quad\mbox{for all
}x\in B_{1/4},\mbox{ and }t,\tau\in(-1/4,0],$

where $\gamma=\frac{\alpha}{ps_{1}-(p-2)\alpha}$. In particular, $u\in
C^{0,\gamma}_{t}(Q_{\frac{1}{4},\frac{1}{4}})$ for any $\gamma<\Gamma$, where
$\Gamma$ is given by (2.6).

Proof.  For fixed $(x_{0},t_{0})\in Q_{1/4,1/4}$, we choose $r\in(0,1/8)$ and
$\theta\in(0,1/8)$. Then, we observe that

$\displaystyle
Q_{r,\theta}(x_{0},t_{0}):=B_{r}(x_{0})\times(t_{0}-\theta,t_{0}]\subset
B_{3/8}\times(-1/2,0].$

Let $\psi\in C^{\infty}_{c}(B_{r/2}(x_{0}))$ be a non-negative cut-off
function such that

$\displaystyle\psi\equiv\|\psi\|_{L^{\infty}(B_{r/2}(x_{0}))}\ \mbox{ in
}B_{r/4}(x_{0}),\quad(\psi)_{r}=1\quad\mbox{and
}\|\nabla\psi\|_{L^{\infty}(B_{r/2}(x_{0}))}\leq C/r,$

for some $C=C(N)>0$ (note that the assumption $(\psi)_{r}$ yields
$\|\psi\|_{L^{\infty}(B_{r/2}(x_{0}))}\leq 4^{N}$). Set

$\displaystyle(u\psi)_{r}(t):=\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}u(x,t)\psi(x)dx.$

Then, on account of [6, Equ. (6.3)], we observe that

$\displaystyle\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}|u(x,t)-(u)_{r,\theta}|dxdt$
$\displaystyle\leq 2\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}|u(x,t)-(u\psi)_{r}(t)|dxdt$
$\displaystyle\quad+2\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}|(u\psi)_{r}(t)-(u\psi)_{r,\theta}|dxdt$
$\displaystyle=:\mathcal{O}_{1}+\mathcal{O}_{2}.$

Using the Poincaré inequality of [6, Lemma 6.1], we obtain

$\displaystyle\mathcal{O}_{1}$ $\displaystyle\leq
2\Bigg{(}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}|u(x,t)-(u\psi)_{r}(t)|^{p}dxdt\Bigg{)}^{1/p}$
$\displaystyle\leq
C\Bigg{(}\frac{r^{s_{1}p}}{|Q_{r,\theta}(x_{0},t_{0})|}\int_{t_{0}-\theta}^{t_{0}}\int_{B_{r}(x_{0})}\int_{B_{r}(x_{0})}\frac{|u(x,t)-u(y,t)|^{p}}{|x-y|^{N+ps_{1}}}dxdydt\Bigg{)}^{1/p},$

for some $C=C(N,p,s_{1})>0$. Note that for $x\in B_{r}(x_{0})$, we have
$B_{r}(x_{0})\subset B_{2r}(x)\subset B_{1/2}$. Thus, using (4.17), we get
(note that $\alpha>s_{1}$)

$\displaystyle\int_{B_{r}(x_{0})}\int_{B_{r}(x_{0})}\frac{|u(x,t)-u(y,t)|^{p}}{|x-y|^{N+ps_{1}}}dxdy$
$\displaystyle\leq
K_{\alpha}^{p}\int_{B_{r}(x_{0})}\left(\int_{B_{r}(x_{0})}|x-y|^{-N+p(\alpha-
s_{1})}dy\right)dx$ $\displaystyle\leq CK_{\alpha}^{p}r^{(\alpha-s_{1})p},$

where $C=C(N,s_{1},p)$. Consequently,

$\displaystyle\mathcal{O}_{1}\leq CK_{\alpha}r^{\alpha}.$

To estimate $\mathcal{O}_{2}$, we observe that

$\displaystyle\mathcal{O}_{2}\leq
2\sup_{T_{0},T_{1}\in(t_{0}-\theta,t_{0}]}|(u\psi)_{r}(T_{0})-(u\psi)_{r}(T_{1})|.$

Thus, for $T_{0},T_{1}\in(t_{0}-\theta,t_{0}]$ with $T_{0}<T_{1}$, from the
weak formulation of problem (${\mathcal{P}}_{t}$), we obtain

$\displaystyle|(u\psi)_{r}(T_{0})-(u\psi)_{r}(T_{1})|$
$\displaystyle=\frac{1}{|B_{r}(x_{0})|}\bigg{|}\int_{B_{r}(x_{0})}u(x,T_{0})\psi(x)dx-\int_{B_{r}(x_{0})}u(x,T_{1})\psi(x)dx\bigg{|}$
$\displaystyle=\frac{1}{|B_{r}(x_{0})|}\bigg{|}\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}[u(x,t)-u(y,t)]^{p-1}(\psi(x)-\psi(y))d\mu_{1}dt$
$\displaystyle\
+\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)[u(x,t)-u(y,t)]^{q-1}(\psi(x)-\psi(y))d\mu_{2}dt$
$\displaystyle\
-\int_{T_{0}}^{T_{1}}\int_{B_{r}(x_{0})}f(x,t)\psi(x)dxdt\bigg{|}$
$\displaystyle\leq\bigg{|}\sum_{(\ell,s)}\int_{T_{0}}^{T_{1}}\int_{B_{r}(x_{0})}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}a_{\ell}[u(x,t)-u(y,t)]^{\ell-1}(\psi(x)-\psi(y))d\mu
dt\bigg{|}$ $\displaystyle\
+2\bigg{|}\sum_{(\ell,s)}\int_{T_{0}}^{T_{1}}\int_{\mathbb{R}^{N}\setminus
B_{r}(x_{0})}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}a_{\ell}[u(x,t)-u(y,t)]^{\ell-1}\psi(x)d\mu
dt\bigg{|}$ $\displaystyle\
+\bigg{|}\int_{T_{0}}^{T_{1}}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}f(x,t)\psi(x)dxdt\bigg{|}$
$\displaystyle=:\mathfrak{J}_{1}+\mathfrak{J}_{2}+\mathfrak{J}_{3}.$

Using Hölder’s inequality, (4.17) and the bounds on $\psi$, we get (noting
$\alpha>\max\\{s_{1},s_{2}\\}$)

$\displaystyle\mathfrak{J}_{1}$
$\displaystyle\leq\frac{\|a_{\ell}\|_{\infty}}{|B_{r}(x_{0})|}[\psi]_{W^{s,\ell}(B_{r}(x_{0}))}\int_{t_{0}-\theta}^{t_{0}}\Bigg{(}\int_{B_{r}(x_{0})}\int_{B_{r}(x_{0})}|u(x,t)-u(y,t)|^{\ell}d\mu\Bigg{)}^{(\ell-1)/\ell}dt$
$\displaystyle\leq
CK_{\alpha}^{\ell-1}r^{-s}\int_{t_{0}-\theta}^{t_{0}}\Bigg{(}\int_{B_{r}(x_{0})}\int_{B_{r}(x_{0})}|x-y|^{\alpha\ell-
N-s\ell}dxdy\Bigg{)}^{(\ell-1)/\ell}dt$ $\displaystyle\leq
C\;K_{\alpha}^{\ell-1}\;\theta\;r^{\alpha(\ell-1)-\ell s},$

where $C=C(N,p,s_{1},\alpha,\|a\|_{\infty})>0$. To estimate
$\mathfrak{J}_{2}$, we observe that for $x\in B_{r}(x_{0})$, $t\in[-1/2,0]$
and $y\in\mathbb{R}^{N}$, we have

$\displaystyle|u(x,t)-u(y,t)|\leq C|x-y|^{\alpha}\quad\mbox{for some
}C=C(K_{\alpha},\alpha)>0.$

Indeed, for $y\in B_{1/2}$ it follows from (4.17) that if
$y\in\mathbb{R}^{N}\setminus B_{1/2}$, then

$\displaystyle|u(x,t)-u(y,t)|\leq
2\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[-1,0])}\leq
2\leq\frac{2}{8^{\alpha}}|x-y|^{\alpha}.$

Therefore (noticing that ${\rm supp}(\psi)\subset B_{r/2}(x_{0})$),

$\displaystyle\mathfrak{J}_{2}$ $\displaystyle\leq
2(T_{1}-T_{0})\|a\|_{\infty}\|\psi\|_{L^{\infty}(B_{r/2}(x_{0}))}\int_{\mathbb{R}^{N}\setminus
B_{r}(x_{0})}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r/2}(x_{0})}\frac{|u(x,t)-u(y,t)|^{\ell-1}}{|x-y|^{N+s\ell}}dxdy$
$\displaystyle\leq C\theta\int_{\mathbb{R}^{N}\setminus
B_{r}(x_{0})}\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{B_{r}(x_{0})}|x-y|^{\alpha(\ell-1)-N-s\ell}dxdy$
$\displaystyle\leq C\;\theta\;r^{\alpha(\ell-1)-\ell s},$

where $C=C(N,p,s_{1},\alpha,\|a\|_{\infty})>0$. Finally, noting the fact that
$r\in(0,1/8)$ and $\alpha(p-1)-ps_{1}<0$, we have

$\displaystyle\mathfrak{J}_{3}\leq
C\|f\|_{L^{\infty}(Q_{r,\theta}(x_{0},t_{0}))}\theta\leq
C\|f\|_{L^{\infty}(Q_{r,\theta}(x_{0},t_{0}))}\theta r^{\alpha(p-1)-ps_{1}}.$

Consequently,

$\displaystyle\mathcal{O}_{2}\leq C\;\theta\;r^{\alpha(\ell-1)-\ell s}.$

Finally, on account of the relation $\alpha(p-1)-ps_{1}\leq\alpha(q-1)-qs_{2}$
and the fact $r\in(0,1/8)$, we have

$\displaystyle\mathchoice{{\vbox{\hbox{$\textstyle-$
}}\kern-7.83337pt}}{{\vbox{\hbox{$\scriptstyle-$
}}\kern-5.90005pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.75003pt}}{{\vbox{\hbox{$\scriptscriptstyle-$
}}\kern-4.25003pt}}\\!\int_{Q_{r,\theta}(x_{0},t_{0})}|u(x,t)-(u)_{r,\theta}|dxdt\leq
Cr^{\alpha}+C\;\theta\;r^{\alpha(p-1)-ps_{1}},$

where
$C=C(N,p,s_{1},\alpha,K_{\alpha},\|a\|_{\infty},\|f\|_{L^{\infty}(Q_{r,\theta}(x_{0},t_{0}))})>0$.
Now, the rest of the proof can be completed as in the proof of [6, Proposition
6.2]. $\square$

######  Remark 2.

It is clear from the proof that when $u\in
L^{\infty}(\mathbb{R}^{N}\times[-1,0])$ is not assumed, then only the proof of
$\mathfrak{J}_{2}$ is changed. Indeed, it can be shown that

$\displaystyle\mathfrak{J}_{2}\leq C\;\theta\;r^{-\ell s}.$

Then the proof of the time regularity can be completed by proceeding similarly
and choosing $\theta$ accordingly.

Proof of Theorem 2.1: The proof of the theorem follows from Theorem 4.5 and
Proposition 4.6 after rescaling (see [6, Section 7]). Additionally, by noting
$u\in L^{\infty}(\mathbb{R}^{N}\times[\hat{t}_{0}-R^{s_{1}p},\hat{t}_{0}])$
and using Corollary 3.2, we see that

$\displaystyle\mathcal{K}_{R}\leq
C\Big{(}\|u\|_{L^{\infty}(\mathbb{R}^{N}\times[\hat{t}_{0}-R^{s_{1}p},\hat{t}_{0}])}+1+\|f\|^{\frac{1}{p-1}}_{L^{\infty}(Q_{R,R^{s_{1}p}})}\Big{)}^{q/p}.$

This completes the proof of the theorem. $\square$

##  5 A further investigation of the stationary case

In this section, we obtain higher Hölder continuity results for weak solutions
to problem (${\mathcal{S}}$) with more general modulating coefficients and
non-homogeneous terms.

###  5.1 Interior regularity

We first start with the following stability result.

######  Lemma 5.1.

Suppose that $p\geq 2$, $q<p^{*}_{s_{1}}$ and $qs_{2}\leq ps_{1}$. Let $u$ be
a weak solution to problem (${\mathcal{S}}$) in $B_{2}\equiv
B_{2}(x_{0})\Subset\Omega$, such that

$\displaystyle\|u\|_{L^{\infty}(B_{{3/2}})}+\sum_{(\ell,s)}T^{\ell}_{\ell
s}(u;x_{0},{3/2})\leq M,$

for some $M\geq 1$. Let $\delta\in(0,1)$ be such that

$\displaystyle\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta\quad\mbox{and}\quad\|f\|_{L^{\Bbbk}(B_{3/2})}\leq\delta,$

where $\tilde{a}(\cdot,\cdot)\in
L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ is a non-negative function
satisfying the assumption A.1 in $B_{1}\times B_{1}$ with

$\left\\{\begin{array}[]{rllll}(-\Delta)_{p}^{s_{1}}v+(-\Delta)_{q,\tilde{a}}^{s_{2}}v&=0\quad\mbox{in
}B_{1},\\\ v&=u\quad\mbox{in }\mathbb{R}^{N}\setminus
B_{1}.\end{array}\right.$

Then, there exists
$\varepsilon(\delta)=\varepsilon(\delta,N,p,q,s_{1},\Bbbk,M,\|a\|_{\infty})>0$
such that

$\displaystyle\|u-v\|_{L^{\infty}(B_{3/8})}<\varepsilon(\delta)$

and $\varepsilon(\delta)$ converges to $0$ as $\delta\to 0$.

Proof.  From Lemma 3.4, we see that $u$ and $v$ are locally bounded. Then, the
triangle inequality and the relation
$\|a\|_{\infty}-\delta\leq\|\tilde{a}\|_{\infty}\leq\|a\|_{\infty}+\delta$
prove the existence of $\varepsilon(\delta)$. To prove the convergence of
$\varepsilon(\delta)$ with respect to $\delta$, we proceed by the method of
contradiction. Suppose that there exist sequences $\\{a_{m}\\}$,
$\\{\tilde{a}_{m}\\}$, $\\{u_{m}\\}$, $\\{f_{m}\\}$ and $\\{v_{m}\\}$ with the
following properties:

$\displaystyle(-\Delta)_{p}^{s_{1}}u_{m}+(-\Delta)_{q,a_{m}}^{s_{2}}u_{m}=f_{m}\quad\mbox{in
}B_{2},$ (5.1)
$\displaystyle\|u_{m}\|_{L^{\infty}(B_{{3/2}})}+\sum_{(\ell,s)}T^{\ell}_{\ell
s}(u_{m};x_{0},{3/2})\leq M,$ (5.2)
$\displaystyle\|a_{m}-\tilde{a}_{m}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq
1/m,\quad\|f_{m}\|_{L^{\Bbbk}(B_{3/2})}\leq 1/m$ (5.3)

and

$\begin{array}[]{rllll}(-\Delta)_{p}^{s_{1}}v_{m}+(-\Delta)_{q,\tilde{a}_{m}}^{s_{2}}v_{m}=0\quad\mbox{in
}B_{1};\quad v_{m}=u_{m}\quad\mbox{ in }\mathbb{R}^{N}\setminus
B_{1},\end{array}$ (5.4)

but

$\displaystyle\liminf_{m\to\infty}\|u_{m}-v_{m}\|_{L^{\infty}(B_{3/8})}>0.$
(5.5)

Note that on account of (5.2) and $qs_{2}\leq ps_{1}$, we have

$\displaystyle\int_{B_{\frac{3}{2}}}\int_{B_{\frac{3}{2}}}\tilde{a}_{m}\frac{|u_{m}(x)-u_{m}(y)|^{q}}{|x-y|^{N+qs_{2}}}dxdy\leq
c\|u_{m}\|_{L^{\infty}(B_{\frac{3}{2}})}^{q-p}\int_{B_{\frac{3}{2}}}\int_{B_{\frac{3}{2}}}\frac{|u_{m}(x)-u_{m}(y)|^{p}}{|x-y|^{N+ps_{1}}}dxdy<\infty.$

Consequently $u_{m}\in X_{u_{m},\tilde{a}_{m}}(B_{1},B_{3/2})$ and
$w_{m}:=u_{m}-v_{m}\in X_{0,\tilde{a}_{m}}(B_{1},B_{3/2})$. Additionally,
since $v_{m}\in\mathcal{W}_{\tilde{a}_{m}}(B_{3/2})$ and $u_{m}\in
L^{\infty}(B_{3/2})$, proceeding similar to the proof of Lemma 3.4, we can
show that $v_{m}\in L^{\infty}(B_{1})$ (however the bound may depend on $m$).
As before, we see that $w_{m}\in\mathcal{W}_{a_{m}}(B_{3/2})$ with support
contained in $B_{1}$. Therefore taking $w_{m}=u_{m}-v_{m}$ as a test function
in the weak formulation of (5.1) and (5.4), and using the monotonocity of the
fractional operators (see (2.4)), we deduce that

$\displaystyle
c_{p}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}|w_{m}(x)-w_{m}(y)|^{p}d\mu_{1}+c_{q}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\tilde{a}_{m}(x,y)|w_{m}(x)-w_{m}(y)|^{q}d\mu_{2}$
$\displaystyle\leq\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}([u_{m}(x)-u_{m}(y)]^{p-1}-[v_{m}(x)-v_{m}(y)]^{p-1})(w_{m}(x)-w_{m}(y))d\mu_{1}$
$\displaystyle\
+\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\tilde{a}_{m}(x,y)\Big{(}[u_{m}(x)-u_{m}(y)]^{q-1}-[v_{m}(x)-v_{m}(y)]^{q-1}\Big{)}(w_{m}(x)-w_{m}(y))d\mu_{2}$
$\displaystyle=\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}(\tilde{a}_{m}(x,y)-a_{m}(x,y))[u_{m}(x)-u_{m}(y)]^{q-1}(w_{m}(x)-w_{m}(y))d\mu_{2}+\int_{B_{2}}f_{m}w_{m}dx$
$\displaystyle\leq\|a_{m}-\tilde{a}_{m}\|_{\infty}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}|u_{m}(x)-u_{m}(y)|^{q-1}|w_{m}(x)-w_{m}(y)|d\mu_{2}+\int_{B_{2}}|f_{m}(x)w_{m}(x)|dx.$
(5.6)

To estimate the above integrals, we set

$\displaystyle\mathbb{J}_{1}:=\int_{B_{\frac{6}{5}}}\int_{B_{\frac{6}{5}}}|u_{m}(x)-u_{m}(y)|^{q-1}|w_{m}(x)-w_{m}(y)|d\mu_{2},$
$\displaystyle\mathbb{J}_{2}:=\int_{B_{\frac{6}{5}}}\int_{\mathbb{R}^{N}\setminus
B_{\frac{6}{5}}}|u_{m}(x)-u_{m}(y)|^{q-1}|w_{m}(x)|d\mu_{2}\quad\mbox{and}\quad\mathbb{J}_{3}:=\int_{B_{2}}f_{m}(x)w_{m}(x)dx.$

Using (5.2), the fact $qs_{2}\leq ps_{1}$ and Hölder’s inequality, we get

$\displaystyle\mathbb{J}_{1}$ $\displaystyle\leq
C_{q}\|u_{m}\|_{L^{\infty}(B_{\frac{6}{5}})}^{q-p}\int_{B_{\frac{6}{5}}}\int_{B_{\frac{6}{5}}}|u_{m}(x)-u_{m}(y)|^{p-1}|w_{m}(x)-w_{m}(y)|\frac{dxdy}{|x-y|^{N+s_{1}p}}$
$\displaystyle\leq
CM^{q-p}[u_{m}]^{p-1}_{W^{s_{1},p}(B_{\frac{6}{5}})}\Bigg{(}\int_{B_{\frac{6}{5}}}\int_{B_{\frac{6}{5}}}\frac{|w_{m}(x)-w_{m}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy\Bigg{)}^{\frac{1}{p}}.$

Moreover, using Caccioppoli type relation of Lemma 3.3, for the choice
$\psi\in C^{\infty}_{c}(B_{7/5})$ satisfying $0\leq\psi\leq 1$ and $\psi\equiv
1$ in $B_{6/5}$, we see that

$\displaystyle[u_{m}\psi^{\frac{q}{p}}]^{p}_{W^{s_{1},p}(B_{\frac{7}{5}})}$
$\displaystyle\leq
C\Bigg{[}\int_{B_{\frac{3}{2}}}\int_{B_{\frac{3}{2}}}\frac{|\psi(x)-\psi(y)|^{p}}{|x-y|^{N+ps_{1}}}dxdy+\int_{\mathbb{R}^{N}\setminus
B_{\frac{3}{2}}}\frac{|u_{m}(y)|^{p-1}+1}{|x_{0}-y|^{N+ps_{1}}}dy\Bigg{]}+C\Big{(}\int_{B_{\frac{7}{5}}}|f_{m}|^{\Bbbk}\Big{)}^{\frac{p}{(p-1)\Bbbk}}$
$\displaystyle+C\Bigg{[}\int_{B_{\frac{3}{2}}}\int_{B_{\frac{3}{2}}}a_{m}(x,y)\frac{|\psi(x)-\psi(y)|^{q}}{|x-y|^{N+qs_{2}}}dxdy+\int_{\mathbb{R}^{N}\setminus
B_{\frac{3}{2}}}a_{m,2}(y)\frac{|u_{m}(y)|^{q-1}+1}{|x_{0}-y|^{N+qs_{2}}}dy\Bigg{]}$
$\displaystyle\leq CM(1+\|a\|_{\infty}),$

where we have set $a_{m,2}(y):=\sup_{x\in B_{2}}a_{m}(x,y)$ and used (5.2)
along with Hölder’s and Young’s inequalities for the term involving $f_{m}$
(see the estimate of $\mathbb{J}_{3}$ for details). This yields that
$[u_{m}]_{W^{s_{1},p}(B_{6/5})}$ is bounded independent of $m$. Consequently,

$\displaystyle\mathbb{J}_{1}\leq
C\Bigg{(}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{|w_{m}(x)-w_{m}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy\Bigg{)}^{\frac{1}{p}}.$
(5.7)

Furthermore, since $w_{m}=0$ in $\mathbb{R}^{N}\setminus B_{1}$, (5.2) and the
fractional Poincaré inequality imply

$\displaystyle\mathbb{J}_{2}$ $\displaystyle\leq
C\Bigg{[}\int_{\mathbb{R}^{N}\setminus
B_{\frac{6}{5}}}\frac{dy}{|x-y|^{N+qs_{2}}}+\int_{\mathbb{R}^{N}\setminus
B_{\frac{6}{5}}}\frac{|u_{m}(y)|^{q-1}}{|x-y|^{N+qs_{2}}}dy\Bigg{]}\int_{B_{1}}|w_{m}(x)|dx$
$\displaystyle\leq
CM\Bigg{(}\int_{B_{1}}|w_{m}(x)|^{p}dx\Bigg{)}^{\frac{1}{p}}\leq
C\Bigg{(}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{|w_{m}(x)-w_{m}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy\Bigg{)}^{\frac{1}{p}}.$
(5.8)

Similarly, noting $(p^{*}_{s_{1}})^{\prime}<\Bbbk$ and using (5.3) together
with the Sobolev embedding, we have

$\displaystyle\mathbb{J}_{3}\leq\int_{B_{\frac{3}{2}R}}|f_{m}w_{m}|dx$
$\displaystyle\leq
C\Big{(}\int_{B_{R}}|f_{m}(x)|^{(p^{*}_{s_{1}})^{\prime}}dx\Big{)}^{\frac{1}{(p^{*}_{s_{1}})^{\prime}}}\Big{(}\int_{B_{R}}|w_{m}(x)|^{p^{*}_{s_{1}}}dx\Big{)}^{\frac{1}{p^{*}_{s_{1}}}}$
$\displaystyle\leq\frac{C}{m}\Bigg{(}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{|w_{m}(x)-w_{m}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy\Bigg{)}^{\frac{1}{p}}.$
(5.9)

Now combining (5.7), (5.1), (5.1) and (5.3) with (5.1), we obtain

$\displaystyle\Bigg{(}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{|w_{m}(x)-w_{m}(y)|^{p}}{|x-y|^{N+s_{1}p}}dxdy\Bigg{)}^{\frac{p-1}{p}}\leq
C\frac{1}{m}.$

Then, by the fractional Sobolev inequality, we have

$\displaystyle\int_{B_{1}}|w_{m}(x)|^{p^{*}_{s_{1}}}dx\to 0\quad\mbox{as
}m\to\infty.$ (5.10)

Invoking Proposition 3.6 and (5.2), we get that the sequence $\\{u_{m}\\}$ is
uniformly bounded in $C^{0,\alpha_{1}}$ in $B_{1/2}$, for some
$\alpha_{1}\in(0,s_{1})$. We next note that on account of the local
boundedness result of Lemma 3.4 and Caccioppoli type result of Lemma 3.3 for
$v_{m}$ and estimating the tail term as in the proof of Proposition 3.6 with
the help of (5.10) and (5.2), there holds

$\displaystyle\|v_{m}\|_{L^{\infty}(B_{3/4})}+[v_{m}]_{W^{s_{1},p}(B_{3/4})}+\sum_{(\ell,s)}T^{\ell}_{\ell
s}(v_{m};0,\frac{3}{4})^{\frac{1}{\ell-1}}$ $\displaystyle\leq
C\Big{(}\|v_{m}\|^{\frac{q}{p}}_{L^{\vartheta}(B_{\frac{3}{4}})}+\|w_{m}\|^{\frac{q}{p}}_{L^{\vartheta}(B_{\frac{3}{4}})}+\sum_{(\ell,s)}T^{\ell}_{\ell
s}(v_{m};0,\frac{3}{4})^{\frac{1}{\ell-1}}+1+\|f\|_{L^{\Bbbk}(B_{1})}^{\frac{1}{p-1}}\Big{)}\leq
M_{1},$

where $\vartheta<p^{*}_{s_{1}}$ is given by Lemma 3.4 and $C,M_{1}>0$ are
independent of $m$ and are non-decreasing functions of $\|a\|_{\infty}$.
Consequently, from the stationary counterpart of Theorem 4.5 and (4.16), we
have

$\displaystyle\|v_{m}\|_{C^{0,\alpha_{1}}(B_{3/8})}\leq CM_{1},$ (5.11)

where $CM_{1}$ is independent of $m$ and is a non-decreasing function of
$\|a\|_{\infty}$. Thus, we get that the sequence $\\{v_{m}\\}$ is also
uniformly bounded in $C^{0,\alpha_{1}}$ norm in $B_{3/8}$. Therefore, by the
Arzela-Ascoli theorem, up to a subsequence, $\\{u_{m}-v_{m}\\}$ converges
uniformly to some function $w$ in $B_{3/8}$. This together with (5.10) implies
that $w\equiv 0$, that is, (up to a subsequence),

$\displaystyle\lim_{m\to\infty}\|u_{m}-v_{m}\|_{L^{\infty}(B_{3/8})}=0,$

which contradicts (5.5). This completes the proof of the lemma. $\square$

######  Proposition 5.2.

Suppose that $p\geq 2$, $q<p^{*}_{s_{1}}$ and $qs_{2}\leq ps_{1}$. Let $u$ be
a local weak solution to problem (${\mathcal{S}}$) in $B_{2}\Subset\Omega$,
such that

$\displaystyle\|u\|_{L^{\infty}(B_{1})}+\sum_{(\ell,s)}\int_{\mathbb{R}^{N}\setminus
B_{1}}\frac{|u(x)|^{\ell-1}}{|x|^{N+s\ell}}dx\leq M,$

with $M\geq 1$. For any $\epsilon>0$, there exists
$\delta_{1}=\delta_{1}(\epsilon,N,p,s_{1},q,\Bbbk,M,\|a\|_{\infty})>0$ such
that if

$\displaystyle\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta_{1}\quad\mbox{and}\quad\|f\|_{L^{\Bbbk}(B_{1})}\leq\delta_{1},$
(5.12)

for some $0\leq\tilde{a}(\cdot,\cdot)\in
L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ satisfying the assumption A.1,
then $u\in C^{0,{\Theta_{1}}-\epsilon}_{\rm loc}(\Omega)$, where
${\Theta_{1}}$ is given by (2.7).

Proof.  We provide the proof in two steps.
Step 1: Regularity at the origin. In this step we will prove that for every
$\epsilon>0$ and $r\in(0,1/2)$, there exists $\delta>0$ such that for the data
of the proposition,

$\displaystyle\sup_{x\in B_{r}}|u(x)-u(0)|\leq
C_{1}r^{{\Theta_{1}}-\epsilon},$

where $C_{1}=C_{1}(N,p,s_{1},s_{2},q,\|a\|_{\infty},M,\epsilon)>0$ is a
constant. Without loss of generality, we may assume that $u(0)=0$ (otherwise,
we consider $\bar{u}=u-u(0)$ as the weak solution). For fixed
$\epsilon\in(0,{\Theta_{1}})$, it suffices to show that there exist $\rho<1/8$
and $\delta>0$ such that if $f$, $a$, $\tilde{a}$ and $u$ are as in the
proposition, then for any $n\in\mathbb{N}\cup\\{0\\}$, we have

$\displaystyle\sup_{x\in B_{\rho^{n}}}|u(x)|\leq
M\rho^{n({\Theta_{1}}-\epsilon)}\quad\mbox{and}\quad\sum_{(\ell,s)}\int_{\mathbb{R}^{N}\setminus
B_{1}}\bigg{|}\frac{u(\rho^{n}x)}{\rho^{n({\Theta_{1}}-\epsilon)}}\bigg{|}^{\ell-1}\frac{dx}{|x|^{N+s\ell}}\leq
M.$ (5.13)

The claim is true for $n=0$ by the assumption. Let us assume that (5.13) holds
for some $n\in\mathbb{N}$, then we need to prove it for $n+1$ provided (5.12)
holds for $\delta>0$ small. For $\rho\in(0,1)$, we define

$\displaystyle
w_{n}(x):=\frac{u(\rho^{n}x)}{\rho^{n({\Theta_{1}}-\epsilon)}}.$

By the scaling property of the fractional operators, it is easy to verify that
$w_{n}$ satisfies the following equation:

$\displaystyle(-\Delta)_{p}^{s_{1}}w_{n}+(-\Delta)_{q,a_{n}}^{s_{2}}w_{n}=f_{n},$

where

$\displaystyle
a_{n}(x,y)=\rho^{n(ps_{1}-qs_{2}+(q-p)({\Theta_{1}}-\epsilon))}a(\rho^{n}x,\rho^{n}y)\quad\mbox{and}\quad
f_{n}(x)=\rho^{n(ps_{1}-({\Theta_{1}}-\epsilon)(p-1))}f(\rho^{n}x).$

Moreover, by the induction hypothesis and the fact $\rho<1$, we have

$\displaystyle\|w_{n}\|_{L^{\infty}(B_{1})}+\sum_{(\ell,s)}\int_{\mathbb{R}^{N}\setminus
B_{1}}\frac{|w_{n}(x)|^{\ell-1}}{|x|^{N+s\ell}}dx\leq M$ (5.14)

and since $\Bbbk>N/(ps_{1})$, using (2.7), we also have

$\displaystyle\|f_{n}\|_{L^{\Bbbk}(B_{1})}\leq\rho^{n(ps_{1}-({\Theta_{1}}-\epsilon)(p-1)-N/\Bbbk)}\|f\|_{L^{\Bbbk}(B_{\rho^{n}})}\leq\tilde{\delta}.$
(5.15)

Let $v_{n}$ be the solution to the following problem:

$(-\Delta)_{p}^{s_{1}}v_{n}+(-\Delta)_{q,\tilde{a}_{n}}^{s_{2}}v_{n}=0\quad\mbox{in
}B_{2/3};\quad v_{n}=w_{n}\quad\mbox{in }\mathbb{R}^{N}\setminus B_{2/3},$

where we have set

$\displaystyle\tilde{a}_{n}(x,y)=\rho^{n(ps_{1}-qs_{2}+(q-p)({\Theta_{1}}-\epsilon))}\tilde{a}(\rho^{n}x,\rho^{n}y).$

On account of (5.12) and the assumption $qs_{2}\leq ps_{1}$, we have

$\displaystyle\|a_{n}-\tilde{a}_{n}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta.$
(5.16)

Therefore, Lemma 5.1 yields

$\displaystyle\|v_{n}-w_{n}\|_{L^{\infty}(B_{1/4})}\leq\varepsilon(\delta),$
(5.17)

where $\varepsilon(\delta)$ is independent of $n$, provided (5.12) holds for
small $\delta>0$. This together with the fact $w_{n}(0)=0$, implies that

$\displaystyle|w_{n}(x)|$
$\displaystyle\leq|w_{n}(x)-v_{n}(x)|+|v_{n}(0)-w_{n}(0)|+|v_{n}(x)-v_{n}(0)|$
$\displaystyle\leq
2\varepsilon(\delta)+[v_{n}]_{C^{0,{\Theta_{1}}-\epsilon/2}(B_{1/8})}|x|^{{\Theta_{1}}-\epsilon/2},\quad\mbox{for
all }x\in B_{1/8},$ (5.18)

where we have used the $C^{0,{\Theta_{1}}-\epsilon/2}$ regularity of $v_{n}$
from Theorem 4.5. Next, from the stationary counterpart of (4.16), we have

$\displaystyle[v_{n}]_{C^{0,\Theta_{1}-\frac{\epsilon}{2}}(B_{\frac{1}{8}})}\leq
C_{4}\Big{(}\|v_{n}\|_{L^{\infty}(B_{\frac{1}{4}})}+[v_{n}]_{W^{s_{1},p}(B_{\frac{1}{4}})}+1+\sum_{(\ell,s)}T^{\ell}_{\ell
s}(v_{n};0,\frac{1}{4})^{\frac{1}{\ell-1}}\Big{)}^{j_{\infty}},$

where $C_{4}=C_{4}(N,p,s_{1},q,\epsilon,M,\|a\|_{\infty})>0$ is a constant.
Using (5.14) and (5.17), we obtain

$\displaystyle\|v_{n}\|_{L^{\infty}(B_{1/4})}\leq\|v_{n}-w_{n}\|_{L^{\infty}(B_{1/4})}+\|w_{n}\|_{L^{\infty}(B_{1/4})}\leq\varepsilon+2M\leq
3M,$

provided (5.15) holds. For the tail terms, noting $q\leq p^{*}$, similar to
(5.10), we have

$\displaystyle\|v_{n}-w_{n}\|_{L^{\ell}(B_{1/2})}\leq C_{5},$

where $C_{5}>0$ is independent of $n$. Therefore, on account of (5.14) and
Hölder’s inequality, we deduce that

$\displaystyle\int_{B_{1/4}^{c}}\frac{|v_{n}(y)|^{\ell-1}}{|y|^{N+s\ell}}dy$
$\displaystyle\leq
C_{6}\Big{(}\int_{B_{1/2}^{c}}\frac{|w_{n}(y)|^{\ell-1}}{|y|^{N+s\ell}}dy+\|w_{n}\|^{\ell-1}_{L^{\infty}(B_{1/2})}+\int_{B_{1/2}}|v_{n}(x)-w_{n}(x)|^{\ell-1}dx\Big{)}$
$\displaystyle\leq C_{7}.$

Similar to Lemma 5.1, using the Caccioppoli inequality of Lemma 3.3, we can
show that the Sobolev seminorm of $v_{n}$ is bounded. Consequently, from
(5.1), we get

$\displaystyle|w_{n}(x)|\leq
2\varepsilon(\delta)+C_{9}|x|^{{\Theta_{1}}-\epsilon/2}\quad\mbox{in
}B_{1/9},$ (5.19)

where $C_{9}$ is a constant independent of $n$. Next, set

$\displaystyle w_{n+1}(x):=\frac{w_{n}(\rho
x)}{\rho^{{\Theta_{1}}-\epsilon}}=\frac{u(\rho^{n+1}x)}{\rho^{(n+1)({\Theta_{1}}-\epsilon)}}.$

By choosing $\delta>0$ small enough such that
$2\varepsilon(\delta)<\rho^{\Theta_{1}}$, from (5.19), we obtain

$\displaystyle|w_{n+1}(x)|\leq
2\varepsilon\rho^{\epsilon-{\Theta_{1}}}+C_{9}\rho^{\epsilon-{\Theta_{1}}}|\rho
x|^{{\Theta_{1}}-\epsilon/2}\leq(1+C_{9}|x|^{{\Theta_{1}}-\epsilon/2})\rho^{\epsilon/2}\quad\mbox{in
}B_{1/(8\rho)}.$

In particular, choosing $\rho>0$ small enough, we get
$\|w_{n+1}\|_{L^{\infty}(B_{1/(8\rho)})}\leq 1$. This on using the definition
of $w_{n+1}$ and the fact that $8\rho<1$ proves the first part of (5.13). For
the second part, by noting the bound on ${\Theta_{1}}$ and proceeding similar
to [5, Proposition 6.2, pp. 838-839], we can show that

$\displaystyle\int_{\mathbb{R}^{N}\setminus
B_{1}}\frac{|w_{n+1}(x)|^{\ell-1}}{|x|^{N+s\ell}}dx\leq\Big{(}\frac{C_{2}}{\epsilon(\ell-1)}+C_{3}+1\Big{)}M\rho^{\varepsilon(\ell-1)/2}\leq
M,$

for appropriate choice of $\rho>0$. This completes the iteration process.
Here, we stress that the constant $\delta>0$, small enough, is chosen such
that $2\varepsilon<\rho^{\Theta_{1}}$ and $\rho$ is also chosen small enough
depending only on $N,p,s_{1},s_{2},q,\Bbbk,M,\epsilon$ and the constant $S$
satisfying $\|a\|_{\infty}\leq S$.
Step 2: Regularity in the ball. We now show the Hölder continuity of $u$ in
the whole ball $B_{1/8}$. Let $z\in B_{1/2}$ and define

$\displaystyle v(x):=L^{-\frac{1}{p-1}}u({x/2}+z),\quad
x\in\mathbb{R}^{N}\quad\mbox{with }L=2^{N+2}(1+|B_{1}|).$

Observe that $v$ is a weak solution of

$\displaystyle(-\Delta)_{p}^{s_{1}}v+(-\Delta)_{q,\hat{a}}^{s_{2}}v=\widehat{f}\quad\mbox{in
}B_{1},$

where

$\displaystyle\widehat{a}(x,y)=L^{\frac{q-p}{p-1}}2^{qs_{2}-ps_{1}}a({x/2}+z,{y/2}+z)\quad\mbox{and
}\widehat{f}(x)=2^{-ps_{1}}L^{-1}f({x/2}+z).$

For $0<\epsilon<{\Theta_{1}}$, let $\delta$ be as obtained in Step 1 for the
choice $\widehat{a}$. Assume further that

$\displaystyle\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq{2^{ps_{1}-qs_{2}}}{L^{\frac{p-q}{p-1}}}\delta.$

Setting
$\widehat{\tilde{a}}(x,y)=L^{\frac{q-p}{p-1}}2^{qs_{2}-ps_{1}}\tilde{a}({x/2}+z,{y/2}+z)$,
we have

$\displaystyle\|\widehat{a}-\widehat{\tilde{a}}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq
L^{\frac{q-p}{p-1}}2^{qs_{2}-ps_{1}}\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta\quad\mbox{and}$
$\displaystyle\|\widehat{f}\|_{L^{\Bbbk}(B_{1})}\leq\frac{2^{N/\Bbbk-
ps_{1}}}{L}\|f\|_{L^{\Bbbk}(B_{1})}\leq\delta.$

By construction, we have

$\displaystyle\|v\|_{L^{\infty}(B_{1})}\leq M\quad\mbox{and}$
$\displaystyle\int_{\mathbb{R}^{N}\setminus
B_{1}}\frac{|v(x)|^{\ell-1}}{|x|^{N+s\ell}}dx\leq\frac{1}{L^{\frac{\ell-1}{p-1}}}\int_{\mathbb{R}^{N}\setminus
B_{1}}\frac{|u(x)|^{\ell-1}}{|x|^{N+s\ell}}dx+\frac{2^{N}|B_{1}|}{L^{\frac{\ell-1}{p-1}}}\|u\|_{L^{\infty}(B_{1})}^{\ell-1}\leq
M.$

Applying Step 1 to $v$, we get

$\displaystyle\sup_{x\in B_{r}}|v(x)-v(0)|\leq
C_{1}r^{{\Theta_{1}}-\epsilon},\quad r\in(0,1/2).$

This implies that, for any $z\in B_{1/2}$,

$\displaystyle\sup_{x\in B_{r}}|u(x)-u(z)|\leq
C_{1}L^{\frac{1}{p-1}}r^{{\Theta_{1}}-\epsilon},\quad r\in(0,1/4).$ (5.20)

Next, we fix $x,y\in B_{1/8}$ and set $|x-y|=r$. We observe that $r<1/4$, then
setting $z=(x+y)/2$, we get from (5.20) that

$\displaystyle|u(x)-u(y)|\leq|u(x)-u(z)|+|u(y)-u(z)|\leq 2\sup_{\xi\in
B_{r}(z)}|u(\xi)-u(z)|\leq 2CL^{\frac{1}{p-1}}|x-y|^{{\Theta_{1}}-\epsilon}.$

This proves the proposition for the choice
$\delta_{1}=\min\Big{\\{}\delta,2^{ps_{1}-qs_{2}}(2^{N+2}(1+|B_{1}|))^{\frac{p-q}{p-1}}\delta\Big{\\}}$.
$\square$

Proof of Theorem 2.2: Let $R_{0}\in(0,1)$ be fixed such that $B_{R_{0}}\equiv
B_{R_{0}}(x_{0})\Subset\Omega$. For $\alpha\in(0,{\Theta_{1}})$ and
$\epsilon=\Theta_{1}-\alpha$, let
$\delta_{1}=\delta_{1}(\epsilon,N,p,s_{1},s_{2},q,\|a\|_{\infty})>0$ be given
by Proposition 5.2. We set

$\displaystyle\mathcal{M}=\mathcal{M}(R_{0}):=\mathcal{C}\|u\|_{L^{\infty}(B_{R_{0}})}+1+\Big{(}\frac{R_{0}^{s_{1}p-\frac{N}{\Bbbk}}\|f\|_{L^{\Bbbk}(B_{R_{0}})}}{\delta_{1}}\Big{)}^{\frac{1}{p-1}}+\sum_{(\ell,s)}2^{\frac{N+s\ell}{\ell-1}}T^{\ell}_{\ell
s}(u;x_{0},R_{0})^{\frac{1}{\ell-1}},$

where
$\mathcal{C}:=\max\Big{\\{}1,\Big{(}\frac{N|B_{1}|}{ps_{1}}\Big{)}^{1/(p-1)},\Big{(}\frac{N|B_{1}|}{qs_{2}}\Big{)}^{1/(q-1)}\Big{\\}}$.
We now choose $R\in(0,R_{0}/2)$ such that

$\displaystyle\mathcal{M}^{q-p}R^{ps_{1}-qs_{2}}\leq 1,$ (5.21)

(thanks to the condition $ps_{1}>qs_{2}$). By the assumption A.2, there exist
$\varrho,h_{\delta_{1}}>0$ such that

$\displaystyle\sup_{\begin{subarray}{c}x,y\in B_{R}(x_{0})\\\
|x-y|\leq\varrho\end{subarray}}|a(x+h,y+h)-a(x,y)|\leq\delta_{1}\quad\mbox{for
all }h\in B_{h_{\delta_{1}}}(0).$ (5.22)

Fix $z\in B_{R/2}(x_{0})$ and $r_{z}\in(0,1)$ small enough such that
$r_{z}\leq\min\\{R/2,\varrho/2,h_{\delta_{1}}\\}$ and $B_{r_{z}}(z)\subset
B_{R/2}(x_{0})$. Then, for all $x,y\in B_{r_{z}}(z)$, we see that $(z-y)$ and
$(z-x)\in B_{h_{\delta_{1}}}(0)$. Therefore,

$\displaystyle\sup_{x,y\in
B_{r_{z}}(z)}|a(x-y+z,z)-a(x,y)|\leq\delta_{1}\quad\mbox{and }\sup_{x,y\in
B_{r_{z}}(z)}|a(z,y-x+z)-a(x,y)|\leq\delta_{1}.$

Consequently, using the symmetry of $a(\cdot,\cdot)$, we see that

$\displaystyle a_{z}(x,y):=\frac{1}{2}\big{(}a(x-y+z,z)+a(y-x+z,z)\big{)}$

is translation invariant in $B_{r_{z}}(z)\times B_{r_{z}}(z)$ and satisfies

$\displaystyle\|a-a_{z}\|_{L^{\infty}(B_{r_{z}}(z)\times
B_{r_{z}}(z))}\leq\delta_{1}.$

Defining

$\displaystyle\tilde{a}(x,y):=\begin{cases}a_{z}(x,y)&\mbox{if }(x,y)\in
B_{r_{z}}(z)\times B_{r_{z}}(z),\\\ a(x,y)&\mbox{if }(x,y)\notin
B_{r_{z}}(z)\times B_{r_{z}}(z),\end{cases}$

it follows that $\tilde{a}(\cdot,\cdot)$ satisfies the assumption A.1 in
$B_{r_{z}}(z)\times B_{r_{z}}(z)$ and

$\displaystyle\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta_{1}.$
(5.23)

We set

$\displaystyle u_{z}(x):=\frac{1}{\mathcal{M}}u(r_{z}x+z)\quad\mbox{for }x\in
B_{2}.$

By the scaling property of the fractional operators, it is easy to observe
that $u_{z}$ solves

$\displaystyle(-\Delta)_{p}^{s_{1}}u_{z}+(-\Delta)_{q,a_{1}}^{s_{2}}u_{z}=f_{z}\quad\mbox{in
}B_{2},$

where $a_{1}(x,y)=\mathcal{M}^{q-p}r_{z}^{s_{1}p-s_{2}q}a(r_{z}x+z,r_{z}y+z)$
and $f_{z}(x):=\mathcal{M}^{1-p}r_{z}^{ps_{1}}f(r_{z}x+z)$. Moreover, since
$z\in B_{R/2}(x_{0})$, we see that

$\displaystyle|y-z|\geq|y-x_{0}|\Big{(}1-\frac{|z-x_{0}|}{|y-x_{0}|}\Big{)}\geq\frac{|y-x_{0}|}{2}\quad\mbox{for
all }y\in B_{R_{0}/2}(x_{0})^{c}.$

Consequently,

$\displaystyle T_{s\ell}^{\ell}(u;z,r_{z})$ $\displaystyle\leq r_{z}^{\ell
s}2^{N+\ell
s}\int_{B_{R_{0}/2}(x_{0})^{c}}\frac{|u(y)|^{\ell-1}}{|y-x_{0}|^{N+\ell
s}}dy+\|u\|_{L^{\infty}(B_{R_{0}/2}(x_{0}))}^{\ell-1}r_{z}^{\ell
s}\int_{B_{r_{z}}(z)^{c}}\frac{dy}{|y-z|^{N+\ell s}}$ $\displaystyle\leq
2^{N+\ell s}T_{s\ell}^{\ell}(u;x_{0},R_{0}/2)+\Big{(}\frac{N|B_{1}|}{\ell
s}\Big{)}\|u\|_{L^{\infty}(B_{R_{0}/2}(x_{0}))}^{\ell-1}.$ (5.24)

Recalling the definition of $\mathcal{M}$ and using (5.1), we have

$\displaystyle\|u_{z}\|_{L^{\infty}(B_{1})}\leq
1,\quad\sum_{(\ell,s)}T^{\ell}_{\ell s}(u_{z};0,1)^{1/(\ell-1)}\leq
1\quad\mbox{and}\quad\|f_{z}\|_{L^{\Bbbk}(B_{1})}\leq\delta_{1}.$

Furthermore, for
$\tilde{a}_{1}(x,y):=\mathcal{M}^{q-p}r_{z}^{s_{1}p-s_{2}q}\tilde{a}(r_{z}x+z,r_{z}y+z)$,
there holds

$\displaystyle\|a_{1}-\tilde{a}_{1}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}=r_{z}^{s_{1}p-s_{2}q}\mathcal{M}^{q-p}\|a-\tilde{a}\|_{L^{\infty}(\mathbb{R}^{N}\times\mathbb{R}^{N})}\leq\delta_{1},$

where we have used (5.23) and (5.21). Therefore, from Proposition 5.2, we
obtain

$\displaystyle[u_{z}]_{C^{0,{\Theta_{1}}-\epsilon}({B_{1/8}})}\leq
C_{1}(N,p,q,s_{1},s_{2},\Bbbk,\|f\|_{L^{\Bbbk}},\|a\|_{\infty},\epsilon).$

By scaling back to $u$ and since $z\in B_{R/2}(x_{0})$ was arbitrary, by a
standard covering argument, we can complete the proof of the theorem.
$\square$

######  Remark 3.

It is clear from the above proof that for $\Bbbk=\infty$, Theorem 2.2 holds
for the case $ps_{1}=qs_{2}$, too. Indeed, in this case, we approximate only
the modulating coefficient $a(\cdot,\cdot)$ in Lemma 5.1 by keeping $f$ fixed.
In the proof of Theorem 2.2, we proceed with
$\mathcal{M}(R):=\mathcal{C}\|u\|_{L^{\infty}(B_{R})}+1+\sum_{(\ell,s)}2^{\frac{N+s\ell}{\ell-1}}T^{\ell}_{\ell
s}(u;x_{0},R)^{\frac{1}{\ell-1}},$ instead of $\mathcal{M}(R_{0})$ and take
${\delta_{1}\mathcal{M}(R)^{p-q}}$ in place of $\delta_{1}$ in (5.22).

###  5.2 Boundary regularity

We define the distance function as $d(x):={\rm
dist}(x,\mathbb{R}^{N}\setminus\Omega)$. Then, we have the following estimate

######  Proposition 5.3.

Let $\Omega\subset\mathbb{R}^{N}$ be a bounded domain with $C^{1,1}$ boundary
$\partial\Omega$ and for $s_{1}<q^{\prime}s_{2}$, we assume that
$a(\cdot,\cdot)\in C^{0,\beta}(\mathbb{R}^{N}\times\mathbb{R}^{N})$ with
$\beta>s_{2}$. Let $u\in\mathcal{W}(\mathbb{R}^{N})\cap L^{\infty}_{\rm
loc}(\Omega)$ be such that $|\mathcal{L}u|\leq K$, weakly in $\Omega$, for
some $K>0$. Then, for all $0<\sigma<s_{1}$, there exists a constant
$\mathcal{C}_{u}>0$ (depending only on data of the problem and
$\|u\|_{L^{\infty}_{\rm loc}(\Omega)}$) such that

$|u|\leq\mathcal{C}_{u}\;d^{\sigma}\quad\mbox{in }\Omega.$

Proof.  We first obtain a lower bound on the quantity
$(-\Delta)_{q,a}^{s_{2}}d^{\alpha}$, for all $\alpha\in[s_{2},1)$. Indeed,
proceeding as in [23, Lemma 3.12], the case $\alpha>q^{\prime}s_{2}$ is
exactly the same as Case (i) of the proof presented there, that is,

$\displaystyle(-\Delta)_{q,a}^{s_{2}}\ d^{\alpha}=g\in
L^{\infty}(\Omega_{\rho})\quad\mbox{weakly in }\Omega_{\rho},$ (5.25)

where $\Omega_{\rho}:=\\{x\in\Omega\ :\ d(x)<\rho\\}$. For the case
$\alpha<q^{\prime}s_{2}$, a similar procedure to that of [23, Lemma 3.12]
yields that

$\displaystyle(-\Delta)_{q,\tilde{a}(\cdot)}^{s_{2}}\ d^{\alpha}\geq
C_{1}d^{\alpha(q-1)-qs_{2}}\quad\mbox{weakly in }\Omega_{\rho},$ (5.26)

where $\tilde{a}(x)=a(x,x)$ and $C_{1}$ is a positive constant. Furthermore,
proceeding with the same notation as in [23, Lemma 3.12], the integrand
$\tilde{I}_{2,\epsilon}(X)$ is estimated as below:

$\displaystyle\tilde{I}_{2,\epsilon}(X)$
$\displaystyle=\int_{B_{\epsilon}(X)^{c}\cap
B_{2\rho}}a(\Phi(X),\Phi(X))\frac{[(X_{N})_{+}^{\alpha}-(Y_{N})_{+}^{\alpha}]^{q-1}}{|\nabla\Phi(X)(X-Y)|^{N+qs_{2}}}J_{\Phi}(X)dY$
$\displaystyle+\int_{B_{\epsilon}(X)^{c}\cap
B_{2\rho}}\big{[}a(\Phi(X),\Phi(Y))-a(\Phi(X),\Phi(X))\big{]}\frac{[(X_{N})_{+}^{\alpha}-(Y_{N})_{+}^{\alpha}]^{q-1}}{|\nabla\Phi(X)(X-Y)|^{N+qs_{2}}}J_{\Phi}(X)dY$
$\displaystyle=:I_{5,\epsilon}(X)+I_{6,\epsilon}(X).$

Using the $\beta$-Hölder continuity of $a(\cdot,\cdot)$ and the Lipschitz
nature of $\Phi$ and $\nabla\Phi$, we obtain

$\displaystyle|I_{6,\epsilon}(X)|\leq C\int_{B_{\epsilon}(X)^{c}\cap
B_{2\rho}}\frac{1}{|X-Y|^{N+qs_{2}-\alpha(q-1)-\beta}}dY=C_{\rho}<\infty,$

where we have used the fact that $\beta>qs_{2}-\alpha(q-1)$. The term
$I_{5,\epsilon}(X)$ is handled as in the estimate to prove (5.26). This
concludes that, for $\alpha<q^{\prime}s_{2}$,

$\displaystyle(-\Delta)_{q,a}^{s_{2}}\ d^{\alpha}\geq
C_{2}d^{\alpha(q-1)-qs_{2}}\quad\mbox{weakly in }\Omega_{\rho}.$ (5.27)

Also, note that $d^{\alpha}\in W^{s_{1},p}(\Omega^{\prime})\cap
W^{s_{2},q}(\Omega^{\prime})$, for some $\Omega\Subset\Omega^{\prime}$ and for
all $s_{2}-1/q<\alpha<1$. Consequently, proceeding similar to [23, Proposition
3.14] and taking into account the local boundedness of $u$ together with
(5.25) and (5.27), we can complete the proof of the proposition. $\square$

Finally, we have our global Hölder continuity result.
Proof of Corollary 2.3: The proof can be completed by standard covering
arguments and using the interior regularity result of Theorem 2.2 (note that
$\Bbbk=\infty$, here) together with the boundary estimate of Proposition 5.3.
$\square$

##  Appendix A On the existence result

Here, we give the functional setup and a basic framework for the existence
result for the following problem:

$\left\\{\begin{array}[]{rl}\partial_{t}u+(-\Delta)_{p}^{s_{1}}u+(-\Delta)_{q,a}^{s_{2}}u&=f\quad\mbox{in
}\Omega\times I,\\\ u&=g\quad\mbox{in }(\mathbb{R}^{N}\setminus\Omega)\times
I,\\\ u(\cdot,t_{0})&=u_{0}(\cdot)\quad\mbox{in }\Omega,\end{array}\right.$
($GNP$)

where $I=[t_{0},t_{1}]$, $u_{0}\in L^{2}(\Omega)$,
$f\in\big{(}\mathcal{Y}(I,\mathcal{W}(\Omega^{\prime}))\big{)}^{*}$,
$g\in\mathcal{Y}(I,\mathcal{W}(\Omega^{\prime}))\cap
L^{p-1}(I;L^{p-1}_{s_{1}p}(\mathbb{R}^{N}))\cap
L^{q-1}(I;L^{q-1}_{s_{2}q,a}(\mathbb{R}^{N}))$ and
$\partial_{t}g\in\big{(}\mathcal{Y}(I,\mathcal{W}(\Omega^{\prime}))\big{)}^{*}$,
for $\Omega\Subset\Omega^{\prime}\subset\mathbb{R}^{N}$.

######  Definition A.1.

A function $u\in\mathcal{Y}(I,\mathcal{W}(\Omega^{\prime}))\cap
L^{p-1}(I;L^{p-1}_{s_{1}p}(\mathbb{R}^{N}))\cap
L^{q-1}(I;L^{q-1}_{s_{1}q,a}(\mathbb{R}^{N}))\cap C(I;L^{2}(\Omega))$ is said
to be a weak solution to problem ($GNP$) if the following hold:

  * (i)

$u\in X_{\rm g(t)}(\Omega,\Omega^{\prime})$ for a.e. $t\in I$, where ${\rm
g(t)}(x)=g(x,t)$;

  * (ii)

$u(\cdot,t)\to u_{0}$ in $L^{2}(\Omega)$, as $t\to t_{0}$;

  * (iii)

for all $J=[T_{0},T_{1}]$ and every
$\phi\in\mathcal{Y}(J;X_{0}(\Omega,\Omega^{\prime}))\cap
C^{1}(J;L^{2}(\Omega))$

$\displaystyle-\int_{J}\int_{\Omega}u(x,t)\partial_{t}\phi(x,t)dxdt+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}\frac{[u(x,t)-u(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi(x,t)-\phi(y,t))dxdydt$
$\displaystyle\quad+\int_{J}\int_{\mathbb{R}^{N}}\int_{\mathbb{R}^{N}}a(x,y)\frac{[u(x,t)-u(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi(x,t)-\phi(y,t))dxdydt$
$\displaystyle=\int_{\Omega}u(x,T_{0})\phi(x,T_{0})dx-\int_{\Omega}u(x,T_{1})\phi(x,T_{1})dx+\langle
f(\cdot,\cdot),\phi(\cdot,\cdot)\rangle_{\mathcal{Y},\mathcal{Y}^{*}}.$

Proceeding similar to [39, Proposition 1.2, Chapter III], we can show that

$\displaystyle\mathbb{W}(I):=\\{v\in\mathcal{Y}(I;\mathcal{W}(\Omega))\ :\
u^{\prime}\in\big{(}\mathcal{Y}(I,\mathcal{W}(\Omega))\big{)}^{*}\\}\subset
C(I;L^{2}(\Omega)).$

Additionally, for $v\in\mathbb{W}(I)$,

$\displaystyle t\mapsto\|v(t)\|_{L^{2}(\Omega)}^{2}\mbox{ is absolutely
continuous and }\frac{d}{dt}\|v(t)\|_{L^{2}(\Omega)}^{2}=2\langle
v^{\prime}(t),v(t)\rangle.$

We have the following existence result.

######  Theorem A.1.

Let $p\geq 2$ and suppose that

$\displaystyle\lim_{t\to
t_{0}}\|g(\cdot,t)-g_{0}\|_{L^{2}(\Omega)}=0\quad\mbox{for some }g_{0}\in
L^{2}(\Omega).$

Then, there exists a unique solution to problem ($GNP$).

Proof.  The proof is similar to [6, Theorem A.3]. Here, we highlight only the
key steps, for instance, boundedness and coercivity of the operators involved
in the proof. For a.e. $t\in I$, we define $\mathcal{A}_{t}:X_{{\rm
g}(t)}(\Omega,\Omega^{\prime})\to(\mathcal{W}(\Omega^{\prime}))^{*}$ by

$\displaystyle\langle\mathcal{A}_{t}(v),\phi\rangle=$
$\displaystyle\sum_{(\ell,s)}\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a_{\ell}(x,y)\frac{[v(x)-v(y)]^{\ell-1}}{|x-y|^{N+s\ell}}(\phi(x)-\phi(y))dxdy$
$\displaystyle+2\sum_{(\ell,s)}\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}a_{\ell}(x,y)\frac{[v(x)-g(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}\phi(x)dxdy$
$\displaystyle=:\langle\mathcal{A}^{p}_{t}(v),\phi\rangle+\langle\mathcal{A}^{q}_{t}(v),\phi\rangle.$

Applying Hölder’s inequality and recalling the definition of $W$ (stated as in
(2.2)), we obtain

$\displaystyle|\langle\mathcal{A}^{q}_{t}(v),\phi\rangle|$
$\displaystyle\leq\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)\frac{|v(x)-v(y)|^{q-1}}{|x-y|^{N+s_{2}q}}|\phi(x)-\phi(y)|dxdy$
$\displaystyle\quad+C(q)\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}a(x,y)\frac{|v(x)|^{q-1}+|g(y,t)|^{q-1}}{|x-y|^{N+s_{2}q}}|\phi(x)|dxdy$
$\displaystyle\leq
C[v]^{q-1}_{W^{s_{2},q}_{a}(\Omega^{\prime})}[\phi]_{W^{s_{2},q}_{a}(\Omega^{\prime})}+C\Big{(}\int_{\Omega}W(x)|v(x)|^{q}dx\Big{)}^{\frac{q-1}{q}}\Big{(}\int_{\Omega}W(x)|\phi(x)|^{q}dx\Big{)}^{\frac{1}{q}}$
$\displaystyle\quad+C\Big{(}\int_{\Omega}|\phi(x)|^{p}dx\Big{)}^{\frac{1}{p}}T^{q}_{qs_{2},a}(g;z,r)^{p/(p-1)}$
$\displaystyle\leq
C\big{(}\|v\|^{q-1}_{\mathcal{W}(\Omega^{\prime})}+T^{q}_{qs_{2},a}(g;z,r)^{p/(p-1)}\big{)}\|\phi\|_{\mathcal{W}(\Omega^{\prime})},$
(E.1)

where $z\in\Omega$ and $r={\rm dist}(\Omega,\partial\Omega^{\prime})>0$. An
analogous result holds for the $p$-term, too. Consequently, $\mathcal{A}_{t}$
is a well defined operator. Next, we define
$\mathcal{A}:X_{0}(\Omega,\Omega^{\prime})\to(\mathcal{W}(\Omega^{\prime}))^{*}$
by

$\displaystyle\mathcal{A}(v,t)=\mathcal{A}_{t}(v+{\rm g}(t)).$

Then, the operator $\mathcal{A}$ is monotone and hemicontinuous (see [31,
Lemma 3]). On account of Hölder’s inequality and the definition of norm on
$\mathcal{W}(\Omega^{\prime})$, it is easy to observe that

$\displaystyle\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}\frac{[v(x)+g(x,t)-v(y)-g(y,t)]^{p-1}}{|x-y|^{N+s_{1}p}}(\phi(x)-\phi(y))dxdy$
$\displaystyle+\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)\frac{[v(x)+g(x,t)-v(y)-g(y,t)]^{q-1}}{|x-y|^{N+s_{2}q}}(\phi(x)-\phi(y))dxdy$
$\displaystyle\leq C\Big{(}\|v\|_{W^{s_{1},p}(\Omega^{\prime})}^{p-1}+\|{\rm
g}(t)\|_{W^{s_{1},p}(\Omega^{\prime})}^{p-1}+\|v\|_{W^{s_{2},q}_{a}(\Omega^{\prime})}^{q-1}+\|{\rm
g}(t)\|_{W^{s_{2},q}_{a}(\Omega^{\prime})}^{q-1}\Big{)}\|\phi\|_{\mathcal{W}(\Omega^{\prime})}.$
(E.2)

Moreover, proceeding similar to (A), we can prove that

$\displaystyle\sum_{(\ell,s)}\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}a_{\ell}(x,y)\frac{[v(x)+g(x,t)-g(y,t)]^{\ell-1}}{|x-y|^{N+s\ell}}\phi(x)dxdy$
$\displaystyle\leq C\Big{(}\|v\|_{\mathcal{W}(\Omega^{\prime})}^{p-1}+\|{\rm
g}(t)\|_{\mathcal{W}(\Omega^{\prime})}^{p-1}+\|{\rm
g}(t)\|_{L^{p-1}_{s_{1}p}(\mathbb{R}^{N})}^{p-1}\Big{)}\|\phi\|_{\mathcal{W}(\Omega^{\prime})}$
$\displaystyle\quad+C\Big{(}\|v\|_{\mathcal{W}(\Omega^{\prime})}^{q-1}+\|{\rm
g}(t)\|_{\mathcal{W}(\Omega^{\prime})}^{q-1}+\|{\rm
g}(t)\|_{L^{q-1}_{s_{2}q,a}(\mathbb{R}^{N})}^{q-1}\Big{)}\|\phi\|_{\mathcal{W}(\Omega^{\prime})}.$
(E.3)

Coupling (A) and (A), we obtain

$\displaystyle|\langle\mathcal{A}(v,t),\phi\rangle|$ $\displaystyle\leq
C\Big{(}\|v\|_{\mathcal{W}(\Omega^{\prime})}^{p-1}+\|{\rm
g}(t)\|_{\mathcal{W}(\Omega^{\prime})}^{p-1}+\|{\rm
g}(t)\|_{L^{p-1}_{s_{1}p}(\mathbb{R}^{N})}^{p-1}$
$\displaystyle\quad+\|v\|_{\mathcal{W}(\Omega^{\prime})}^{q-1}+\|{\rm
g}(t)\|_{\mathcal{W}(\Omega^{\prime})}^{q-1}+\|{\rm
g}(t)\|_{L^{q-1}_{s_{2}q,a}(\mathbb{R}^{N})}^{q-1}\Big{)}\|\phi\|_{\mathcal{W}(\Omega^{\prime})},$

which proves the boundedness of $\mathcal{A}$. Next, to obtain coercivity of
the operator $\mathcal{A}$, for any $v\in X_{0}(\Omega,\Omega^{\prime})$,
using (2.4) and Young’s inequality, we first see that

$\displaystyle\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)[v(x)+g(x,t)-v(y)-g(y,t)]^{q-1}(v(x)-v(y))d\mu_{2}$
$\displaystyle=\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)\Big{(}[v(x)+g(x,t)-v(y)+g(y,t)]^{q-1}-[g(x,t)-g(y,t)]^{q-1}\Big{)}(v(x)-v(y))d\mu_{2}$
$\displaystyle\
+\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)[g(x,t)-g(y,t)]^{q-1}(v(x)-v(y))d\mu_{2}$
$\displaystyle\geq
c_{q}\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)|v(x)-v(y)|^{q}d\mu_{2}-c_{q}\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}a(x,y)|g(x,t)-g(y,t)|^{q-1}|v(x)-v(y)|d\mu_{2}$
$\displaystyle\geq
c\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}|v(x)-v(y)|^{q}d\mu_{2}-C\int_{\Omega^{\prime}}\int_{\Omega^{\prime}}|g(x,t)-g(y,t)|^{q}d\mu_{2}.$

Proceeding similarly, and using the definition of $W$, we observe that

$\displaystyle\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}a(x,y)[v(x)+g(x,t)-g(y,t)]^{q-1}v(x)d\mu_{2}$
$\displaystyle\geq
c(q)\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}\frac{a(x,y)}{|x-y|^{N+qt}}|v(x)|^{q}dxdy$
$\displaystyle\
-c(q)\int_{\Omega}\int_{\mathbb{R}^{N}\setminus\Omega^{\prime}}\frac{a(x,y)}{|x-y|^{N+qs_{2}}}(|g(x,t)|^{q-1}+|g(y,t)|^{q-1})|v(x)|dxdy$
$\displaystyle\geq
c\int_{\Omega^{\prime}}W(x)|v(x)|^{q}dx-C\int_{\Omega^{\prime}}W(x)|g(x,t)|^{q}dx-\epsilon\int_{\Omega^{\prime}}|v(x)|^{p}dx-
CT^{q}_{qs_{2},a}(g;z,r)^{\frac{p}{p-1}},$

where we have also used the fact that $v=0$ in
$\Omega^{\prime}\setminus\Omega$ and Young’s inequality on the last line.
Analogous results hold for the corresponding $p$-counterparts, too.
Consequently, recalling the definition of the norm on
$\mathcal{W}(\Omega^{\prime})$, we obtain