# Harnack inequality for nonlocal problems with non-standard growth

Jamil Chaker  Fakultät für Mathematik, Universität Bielefeld, 33615 Bielefeld,
Germany<EMAIL_ADDRESS>,  Minhyun Kim  Fakultät für
Mathematik, Universität Bielefeld, 33615 Bielefeld, Germany  minhyun.kim@uni-
bielefeld.de and  Marvin Weidner  Fakultät für Mathematik, Universität
Bielefeld, 33615 Bielefeld, Germany<EMAIL_ADDRESS>
###### Abstract.

We prove a full Harnack inequality for local minimizers, as well as weak
solutions to nonlocal problems with non-standard growth. The main auxiliary
results are local boundedness and a weak Harnack inequality for functions in a
corresponding De Giorgi class. This paper builds upon a recent work on
regularity estimates for such nonlocal problems by the same authors.

###### Key words and phrases:

Harnack inequality, local boundedness, nonlocal problem, non-standard growth,
minimizer, weak solution, De Giorgi class

###### 2020 Mathematics Subject Classification:

35B65, 47G20, 35D30, 35B45, 35A15

Jamil Chaker gratefully acknowledges funding by the Deutsche
Forschungsgemeinschaft (SFB 1283/2 2021 - 317210226). Minhyun Kim and Marvin
Weidner gratefully acknowledge funding by the Deutsche Forschungsgemeinschaft
(GRK 2235/2 2021 - 282638148).

##  1\. Introduction

The goal of the present work is to prove a full Harnack inequality for local
minimizers and weak solutions to a class of nonlocal problems which exhibit
non-standard growth. This article builds upon the recent paper [CKW21], in
which we study regularity properties for local minimizers of nonlocal energy
functionals, as well as weak solutions to nonlocal equations with non-standard
growth. We prove that these objects satisfy a suitable fractional Caccioppoli
inequality and therefore belong to corresponding De Giorgi classes. In this
work, we show that any function in such De Giorgi class satisfies a full
Harnack inequality. As a consequence, we obtain the full Harnack inequality
for local minimizers and weak solutions.

Before we state the main result of this paper, let us formulate the main
assumptions and briefly present the energy functionals respectively the
nonlocal operators considered in this work. We point out that the setup of
this article is in align with [CKW21].

Let $\Omega\subset\mathbb{R}^{d}$ be open, $s\in(0,1)$ and $1\leq p\leq q$.
Throughout the paper, let $f:[0,\infty)\to[0,\infty)$ be convex, strictly
increasing and differentiable with $f(0)=0$ and $f(1)=1$. We say that $f$
satisfies (1.1) if for all $t\geq 0$:

($f_{p}$) $\displaystyle pf(t)\leq~{}$ $\displaystyle tf^{\prime}(t),$
($f^{q}$) $\displaystyle tf^{\prime}(t)\leq qf(t).$

The growth function $f$ is naturally associated with nonlocal energy
functionals and nonlocal operators. On the one hand, consider

(1.2)
$u\mapsto\mathcal{I}_{f}(u)=(1-s)\iint_{(\Omega^{c}\times\Omega^{c})^{c}}f\left(\frac{|u(x)-u(y)|}{|x-y|^{s}}\right)\frac{k(x,y)}{|x-y|^{d}}\,\mathrm{d}y\,\mathrm{d}x,$

where $k:\mathbb{R}^{d}\times\mathbb{R}^{d}\to\mathbb{R}$ is a measurable
function satisfying

($k$) $k(x,y)=k(y,x)\quad\text{and}\quad\Lambda^{-1}\leq
k(x,y)\leq\Lambda\quad\text{for a.e. }x,y\in\mathbb{R}^{d}$

for some $\Lambda\geq 1$. In [CKW21, Theorem 6.2], we prove that local
minimizers of $\mathcal{I}_{f}$ belong to the De Giorgi class
$G(\Omega;q,c,s,f)$ for some constant $c=c(d,q,\Lambda)>0$ if $f$ satisfies
($f^{q}$) for some $q>1$. For the precise definition of the De Giorgi class,
see Definition 2.8.

On the other hand, we consider weak solutions to

(1.3) $\mathcal{L}_{h}u=0\quad\text{in }\Omega,$

where $\mathcal{L}_{h}$ is a nonlocal operator of the form

$\mathcal{L}_{h}u(x)=(1-s)\text{p.v.}\int_{\mathbb{R}^{d}}h\left(x,y,\frac{u(x)-u(y)}{|x-y|^{s}}\right)\frac{\,\mathrm{d}y}{|x-y|^{d+s}}.$

Here, $h:\mathbb{R}^{d}\times\mathbb{R}^{d}\times\mathbb{R}\to\mathbb{R}$ is a
measurable function satisfying the structure condition

($h$)
$h(x,y,t)=h(y,x,t),\quad\text{sign}(t)\frac{1}{\Lambda}f^{\prime}(|t|)\leq
h(x,y,t)\leq\Lambda f^{\prime}(|t|)$

for a.e. $x,y\in\mathbb{R}^{d}$ and for all $t\in\mathbb{R}$. We show in
[CKW21, Theorem 7.3] that weak solutions to (1.3) are in $G(\Omega;q,c,s,f)$
for some constant $c=c(d,q,\Lambda)>0$ if ($f^{q}$) holds true for $q>1$.

The goal of this article is to prove a full Harnack inequality of the
following form.

######  Theorem 1.1.

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $0<s_{0}\leq s<1$,
$1<p\leq q$, $c>0$ and assume that $f$ satisfies (1.1). There exists a
constant $C>0$ such that if $u\in G(\Omega;q,c,s,f)$ is nonnegative in
$B_{R}(x_{0})\subset\Omega$, then

(1.4) $\sup_{B_{R/2}(x_{0})}u\leq
C\left(\inf_{B_{R/2}(x_{0})}u+\mathrm{Tail}_{f^{\prime}}(u_{-};x_{0},R)\right).$

The constant $C$ depends only on $d$, $s_{0}$, $p$, $q$ and $c$.

The Harnack inequality was originally proved for harmonic functions and later
obtained for several elliptic and parabolic local operators. It is known to
have important consequences such as a priori estimates in Hölder spaces or
convergence theorems. Therefore, it plays an important role in several
mathematical fields such as geometric analysis, probability or analysis of
partial differential equations. For an introduction to Harnack inequalities,
their history and consequences, we refer the reader to the article by Kassmann
[Kas07b]. The appearance of the tail term on the right-hand side of (1.4) is a
purely nonlocal phenomenon. It is shown in [Kas07a] that the classical Harnack
inequality fails for $s$-harmonic functions if nonnegativity of the function
is assumed in the solution domain only. In [Kas11], a new formulation of the
Harnack inequality is introduced. It involves a nonlocal tail term as in (1.4)
which captures the negative values of the $s$-harmonic function outside the
solution domain. In our setup the nonlocal tail is of the following form

$\mathrm{Tail}_{f^{\prime}}(u;x_{0},R)=R^{s}(f^{\prime})^{-1}\left((1-s)R^{s}\int_{\mathbb{R}^{d}\setminus
B_{R}(x_{0})}f^{\prime}\left(\frac{|u(y)|}{|y-x_{0}|^{s}}\right)\frac{\textnormal{d}y}{|y-x_{0}|^{d+s}}\right),$

see Section 2.2 for details.
Further important contributions to investigation of Harnack inequalities for
nonlocal operators are, among others, the articles [BL02, CK03, SV04, BK05,
BS05, CS07, CS09, BBCK09, SZ13, DCKP14, Coz17, Str19, CKW19] and the
references therein.

Since local minimizers of (1.2) belong to the De Giorgi class, we have the
following corollary of Theorem 1.1, that is the full Harnack inequality for
local minimizers.

######  Corollary 1.2.

Let $s_{0}\in(0,1)$, $1<p\leq q$, $\Lambda\geq 1$ and assume $s\in[s_{0},1)$.
Assume that $f$ satisfies (1.1) and let
$k:\mathbb{R}^{d}\times\mathbb{R}^{d}\to\mathbb{R}$ be a measurable function
satisfying ($k$). There exists a constant $C>0$, depending only on $d$,
$s_{0}$, $p$, $q$ and $\Lambda$, such that if $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$ is a local minimizer of (1.2) that is
nonnegative in $B_{R}(x_{0})\subset\Omega$, then the full Harnack inequality
(1.4) holds true for $u$.

Another direct consequence of Theorem 1.1, together with the observation that
weak solutions belong to the De Giorgi class, is the full Harnack inequality
for weak solutions.

######  Corollary 1.3.

Let $s_{0}\in(0,1)$, $1<p\leq q$, $\Lambda\geq 1$ and assume $s\in[s_{0},1)$.
Assume that $f$ satisfies (1.1) and let
$h:\mathbb{R}^{d}\times\mathbb{R}^{d}\times\mathbb{R}\to\mathbb{R}$ be a
measurable function satisfying ($h$). There exists a constant $C>0$, depending
only on $d$, $s_{0}$, $p$, $q$ and $\Lambda$, such that if $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$ is a weak solution to (1.3) that is
nonnegative in $B_{R}(x_{0})\subset\Omega$, then the full Harnack inequality
(1.4) holds true for $u$.

The proof of the main result Theorem 1.1 follows from a weak Harnack
inequality together with the local boundedness of functions in the De Giorgi
class. Our approach roughly follows the ideas of Mascolo and Papi [MP96],
where they establish a Harnack inequality for minimizers of functionals with
non-standard growth in the local case.

We would like to point out that all results in the present paper are robust in
the sense that the constants stay uniform as $s\to 1^{-}$, since they depend
on $s_{0}$ and are independent of the actual order of differentiability $s$.
Since the tail contribution vanishes as $s\to 1^{-}$, we recover purely local
estimates in the limit case. Note that the nonlocal energy functional is known
to converge to a local energy form as $s\to 1^{-}$, see for instance [FBS19,
ACPS21].

The family of operators studied in this paper exhibits non-standard growth
behavior. The regularity theory has been intensively studied in recent years.
In [CKW21] and [BKO21], local boundedness and local Hölder regularity are
established for such operators using two different approaches and under
slightly different conditions on the growth function. However, the condition
(1.1) is the main assumption in both articles. For a deeper discussion on the
literature about nonlocal operators with different types of non-standard
growth behavior and their regularity theory, we refer the reader to the
references given in those two articles. See also [DFP19, GKS20, BSV20, CK21,
Ok21, GKS21, BOS21, FZ21, GKS22] and the references therein.

Moreover, we want to point out that the local boundedness estimate in this
paper, see Theorem 3.1, significantly improves the previous result from
[CKW21] and is proved without the condition $q<p^{\ast}$. Note that a similar
result is established in [BKO21].

Recently, Fang and Zhang have investigated Harnack inequalities for nonlocal
operators with general growth [FZ22]. In comparison to our setup, they impose
more restrictive structural assumptions on the growth function $f$. Similar to
the approach in this article, they derive local boundedness and a tail
estimate as in Theorem 3.1, Lemma 5.1, as well as a weak Harnack inequality.
By combining these results, [FZ22] derive an upper estimate for $\sup u$ in
terms of $\inf u$ and a nonlocal tail term. However, for $p<q$, this result is
not optimal due to the appearance of an additional power $\iota=q/p$ in the
Harnack inequality. In this article, we prove a different version of a weak
Harnack inequality taking into account the growth function $f$, see Theorem
4.1. This allows us to deduce a full Harnack inequality in the classical form
(1.4).

### Notation

Throughout the paper, we will denote by $C>0$ a universal constant, which may
be different from line to line.

### Outline

This article is structured as follows. In Section 2 we collect several
auxiliary results for the growth functions under consideration and provide
definitions of related function spaces and De Giorgi classes. Section 3 and
Section 4 are devoted to the proof of local boundedness and a weak Harnack
inequality for functions $u$ in appropriate De Giorgi classes. Finally, the
proof of the main result Theorem 1.1 is provided in Section 5.

##  2\. Preliminaries

This section contains several auxiliary results on the growth function $f$ and
introduces the function spaces related to our setup.

###  2.1. Properties of growth functions

We collect several properties of growth functions $f:[0,\infty)\to[0,\infty)$
which were proved in [CKW21] and will be used in the course of this article.
Recall that we will assume throughout this paper that $f$ is convex, strictly
increasing and differentiable with $f(0)=0$ and $f(1)=1$.

######  Lemma 2.1.

[CKW21, Lemma 2.1] Let $q\geq 1$. Then the following are equivalent:

  1. (i)

($f^{q}$),

  2. (ii)

$t\mapsto t^{-q}f(t)$ is decreasing,

  3. (iii)

$f(\lambda t)\leq\lambda^{q}f(t)$ for all $\lambda\geq 1$,

  4. (iv)

$\lambda^{q}f(t)\leq f(\lambda t)$ for all $\lambda\leq 1$.

######  Lemma 2.2.

[CKW21, Lemma 2.2] Let $p\geq 1$. Then the following are equivalent:

  1. (i)

($f_{p}$),

  2. (ii)

$t\mapsto t^{-p}f(t)$ is increasing,

  3. (iii)

$\lambda^{p}f(t)\leq f(\lambda t)$ for all $\lambda\geq 1$,

  4. (iv)

$f(\lambda t)\leq\lambda^{p}f(t)$ for all $\lambda\leq 1$.

######  Corollary 2.3.

[CKW21, Corollary 2.4] Let $1\leq p\leq q$. Assume that $f$ satisfies (1.1).
Then,

(2.1) $\displaystyle\frac{p}{q}\lambda^{p-1}f^{\prime}(t)\leq
f^{\prime}(\lambda t)\leq\frac{q}{p}\lambda^{q-1}f^{\prime}(t)\quad\text{for
all }\lambda\geq 1,$ (2.2)
$\displaystyle\frac{p}{q}\lambda^{q-1}f^{\prime}(t)\leq f^{\prime}(\lambda
t)\leq\frac{q}{p}\lambda^{p-1}f^{\prime}(t)\quad\text{for all }\lambda\leq 1,$
(2.3) $\displaystyle\frac{1}{2}f^{\prime}(t)+\frac{1}{2}f^{\prime}(s)\leq
f^{\prime}(t+s)\leq\frac{q}{p}2^{q-1}(f^{\prime}(t)+f^{\prime}(s))\quad\text{for
all }t,s\geq 0.$

######  Lemma 2.4.

[CKW21, Lemma 2.5] Let $c>1$ and assume that for some $t,s>0$ it holds that
$f(t)\leq cf(s)$. Then $t\leq cs$.

Note that under the assumptions on $f$, it does not necessarily follow that
$f^{\prime}$ is invertible. Throughout this article, we will work with the
following generalized inverse of $f^{\prime}$:

(2.4) $(f^{\prime})^{-1}(y)=\inf\\{t:f^{\prime}(t)\geq y\\}.$

We collect a few properties of $(f^{\prime})^{-1}$. First, we recall a
proposition from [CKW21].

######  Proposition 2.5.

[CKW21, Proposition 3.1] It holds that

(2.5) $\displaystyle(f^{\prime}\circ(f^{\prime})^{-1})(y)\geq y\quad\text{for
all }y\geq 0,$ (2.6) $\displaystyle((f^{\prime})^{-1}\circ f^{\prime})(t)\leq
t\quad\text{for all }t\geq 0.$

The following are simple consequences of the previous results.

######  Lemma 2.6.

For every $t,s\geq 0$:

$\displaystyle(f^{\prime})^{-1}\left(\frac{t+s}{2}\right)\leq(f^{\prime})^{-1}(t)+(f^{\prime})^{-1}(s).$

###### Proof.

By (2.5), (2.6) and monotonicity of $(f^{\prime})^{-1}$:

$\displaystyle(f^{\prime})^{-1}\left[\frac{t+s}{2}\right]$
$\displaystyle\leq(f^{\prime})^{-1}\left[\frac{f^{\prime}((f^{\prime})^{-1}(t))+f^{\prime}((f^{\prime})^{-1}(s))}{2}\right]\leq(f^{\prime})^{-1}\left[f^{\prime}\left((f^{\prime})^{-1}(t)+(f^{\prime})^{-1}(s)\right)\right]$
$\displaystyle\leq(f^{\prime})^{-1}(t)+(f^{\prime})^{-1}(s).$

∎

######  Lemma 2.7.

Let $1<p\leq q$. Assume that $f$ satisfies (1.1). Then

$\displaystyle(f^{\prime})^{-1}(\lambda t)\leq
c_{\lambda}(f^{\prime})^{-1}(t)\quad\text{for all }\lambda\geq 0,$

where $c_{\lambda}=(q\lambda/p)^{1/(p-1)}$ if $\lambda\geq p/q$ and
$c_{\lambda}=(q\lambda/p)^{1/(q-1)}$ if $\lambda\leq p/q$.

###### Proof.

First, we observe that by (2.1) and (2.2), $\lambda f^{\prime}(t)\leq
f^{\prime}(c_{\lambda}t)$. Therefore, using (2.5), (2.6) and monotonicity of
$(f^{\prime})^{-1}$:

$\displaystyle(f^{\prime})^{-1}\left[\lambda
t\right]\leq(f^{\prime})^{-1}\left[\lambda
f^{\prime}((f^{\prime})^{-1}(t))\right]\leq(f^{\prime})^{-1}\left[f^{\prime}(c_{\lambda}(f^{\prime})^{-1}(t))\right]\leq
c_{\lambda}(f^{\prime})^{-1}(t).$

∎

###  2.2. Function spaces and De Giorgi classes

Let $s\in(0,1)$ and $\Omega\subset\mathbb{R}^{d}$ be open. We define the
Orlicz and Orlicz–Sobolev spaces by

$\begin{split}L^{f}(\Omega)&=\\{u:\Omega\to\mathbb{R}~{}\text{measurable}:\Phi_{L^{f}(\Omega)}(u)<\infty\\},\\\
W^{s,f}(\Omega)&=\\{u\in L^{f}(\Omega):\Phi_{W^{s,f}(\Omega)}(u)<\infty\\},\\\
V^{s,f}(\Omega|\mathbb{R}^{d})&=\\{u\in
L^{f}(\Omega):\Phi_{V^{s,f}(\Omega)}(u)<\infty\\},\end{split}$

where $\Phi_{L^{f}(\Omega)}$, $\Phi_{W^{s,f}(\Omega)}$ and
$\Phi_{V^{s,f}(\Omega|\mathbb{R}^{d})}$ are modulars defined by

$\begin{split}\Phi_{L^{f}(\Omega)}(u)&=\int_{\Omega}f(|u(x)|)\,\mathrm{d}x,\\\
\Phi_{W^{s,f}(\Omega)}(u)&=(1-s)\int_{\Omega}\int_{\Omega}f\left(\frac{|u(x)-u(y)|}{|x-y|^{s}}\right)\frac{\,\mathrm{d}y\,\mathrm{d}x}{|x-y|^{d}},\\\
\Phi_{V^{s,f}(\Omega|\mathbb{R}^{d})}(u)&=(1-s)\iint_{(\Omega^{c}\times\Omega^{c})^{c}}f\left(\frac{|u(x)-u(y)|}{|x-y|^{s}}\right)\frac{\,\mathrm{d}y\,\mathrm{d}x}{|x-y|^{d}}.\end{split}$

Next, we define nonlocal tails, which capture the behavior of functions $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$ at large scales. We define the nonlocal
$f^{\prime}$-Tail by

(2.7)
$\mathrm{Tail}_{f^{\prime}}(u;x_{0},R)=R^{s}(f^{\prime})^{-1}\left((1-s)R^{s}\int_{\mathbb{R}^{d}\setminus
B_{R}(x_{0})}f^{\prime}\left(\frac{|u(y)|}{|y-x_{0}|^{s}}\right)\frac{\textnormal{d}y}{|y-x_{0}|^{d+s}}\right).$

Recall that the function $f^{\prime}$ does not have to be invertible. Here
$(f^{\prime})^{-1}$ denotes the generalized inverse, see (2.4). In our
previous work [CKW21], we prove that this expression is finite if $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$ for $B_{R}(x_{0})\subset\Omega$.

We are now ready to provide the definition of De Giorgi classes.

######  Definition 2.8 (De Giorgi classes).

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $s\in(0,1)$, $q>1$ and
$c>0$. We say that $u\in G_{+}(\Omega;q,c,s,f)$ if $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$ and if for every $x_{0}\in\Omega$, $0<r<R\leq
d(x_{0},\partial\Omega)$ and $k\in\mathbb{R}$, it holds that

(2.8)
$\begin{split}&\Phi_{W^{s,f}(B_{r}(x_{0}))}(w_{+})+(1-s)\int_{B_{r}(x_{0})}\int_{A_{k}^{-}}f^{\prime}\left(\frac{w_{-}(y)}{|x-y|^{s}}\right)\frac{w_{+}(x)}{|x-y|^{s}}\frac{\textnormal{d}y\,\mathrm{d}x}{|x-y|^{d}}\\\
&\leq
c\left(\frac{R}{R-r}\right)^{q}\Phi_{L^{f}(B_{R}(x_{0}))}\left(\frac{w_{+}}{R^{s}}\right)\\\
&\quad+c\left(\frac{R}{R-r}\right)^{d+sq}\|w_{+}\|_{L^{1}(B_{R}(x_{0}))}(1-s)\int_{\mathbb{R}^{d}\setminus
B_{r}(x_{0})}f^{\prime}\left(\frac{w_{+}(y)}{|y-x_{0}|^{s}}\right)\frac{\textnormal{d}y}{|y-x_{0}|^{d+s}},\end{split}$

where $w_{\pm}(x)=(u(x)-k)_{\pm}$ and
$A_{k}^{-}=\\{y\in\mathbb{R}^{d}:u(y)<k\\}$. We say that $u\in
G_{-}(\Omega;q,c,s,f)$ if (2.8) holds with $w_{+}$, $w_{-}$ and $A_{k}^{-}$
replaced by $w_{-}$, $w_{+}$ and $A_{k}^{+}=\\{y\in\mathbb{R}^{d}:u(y)>k\\}$,
respectively. Moreover, we denote by
$G(\Omega;q,c,s,f)=G_{+}(\Omega;q,c,s,f)\cap G_{-}(\Omega;q,c,s,f)$.

##  3\. Local boundedness

The goal of this section is to prove local boundedness of functions $u\in
G_{+}(\Omega;q,c,s,f)$ under ($f^{q}$), see Theorem 3.1. This result
significantly improves [CKW21, Theorem 5.1]. Let us mention that a similar
estimate has been obtained in [BKO21] using a different technique based on a
Poincaré–Sobolev-type inequality for nonlocal Orlicz–Sobolev spaces. Our proof
solely relies on the classical fractional Sobolev embedding and our estimate
is robust for $s\to 1^{-}$.

######  Theorem 3.1.

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $0<s_{0}\leq s<1$,
$q>1$, $c>0$ and assume that $f$ satisfies ($f^{q}$). There exists a constant
$C>0$ such that if $u\in G_{+}(\Omega;q,c,s,f)$, then for any
$B_{R}(x_{0})\subset\Omega$, $1/2\leq\rho<\tau\leq 1$ and $\delta\in(0,1)$, it
holds that

$\begin{split}f\left(\sup_{B_{\rho R}(x_{0})}\frac{u_{+}}{R^{s}}\right)&\leq
C\frac{\delta^{(1-q)2d/s_{0}}}{(\tau-\rho)^{\gamma}}\fint_{B_{\tau
R}(x_{0})}f\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x+\delta
f\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};x_{0},R/2)}{(R/2)^{s}}\right),\end{split}$

where $\gamma=2d(d+q)/s_{0}$. The constant $C$ depends only on $d$, $s_{0}$,
$q$ and $c$.

######  Remark 3.2.

In particular, Theorem 3.1 implies that functions $u\in G_{+}(\Omega;q,c,s,f)$
are locally bounded from above in $\Omega$ under the assumptions of Theorem
3.1. Local boundedness from below for functions $u\in G_{-}(\Omega;q,c,s,f)$
can be proved in the same way. Finiteness of
$\mathrm{Tail}_{f^{\prime}}(u_{+};x_{0},R/2)$ is a consequence of $u\in
V^{s,f}(\Omega|\mathbb{R}^{d})$, see [CKW21, Proposition 3.2].

###### Proof.

We may assume that $x_{0}=0$. For $j\geq 0$, we set

$\begin{split}&R_{j}=\rho R+2^{-j}(\tau-\rho)R,\quad B_{j}=B_{R_{j}},\\\
&k_{j}=(1-2^{-j})k,\quad\tilde{k}_{j}=(k_{j}+k_{j+1})/2,\\\
&w_{j}=(u-k_{j})_{+},\quad\tilde{w}_{j}=(u-\tilde{k}_{j})_{+},\end{split}$

where $k$ is a positive number that will be determined later. Note that
$R_{j+1}<R_{j}\leq 2R_{j+1}$, $k_{j}\leq\tilde{k}_{j}\leq k_{j+1}$ and
$w_{j+1}\leq\tilde{w}_{j}\leq w_{j}$. We denote by $A^{+}_{h,r}$ the set
$\\{x\in B_{r}:u(x)>h\\}$.

Let $\sigma=\max\\{s_{0}/2,(3s-1)/2\\}\in(0,s)$. Then, it is easy to check
that

(3.1) $1-\sigma\leq\frac{3}{2}(1-s)\quad\text{and}\quad
s-\sigma\geq\min\\{s_{0}/2,(1-\sigma)/3\\}.$

First, by Hölder’s inequality we have

$\begin{split}\fint_{B_{j+1}}f\left(\frac{w_{j+1}(x)}{R^{s}}\right)\,\mathrm{d}x&\leq\frac{1}{|B_{j+1}|}\int_{A_{\tilde{k}_{j},R_{j+1}}^{+}}f\left(\frac{\tilde{w}_{j}(x)}{R^{s}}\right)\,\mathrm{d}x\\\
&\leq\frac{|A_{\tilde{k}_{j},R_{j+1}}^{+}|^{\frac{\sigma}{d}}}{2^{-d}|B_{j}|}\left(\int_{B_{j+1}}f\left(\frac{\tilde{w}_{j}(x)}{R^{s}}\right)^{\frac{d}{d-\sigma}}\,\mathrm{d}x\right)^{\frac{d-\sigma}{d}}.\end{split}$

By applying the fractional Sobolev inequality to $\tilde{w}_{j}/R^{s}$ in
$B_{j+1}$, we estimate

$\begin{split}\left(\int_{B_{j+1}}f\left(\frac{\tilde{w}_{j}(x)}{R^{s}}\right)^{\frac{d}{d-\sigma}}\,\mathrm{d}x\right)^{\frac{d-\sigma}{d}}&\leq
C(1-\sigma)\int_{B_{j+1}}\int_{B_{j+1}}\frac{|f(\tilde{w}_{j}(x)/R^{s})-f(\tilde{w}_{j}(y)/R^{s})|}{|x-y|^{d+\sigma}}\,\mathrm{d}y\,\mathrm{d}x\\\
&\quad+CR_{j+1}^{-\sigma}\int_{B_{j+1}}f\left(\frac{\tilde{w}_{j}(x)}{R^{s}}\right)\,\mathrm{d}x.\end{split}$

Note that from monotonicity of $f^{\prime}$ and assumption ($f^{q}$) it
follows:

$\begin{split}\frac{|f(a/R^{s})-f(b/R^{s})|}{|x-y|^{\sigma}}&\leq\max\left\\{f^{\prime}\left(\frac{a}{R^{s}}\right),f^{\prime}\left(\frac{b}{R^{s}}\right)\right\\}\frac{|a-b|}{|x-y|^{s}}R^{-s}|x-y|^{s-\sigma}\\\
&\leq
q\left(\max\left\\{f\left(\frac{a}{R^{s}}\right),f\left(\frac{b}{R^{s}}\right)\right\\}+f\left(\frac{|a-b|}{|x-y|^{s}}\right)\right)R^{-s}|x-y|^{s-\sigma}\end{split}$

for any $a,b\in\mathbb{R}$. This inequality applied to $a=\tilde{w}_{j}(x)$,
$b=\tilde{w}_{j}(y)$, together with (3.1), yields

$\begin{split}&(1-\sigma)\int_{B_{j+1}}\int_{B_{j+1}}\frac{|f(\tilde{w}_{j}(x)/R^{s})-f(\tilde{w}_{j}(y)/R^{s})|}{|x-y|^{d+\sigma}}\,\mathrm{d}y\,\mathrm{d}x\\\
&\leq
C(1-\sigma)\int_{B_{j}}f\left(\frac{\tilde{w}_{j}(x)}{R^{s}}\right)\int_{B_{j}}\frac{R^{-s}|x-y|^{s-\sigma}}{|x-y|^{d}}\,\mathrm{d}y\,\mathrm{d}x+CR^{-s}R_{j+1}^{s-\sigma}\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})\\\
&\leq\frac{C}{R_{j}^{\sigma}}\left(|B_{j}|\fint_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x+\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})\right).\end{split}$

We have thus far obtained

(3.2) $\fint_{B_{j+1}}f\left(\frac{w_{j+1}(x)}{R^{s}}\right)\,\mathrm{d}x\leq
C\left(\frac{|A_{\tilde{k}_{j},R_{j+1}}^{+}|}{|B_{j}|}\right)^{\frac{\sigma}{d}}\left(\fint_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x+\frac{\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})}{|B_{j}|}\right),$

where the constant $C$ depends only on $d$, $s_{0}$ and $q$ at this point.
In order to set up a suitable iteration scheme based on (3.2), it remains to
estimate the quantity $\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})$. Since $u\in
G_{+}(\Omega;q,c,s,f)$, we have

(3.3) $\begin{split}\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})&\leq
c\left(\frac{R_{j}}{R_{j}-R_{j+1}}\right)^{q}\int_{B_{j}}f\left(\frac{\tilde{w}_{j}(x)}{R_{j}^{s}}\right)\,\mathrm{d}x\\\
&\quad+c\left(\frac{R_{j}}{R_{j}-R_{j+1}}\right)^{d+sq}\|\tilde{w}_{j}\|_{L^{1}(B_{j})}\int_{\mathbb{R}^{d}\setminus
B_{j+1}}f^{\prime}\left(\frac{\tilde{w}_{j}(y)}{|y|^{s}}\right)\frac{1-s}{|y|^{d+s}}\,\mathrm{d}y\\\
&=:I_{1}+I_{2}.\end{split}$

For $I_{1}$, we use Lemma 2.1 and the fact that $R_{j}\leq R$ to deduce

(3.4) $I_{1}\leq
c\left(\frac{R_{j}}{R_{j}-R_{j+1}}\right)^{q}\left(\frac{R}{R_{j}}\right)^{sq}\int_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x\leq
C\frac{2^{qj}}{(\tau-\rho)^{q}}\int_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x.$

For $I_{2}$, using monotonicity of $f^{\prime}$ and the assumption ($f^{q}$)
again, we observe that

$\begin{split}\frac{\tilde{w}_{j}}{R^{s}}f^{\prime}\left(\frac{k}{R^{s}}\right)&\leq
q2^{(q-1)(j+2)}\frac{\tilde{w}_{j}}{R^{s}}f^{\prime}\left(\frac{\tilde{k}_{j}-k_{j}}{R^{s}}\right)\\\
&\leq
q2^{(q-1)(j+2)}\frac{w_{j}}{R^{s}}f^{\prime}\left(\frac{w_{j}}{R^{s}}\right)\leq
q^{2}2^{(q-1)(j+2)}f\left(\frac{w_{j}}{R^{s}}\right).\end{split}$

Thus, $I_{2}$ can be estimated as

$I_{2}\leq
C\frac{2^{(d+2q)j}}{(\tau-\rho)^{d+sq}}\left(\int_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x\right)\frac{R^{s}}{f^{\prime}(k/R^{s})}\int_{\mathbb{R}^{d}\setminus
B_{\rho
R}}f^{\prime}\left(\frac{u_{+}(y)}{|y|^{s}}\right)\frac{1-s}{|y|^{d+s}}\,\mathrm{d}y.$

If we assume that for some $\delta\in(0,1)$, $k\geq k_{1}:=\delta
2^{s}\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)$, it follows:

(3.5) $I_{2}\leq
C\frac{\delta^{1-q}}{(\tau-\rho)^{d+sq}}2^{(d+2q)j}\int_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x,$

where we used (2.2), and $C=C(d,s_{0},q,c)$ is a positive constant. Combining
(3.3), (3.4) and (3.5):

(3.6) $\frac{1}{|B_{j}|}\Phi_{W^{s,f}(B_{j+1})}(\tilde{w}_{j})\leq
C\frac{\delta^{1-q}}{(\tau-\rho)^{d+q}}2^{(d+2q)j}\fint_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x.$

Since we have by Lemma 2.1

(3.7)
$f\left(\frac{k}{R^{s}}\right)\frac{|A_{\tilde{k}_{j},R_{j+1}}^{+}|}{|B_{j}|}\leq
q2^{q(j+2)}f\left(\frac{\tilde{k}_{j}-k_{j}}{R^{s}}\right)\frac{|A_{\tilde{k}_{j},R_{j}}^{+}|}{|B_{j}|}\leq
q2^{q(j+2)}\fint_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x,$

it follows from (3.2), (3.6) and (3.7) that

$Y_{j+1}\leq C_{0}\frac{\delta^{1-q}}{(\tau-\rho)^{d+q}}b^{j}Y_{j}^{1+\beta},$

where $b=2^{d+3q}>1$, $\beta=\sigma/d>0$, $C_{0}=C_{0}(d,s_{0},q,c)>1$ and

$Y_{j}=\frac{1}{f(k/R^{s})}\fint_{B_{j}}f\left(\frac{w_{j}(x)}{R^{s}}\right)\,\mathrm{d}x.$

Let us take

$k=R^{s}f^{-1}\left(C_{1}\frac{\delta^{(1-q)2d/s_{0}}}{(\tau-\rho)^{\gamma}}\fint_{B_{\tau
R}}f\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x\right)+k_{1},$

where $\gamma=2d(d+q)/s_{0}$, $C_{1}=C_{0}^{2d/s_{0}}b^{4d^{2}/s_{0}^{2}}$ and
$k_{1}$ is as before. This choice provides

$Y_{0}\leq\left(\frac{C_{0}\delta^{1-q}}{(\tau-\rho)^{d+q}}\right)^{-1/\beta}b^{-1/\beta^{2}},$

where we used that $\sigma\geq s_{0}/2$. Therefore, [LU68, Lemma 4.7] shows
that $Y_{j}\to 0$ as $j\to\infty$, which concludes that $u\leq k$ a.e. in
$B_{\rho R}$. By monotonicity of $f$, it follows that $f(u/R^{s})\leq
f(k/R^{s})$ a.e. in $B_{\rho R}$, which implies the desired result due to
Lemma 2.1 and since $f(a+b)\leq 2^{q}(f(a)+f(b))$. ∎

The following result includes Theorem 3.1 as the special case $\varepsilon=1$.
It is a direct consequence of Theorem 3.1 and a classical iteration argument.

######  Corollary 3.3.

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $0<s_{0}\leq s<1$,
$q>1$, $c>0$, $\varepsilon\in(0,1]$ and assume that $f$ satisfies ($f^{q}$).
There exists a constant $C>0$ such that if $u\in G_{+}(\Omega;q,c,s,f)$, then
for any $B_{R}(x_{0})\subset\Omega$ and $\delta\in(0,1)$, it holds that

(3.8)
$f^{\varepsilon}\left(\sup_{B_{R/2}(x_{0})}\frac{u_{+}}{R^{s}}\right)\leq
C\delta^{-\mu}\fint_{B_{R}(x_{0})}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x+\delta
f^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};x_{0},R/2)}{(R/2)^{s}}\right),$

where $\mu=2d(q-1)/(\varepsilon s_{0})$. The constant $C$ depends only on $d$,
$s_{0}$, $q$, $c$ and $\varepsilon$.

###### Proof.

We may assume that $x_{0}=0$. Let $1/2\leq\rho<\tau\leq 1$,
$\delta_{0}\in(0,1)$. By applying Theorem 3.1 with $\delta_{0}$, we have

(3.9) $\begin{split}g(\rho)&\leq
C\frac{\delta_{0}^{(1-q)2d/s_{0}}}{(\tau-\rho)^{\gamma}}\fint_{B_{\tau
R}}f\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x+\delta_{0}f\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)}{(R/2)^{s}}\right)=:I_{1}+I_{2},\end{split}$

where $C=C(d,s_{0},q,c)>0$, $\gamma=2d(d+q)/s_{0}$ and

$g(\rho)=f\left(\sup_{B_{\rho R}}\frac{u_{+}}{R^{s}}\right).$

Using Young’s inequality, we obtain

(3.10) $\begin{split}I_{1}&\leq
C\frac{\delta_{0}^{(1-q)2d/s_{0}}}{(\tau-\rho)^{\gamma}}g(\tau)^{1-\varepsilon}\fint_{B_{\tau
R}}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x\\\
&\leq\frac{1}{2}g(\tau)+C\frac{\delta_{0}^{-\mu}}{(\tau-\rho)^{\gamma/\varepsilon}}\left(\fint_{B_{R}}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x\right)^{1/\varepsilon}\end{split}$

for some $C=C(d,s_{0},q,c,\varepsilon)>0$, where $\mu=2d(q-1)/(\varepsilon
s_{0})$. Combining (3.9) and (3.10):

$g(\rho)\leq\frac{1}{2}g(\tau)+C\frac{\delta_{0}^{-\mu}}{(\tau-\rho)^{\gamma/\varepsilon}}\left(\fint_{B_{R}}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x\right)^{1/\varepsilon}+\delta_{0}f\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)}{(R/2)^{s}}\right)$

for any $1/2\leq\rho<\tau\leq 1$. Therefore, by an iteration lemma, see [GG82,
Lemma 1.1]:

$g(1/2)\leq
C\delta_{0}^{-\mu}\left(\fint_{B_{R}}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x\right)^{1/\varepsilon}+C\delta_{0}f\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)}{(R/2)^{s}}\right).$

Using the inequality $(a+b)^{\varepsilon}\leq
a^{\varepsilon}+b^{\varepsilon}$, we obtain

(3.11) $f^{\varepsilon}\left(\sup_{B_{R/2}}\frac{u_{+}}{R^{s}}\right)\leq
C\delta_{0}^{-\varepsilon\mu}\fint_{B_{R}}f^{\varepsilon}\left(\frac{u_{+}(x)}{R^{s}}\right)\,\mathrm{d}x+C\delta_{0}^{\varepsilon}f^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)}{(R/2)^{s}}\right),$

where $C=C(d,s_{0},q,c,\varepsilon)>1$. For a given $\delta\in(0,1)$, the
inequality (3.8) follows from (3.11) by setting
$\delta_{0}=(\delta/C)^{1/\varepsilon}\in(0,1)$. ∎

##  4\. Weak Harnack inequality

The goal of this section is to prove a weak Harnack inequality for functions
$u\in G_{-}(\Omega;q,c,s,f)$. There exist several possible estimates in the
literature, which go under the name “weak Harnack inequality”. They all differ
in the aspect that $\inf u$ is estimated by different Lebesgue-norms of $u$.
We will prove an estimate of the following type since it allows us to deduce a
full Harnack inequality by combination with Corollary 3.3.

######  Theorem 4.1.

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $0<s_{0}\leq s<1$,
$1<p\leq q$, $c>0$ and assume that $f$ satisfies (1.1). There exist constants
$C>0$ and $\varepsilon\in(0,1)$ such that if $u\in G_{-}(\Omega;q,c,s,f)$ is
nonnegative in $B_{R}(x_{0})\subset\Omega$, then

$\fint_{B_{R}(x_{0})}f^{\varepsilon}\left(\frac{u(x)}{R^{s}}\right)\,\mathrm{d}x\leq
Cf^{\varepsilon}\left(\inf_{B_{R/2}(x_{0})}\frac{u}{R^{s}}\right)+Cf^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{-};x_{0},R)}{R^{s}}\right).$

The constants $C$ and $\varepsilon$ depend only on $d$, $s_{0}$, $q$ and $c$.

Before we give the proof of Theorem 4.1, we recall the following growth lemma
from [CKW21]:

######  Lemma 4.2.

[CKW21, Theorem 4.1] Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let
$1<p\leq q$, $c,H>0$, $R>0$, $s_{0}\in(0,1)$ and assume $s\in[s_{0},1)$.
Assume that $f$ satisfies (1.1). Suppose that
$B_{4R}=B_{4R}(x_{0})\subset\Omega$. Let $u\in G_{-}(\Omega;q,c,s,g)$ satisfy
$u\geq 0$ in $B_{4R}$ and

$|B_{2R}\cap\\{u\geq H\\}|\geq\gamma|B_{2R}|$

for some $\gamma\in(0,1)$. There exists $\delta\in(0,1)$ such that if

$\mathrm{Tail}_{f^{\prime}}(u_{-};x_{0},4R)\leq\delta H,$

then

$u\geq\delta H\quad\text{in}~{}B_{R}.$

The constant $\delta$ depends only on $d$, $s_{0}$, $p$, $q$, $c$ and
$\gamma$.

###### Proof of Theorem 4.1.

Without loss of generality, we assume that $x_{0}=0$. Let us define

$\displaystyle L:=\inf_{B_{R/2}}u+\mathrm{Tail}_{f^{\prime}}(u_{-};0,8R).$

First of all, we claim that for any $H>0$ it holds:

(4.1) $\displaystyle\frac{|A^{+}_{t,R}|}{|B_{R}|}\leq\left(\frac{L}{\delta
t}\right)^{a}.$

Here, $a=\frac{\log\frac{1}{2}}{\log\delta}$, where $\delta\in(0,1)$ is the
constant from Lemma 4.2 applied with $\gamma=\frac{6^{-d}}{2}$ and $H:=t$. The
proof of (4.1) is a well-known consequence of Lemma 4.2 and a covering
argument due to Krylov and Safonov. It is explained in detail in [Coz17, Lemma
6.7 and (6.41)] and can be adapted to our setup without any changes being
necessary.
Let us explain how to deduce the desired result from (4.1). We choose
$\varepsilon=\frac{1}{2}\min(1,\frac{a}{q})$ and compute by Cavalieri’s
principle and performing a change of variables

$\displaystyle\fint_{B_{R}}f^{\varepsilon}\left(\frac{u(x)}{R^{s}}\right)\,\mathrm{d}x$
$\displaystyle=\varepsilon\int_{0}^{\infty}\frac{|B_{R}\cap\\{f(u/R^{s})\geq
t\\}|}{|B_{R}|}t^{\varepsilon-1}\textnormal{d}t=\varepsilon\int_{0}^{\infty}\frac{|A^{+}_{tR^{s},R}|}{|B_{R}|}f^{\varepsilon-1}(t)f^{\prime}(t)\textnormal{d}t$
$\displaystyle\leq\varepsilon\int_{0}^{L/R^{s}}f^{\varepsilon-1}(t)f^{\prime}(t)\textnormal{d}t+\varepsilon\int_{L/R^{s}}^{\infty}\left(\frac{L}{\delta
tR^{s}}\right)^{a}f^{\varepsilon-1}(t)f^{\prime}(t)\textnormal{d}t$
$\displaystyle=:I_{1}+I_{2}.$

For $I_{1}$, by a change of variables,

(4.2) $\displaystyle
I_{1}=\varepsilon\int_{0}^{f(L/R^{s})}t^{\varepsilon-1}\textnormal{d}t=f^{\varepsilon}\left(\frac{L}{R^{s}}\right).$

For $I_{2}$, we apply ($f^{q}$) and obtain

$\displaystyle I_{2}\leq\varepsilon
q\delta^{-a}\left(\frac{L}{R^{s}}\right)^{a}\int_{L/R^{s}}^{\infty}t^{-1-a}f^{\varepsilon}(t)\textnormal{d}t.$

From integration by parts and again ($f^{q}$), we see that

$\displaystyle\int_{L/R^{s}}^{\infty}t^{-1-a}f^{\varepsilon}(t)\textnormal{d}t$
$\displaystyle=\frac{1}{a}\left(\frac{L}{R^{s}}\right)^{-a}f^{\varepsilon}\left(\frac{L}{R^{s}}\right)-\frac{1}{a}\lim_{t\to\infty}t^{-a}f^{\varepsilon}(t)+\frac{\varepsilon}{a}\int_{L/R^{s}}^{\infty}t^{-a}f^{\varepsilon-1}(t)f^{\prime}(t)\textnormal{d}t$
$\displaystyle\leq\frac{1}{a}\left(\frac{L}{R^{s}}\right)^{-a}f^{\varepsilon}\left(\frac{L}{R^{s}}\right)+\frac{1}{2}\int_{L/R^{s}}^{\infty}t^{-1-a}f^{\varepsilon}(t)\textnormal{d}t,$

where we used the definition of $\varepsilon$ and Lemma 2.1(ii) in the last
step. It follows that

(4.3) $\displaystyle
I_{2}\leq\frac{2q}{a}\delta^{-a}f^{\varepsilon}\left(\frac{L}{R^{s}}\right),$

which yields, upon combining (4.2) and (4.3):

$\displaystyle\fint_{B_{R}}f^{\varepsilon}\left(\frac{u(x)}{R^{s}}\right)\,\mathrm{d}x\leq
Cf^{\varepsilon}\left(\frac{L}{R^{s}}\right)\leq
Cf^{\varepsilon}\left(\inf_{B_{R/2}}\frac{u}{R^{s}}\right)+Cf^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{-};0,8R)}{R^{s}}\right),$

where we used that $f(a+b)\leq 2^{q}(f(a)+f(b))$ and $(a+b)^{\varepsilon}\leq
a^{\varepsilon}+b^{\varepsilon}$. The desired result follows by noticing that

$\displaystyle\mathrm{Tail}_{f^{\prime}}(u_{-};0,8R)\leq
C\mathrm{Tail}_{f^{\prime}}(u_{-};0,R),$

which is a direct consequence of Lemma 2.7 applied with $\lambda=8^{s}$. ∎

##  5\. Harnack inequality

In this section, we prove our main result Theorem 1.1. First, we prove the
following estimate for $\mathrm{Tail}_{f^{\prime}}(u_{+};x_{0},R)$.

######  Lemma 5.1.

Let $\Omega$ be an open subset in $\mathbb{R}^{d}$. Let $0<s_{0}\leq s<1$,
$1<p\leq q$, $c>0$ and assume that $f$ satisfies (1.1). There exists a
constant $C>0$ such that if $u\in G_{-}(\Omega;q,c,s,f)$ is nonnegative in
$B_{R}(x_{0})\subset\Omega$, then

$\mathrm{Tail}_{f^{\prime}}(u_{+};x_{0},R/2)\leq
C\left(\sup_{B_{R/2}(x_{0})}u+\mathrm{Tail}_{f^{\prime}}(u_{-};x_{0},R/2)\right).$

The constant $C$ depends on $d$, $s_{0}$, $p$, $q$ and $c$.

###### Proof.

Without loss of generality, we may assume $x_{0}=0$. Let $w=u-2M$, where
$M=\sup_{B_{R/2}}u$. By $u\in G_{-}(\Omega;q,c,s,f)$:

(5.1)
$\begin{split}(1-s)&\int_{B_{R/4}}w_{-}(x)\left(\int_{\mathbb{R}^{d}}f^{\prime}\left(\frac{w_{+}(y)}{|x-y|^{s}}\right)|x-y|^{-d-s}\,\mathrm{d}y\right)\,\mathrm{d}x\\\
&\leq
c\int_{B_{R/2}}f\left(\frac{w_{-}(y)}{R^{s}}\right)\,\mathrm{d}y+c(1-s)\|w_{-}\|_{L^{1}(B_{R/2})}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{w_{-}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y.\end{split}$

Note that due to ($f_{p}$) it holds $f^{\prime}(0)=0$. This allows us to
consider the integral over $\mathbb{R}^{d}$ for the term on the left-hand
side. Since $|x-y|\leq 2|y|$ for every $x\in B_{R/4}$, $y\in B_{R/2}^{c}$ and
by (2.3), we estimate the first term from below by

$\displaystyle
c(1-s)\int_{B_{R/4}}w_{-}(x)\left(\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{w_{+}(y)}{|x-y|^{s}}\right)|x-y|^{-d-s}\,\mathrm{d}y\right)\,\mathrm{d}x$
$\displaystyle\geq
C(1-s)MR^{d}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{u_{+}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y-C(1-s)MR^{d}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{M}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y.$

Note that by monotonicity of $f^{\prime}$ and (2.1)

$\displaystyle(1-s)MR^{d}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{M}{|y|^{s}}\right)|y|^{-s-d}\,\mathrm{d}y\leq
CMR^{d-s}f^{\prime}\left(\frac{M}{R^{s}}\right).$

Furthermore, the terms on the right-hand side of (5.1) can be estimated from
above by:

$CMR^{d-s}f^{\prime}\left(\frac{M}{R^{s}}\right)+C(1-s)MR^{d}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{u_{-}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y$

using similar arguments. Altogether, we obtain

$\displaystyle(1-s)\left(\frac{R}{2}\right)^{s}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{u_{+}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y$
$\displaystyle\leq
Cf^{\prime}\left(\frac{M}{R^{s}}\right)+C(1-s)\left(\frac{R}{2}\right)^{s}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{u_{-}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y$
$\displaystyle\leq\frac{1}{2}\left[f^{\prime}\left(\frac{CM}{R^{s}}\right)+C(1-s)\left(\frac{R}{2}\right)^{s}\int_{B_{R/2}^{c}}f^{\prime}\left(\frac{u_{-}(y)}{|y|^{s}}\right)|y|^{-d-s}\,\mathrm{d}y\right],$

where we used (2.1) in the last step. Next, we apply $(f^{\prime})^{-1}$ on
both sides of the estimate and multiply with $(R/2)^{s}$ to obtain

$\displaystyle\mathrm{Tail}_{f^{\prime}}(u_{+};0,R/2)\leq
CM+C\mathrm{Tail}_{f^{\prime}}(u_{-};0,R/2),$

where we applied Lemma 2.6, (2.6) and Lemma 2.7. ∎

We are now ready to give the proof of Theorem 1.1.

###### Proof of Theorem 1.1.

We may assume that $x_{0}=0$. Let $\varepsilon\in(0,1)$ be the constant from
Theorem 4.1. By Corollary 3.3 and Lemma 5.1, we have

$f^{\varepsilon}\left(\sup_{B_{R/2}}\frac{u}{R^{s}}\right)\leq
C\delta^{-\mu}\fint_{B_{R}}f^{\varepsilon}\left(\frac{u(x)}{R^{s}}\right)\,\mathrm{d}x+\delta
f^{\varepsilon}\left(C\frac{\sup_{B_{R/2}}u+\mathrm{Tail}_{f^{\prime}}(u_{-};0,R/2)}{(R/2)^{s}}\right),$

where $\mu=2d(q-1)/(\varepsilon s_{0})$. Using ($f^{q}$), Lemma 2.1 and $u\geq
0$ in $B_{R}$, we obtain

$\delta
f^{\varepsilon}\left(C\frac{\sup_{B_{R/2}}u+\mathrm{Tail}_{f^{\prime}}(u_{-};0,R/2)}{(R/2)^{s}}\right)\leq
C\delta\left(f^{\varepsilon}\left(\sup_{B_{R/2}}\frac{u}{R^{s}}\right)+f^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{-};0,R)}{R^{s}}\right)\right).$

By taking $\delta$ sufficiently small so that $C\delta<1/2$, we have

$f^{\varepsilon}\left(\sup_{B_{R/2}}\frac{u}{R^{s}}\right)\leq
C\fint_{B_{R}}f^{\varepsilon}\left(\frac{u(x)}{R^{s}}\right)\,\mathrm{d}x+f^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{-};0,R)}{R^{s}}\right).$

Thus, it follows from Theorem 4.1 that

$f^{\varepsilon}\left(\sup_{B_{R/2}}\frac{u}{R^{s}}\right)\leq
Cf^{\varepsilon}\left(\inf_{B_{R/2}}\frac{u}{R^{s}}\right)+Cf^{\varepsilon}\left(\frac{\mathrm{Tail}_{f^{\prime}}(u_{-};0,R)}{R^{s}}\right).$

The desired inequality follows by using $(a+b)^{\frac{1}{\varepsilon}}\leq
2^{\frac{1}{\varepsilon}-1}(a^{\frac{1}{\varepsilon}}+b^{\frac{1}{\varepsilon}})$,
as well as the estimate $f(a)+f(b)\leq 2f(a+b)$ and Lemma 2.4. ∎

## References

  * [ACPS21] Angela Alberico, Andrea Cianchi, Luboš Pick, and Lenka Slavíková.  On fractional Orlicz-Sobolev spaces.  Anal. Math. Phys., 11(2):Paper No. 84, 21, 2021.
  * [BBCK09] Martin T. Barlow, Richard F. Bass, Zhen-Qing Chen, and Moritz Kassmann.  Non-local Dirichlet forms and symmetric jump processes.  Trans. Amer. Math. Soc., 361(4):1963–1999, 2009.
  * [BK05] Richard F. Bass and Moritz Kassmann.  Harnack inequalities for non-local operators of variable order.  Trans. Amer. Math. Soc., 357(2):837–850, 2005.
  * [BKO21] Sun-Sig Byun, Hyojin Kim, and Jihoon Ok.  Local Hölder continuity for fractional nonlocal equations with general growth.  arXiv:2112.13958, 2021.
  * [BL02] Richard F. Bass and David A. Levin.  Harnack inequalities for jump processes.  Potential Anal., 17(4):375–388, 2002.
  * [BOS21] Sun-Sig Byun, Jihoon Ok, and Kyeong Song.  Hölder regularity for weak solutions to nonlocal double phase problems.  arXiv:2108.09623, 2021.
  * [BS05] Krzysztof Bogdan and Paweł Sztonyk.  Harnack’s inequality for stable Lévy processes.  Potential Anal., 22(2):133–150, 2005.
  * [BSV20] Julián Fernández Bonder, Ariel Salort, and Hernán Vivas.  Interior and up to the boundary regularity for the fractional $g$-Laplacian: the convex case.  arXiv:2008.05543, 2020.
  * [CK03] Zhen-Qing Chen and Takashi Kumagai.  Heat kernel estimates for stable-like processes on $d$-sets.  Stochastic Process. Appl., 108(1):27–62, 2003.
  * [CK21] Jamil Chaker and Minhyun Kim.  Local regularity for nonlocal equations with variable exponents.  arXiv:2107.06043, 2021.
  * [CKW19] Zhen-Qing Chen, Takashi Kumagai, and Jian Wang.  Elliptic Harnack inequalities for symmetric non-local Dirichlet forms.  J. Math. Pures Appl. (9), 125:1–42, 2019.
  * [CKW21] Jamil Chaker, Minhyun Kim, and Marvin Weidner.  Regularity for nonlocal problems with non-standard growth.  arXiv:2111.09182, 2021.
  * [Coz17] Matteo Cozzi.  Regularity results and Harnack inequalities for minimizers and solutions of nonlocal problems: a unified approach via fractional De Giorgi classes.  J. Funct. Anal., 272(11):4762–4837, 2017.
  * [CS07] Luis Caffarelli and Luis Silvestre.  An extension problem related to the fractional Laplacian.  Comm. Partial Differential Equations, 32(7-9):1245–1260, 2007.
  * [CS09] Luis Caffarelli and Luis Silvestre.  Regularity theory for fully nonlinear integro-differential equations.  Comm. Pure Appl. Math., 62(5):597–638, 2009.
  * [DCKP14] Agnese Di Castro, Tuomo Kuusi, and Giampiero Palatucci.  Nonlocal Harnack inequalities.  J. Funct. Anal., 267(6):1807–1836, 2014.
  * [DFP19] Cristiana De Filippis and Giampiero Palatucci.  Hölder regularity for nonlocal double phase equations.  J. Differential Equations, 267(1):547–586, 2019.
  * [FBS19] Julián Fernández Bonder and Ariel M. Salort.  Fractional order Orlicz-Sobolev spaces.  J. Funct. Anal., 277(2):333–367, 2019.
  * [FZ21] Yuzhou Fang and Chao Zhang.  On weak and viscosity solutions of nonlocal double phase equations.  Int. Math. Res. Not., 2021.
  * [FZ22] Yuzhou Fang and Chao Zhang.  Harnack inequality for the nonlocal equations with general growth.  arXiv:2201.09495, 2022.
  * [GG82] Mariano Giaquinta and Enrico Giusti.  On the regularity of the minima of variational integrals.  Acta Math., 148:31–46, 1982.
  * [GKS20] Divya Goel, Deepak Kumar, and Konijeti Sreenadh.  Regularity and multiplicity results for fractional $(p,q)$-Laplacian equations.  Commun. Contemp. Math., 22(8):1950065, 37, 2020.
  * [GKS21] Jacques Giacomoni, Deepak Kumar, and Konijeti Sreenadh.  Interior and boundary regularity results for strongly nonhomogeneous $p$, $q$-fractional problems.  Adv. Calc. Var., 2021.
  * [GKS22] Jacques Giacomoni, Deepak Kumar, and Konijeti Sreenadh.  Global regularity results for non-homogeneous growth fractional problems.  J. Geom. Anal., 32(1):Paper No. 36, 41, 2022.
  * [Kas07a] Moritz Kassmann.  The classical Harnack inequality fails for non-local operators.  SFB 611 - Preprint, No. 360., 2007.
  * [Kas07b] Moritz Kassmann.  Harnack inequalities: an introduction.  Bound. Value Probl., pages Art. ID 81415, 21, 2007.
  * [Kas11] Moritz Kassmann.  A new formulation of Harnack’s inequality for nonlocal operators.  C. R. Math. Acad. Sci. Paris, 349(11-12):637–640, 2011.
  * [LU68] Olga A. Ladyzhenskaya and Nina N. Ural’tseva.  Linear and quasilinear elliptic equations.  Academic Press, New York-London, 1968.  Translated from the Russian by Scripta Technica, Inc, Translation editor: Leon Ehrenpreis.
  * [MP96] Elvira Mascolo and Gloria Papi.  Harnack inequality for minimizers of integral functionals with general growth conditions.  NoDEA Nonlinear Differential Equations Appl., 3(2):231–244, 1996\.
  * [Ok21] Jihoon Ok.  Local Hölder regularity for nonlocal equations with variable powers.  arXiv:2107.06611, 2021.
  * [Str19] Martin Strömqvist.  Harnack’s inequality for parabolic nonlocal equations.  Ann. Inst. H. Poincaré Anal. Non Linéaire, 36(6):1709–1745, 2019.
  * [SV04] Renming Song and Zoran Vondraček.  Harnack inequality for some classes of Markov processes.  Math. Z., 246(1-2):177–202, 2004.
  * [SZ13] Pablo Raúl Stinga and Chao Zhang.  Harnack’s inequality for fractional nonlocal equations.  Discrete Contin. Dyn. Syst., 33(7):3153–3170, 2013.


# On Using Deep Learning Proxies as Forward Models in Optimization Problems

Fatima Albreiki, Nidhal Belayouni and Deepak K. Gupta
AIQ, UAE
<EMAIL_ADDRESS><EMAIL_ADDRESS>
Work done as part of the internship at AIQ.

###### Abstract

Physics-based optimization problems are generally very time-consuming,
especially due to the computational complexity associated with the forward
model. Recent works have demonstrated that physics-modelling can be
approximated with neural networks. However, there is always a certain degree
of error associated with this learning, and we study this aspect in this
paper. We demonstrate through experiments on popular mathematical benchmarks,
that neural network approximations (NN-proxies) of such functions when plugged
into the optimization framework, can lead to erroneous results. In particular,
we study the behaviour of particle swarm optimization and genetic algorithm
methods and analyze their stability when coupled with NN-proxies. The
correctness of the approximate model depends on the extent of sampling
conducted in the parameter space, and through numerical experiments, we
demonstrate that caution needs to be taken when constructing this landscape
with neural networks. Further, the NN-proxies are hard to train for higher
dimensional functions, and we present our insights for 4D and 10D problems.
The error is higher for such cases, and we demonstrate that it is sensitive to
the choice of the sampling scheme used to build the NN-proxy. The code is
available at https://github.com/Fa-ti-ma/NN-proxy-in-optimization.

##  1 Introduction

In the past few years, deep learning has led to several advancements in the
different disciplines of science and engineering. Some popular examples
include the development of a rigorous solution to the protein-fold problem [5]
and enhancing computational fluid dynamics simulations [9], among others. One
popular application of deep learning is to develop proxy models that can serve
as an alternative for the computationally-intensive physics-based models. An
interesting example is the deep-learning tomography method [2] which aims at
by-passing the computationally demanding steps of the traditional seismic
tomography processes.

The gain with deep learning can be even higher when such a proxy model is used
as an alternative to a physics-based simulation that serves as a forward model
of an optimization problem. Note that the optimization problem referred is
different from the optimization problem solved for model training in deep
learning. We refer here to the physics-based optimization problems where the
goal is to optimize the parameters associated with the underlying physics of a
certain problem. An example is the problem of near-surface velocity modeling
of the earth, where the goal of the optimization is to identify the
distribution of compression and shear velocities for a predefined model space,
such that the resultant spectrum images of the earth, as acquired by sensors,
can be optimally reconstructed [10]. This optimization involves solving the
forward model and computing its gradient repeatedly for a number of
iterations, and reducing the cost associated with the forward model can speed
up the overall optimization problem by a large margin.

When using deep learning proxy models, it is important to note that the
parameter space of the proxy model might not be the same as the original
forward model. While the physics-model is parameterized with a few control
variables, the proxy NN-model is generally designed using millions of weights.
The choice of parameters for the NN-proxy depends on several factors, such as
the sampling method for picking the training data, extent of convergence of
the model, _etc_. It is well known that the performance of the optimization
methods depends significantly on the complexity of the parameter landscape,
and a changed landscape due to the use of deep learning can also lead to
completely different optimized solutions which are not necessarily optimal.
Clearly, there are two different aspects that need to be studied when creating
such a proxy model for the optimization process: the robustness of the deep
learning approach as well as efficacy of the chosen optimization method.

Contributions. In this work, we present our first observations on how the
choice of data sampling as well as the amount of data samples for training the
NN-proxy can influence the performance of the overall optimization process. We
use the popular benchmark functions from the field of global optimization
(_Rosenbrock, Rastrigin, Ackley,_etc.__) as alternative representations for
the physics models used in NN-proxies. These mathematical functions have long
served as benchmarks for the generic field of optimization, and the exact
solutions are known. For optimization, we use particle swarm optimization
(PSO) [6] method and genetic algorithm (GA) [4] and start with first analyzing
the stability of these methods when the original parameter landscapes are
replaced with those obtained from NN-proxies. Our investigation reveals that
the chosen optimization methods become very sensitive to their respective
initialization schemes when the NN-proxy is employed. Sampling scheme as well
as the amount of training data plays an important role in the construction of
good-NN proxies, and we study this aspect as well qualitative as well as
quantitative assessment of the error in optimization. Lastly, we also study
the influence of NN-proxy error on optimization in higher dimensions and
report results for 4D and 10D functions.

##  2 Deep Learning Proxy Models (NN-Proxies)

Building a deep learning proxy model, also referred as NN-proxy, implies
learning the mapping of the physics-based models in a data-driven manner and
using it as a forward model in an optimization problem. A better understanding
of this can be obtained from Figure 1. Similar to the physics-model, the NN-
proxy computes the desired output which is then evaluated using the objective
function to decide if the model input parameters are to be updated or not. It
is evident from the figure that a large error in the predictions made by the
NN-proxy can drift the whole optimization process off track, and it is very
important that the neural network model is sufficiently good before it can be
plugged into the optimization process as a NN-proxy. We analyze this aspect
through multiple experiments on popular mathematical benchmarks below.

Experimental setup. There are several factors that need to be understood when
using a NN-proxy in a optimization process, and in this paper, we study a few
important ones. In place of the original physics-models, we use popular
mathematical benchmarks to study the efficacy of global optimization methods,
and use neural networks to learn proxy representations of these functions. We
use Rosenbrock, Rastrigin and Ackley functions [8, 7, 1] and experiment with
particle swarm optimization (PSO) [6] and genetic algorithm (GA) [4] as the
two optimization methods. The two chosen optimization methods are well
established to solve complex multimodal global optimization methods, and we
experiment here how they fair when the original models are replaced by their
NN-proxies. We also consider how the error associated with the NN-proxy scales
in higher dimensions, and for this purpose, we consider the cases of 4D and
10D. More details related to the chosen benchmark functions are described in
Appendix A. Further, brief descriptions of the chosen optimization methods are
presented in Appendix B.

Below we present our observations related to a series of experiments aimed at
understanding the efficacy of NN-proxies.

Figure 1: Schematic representation of a conventional optimization pipeline
using conventional physics-based forward models as well as their alternative
representation using NN-proxies.

Choice of the initial model for optimization. Optimization methods are
generally sensitive to the choice of the initial values of the optimization
parameters. Thus, it is of interest to investigate how this behaviour changes
when the NN-proxies replace the actual forward models. Figures 2 and 3 show
the results obtained for Rosenbrock and Rstrigin functions, respectively,
using PSO and GA. For both the functions, we present the solution of PSO and
GA on the actual function, and it is observed that for almost all cases, the
global minima is found. On the contrary, when the true function is replaced by
its NN-proxy, both the optimization methods become very sensitive. For the
case where dense sampling is used for creating the training dataset, the
sensitivity to the initial values is low. However, for sparse sampling, the
error is significantly high with the predicted solutions diverging
significantly from the global minimum. Experimental details are described in
Appendix C and the quantitative results related to this set of experiments as
well as for Ackley function are presented in Table 1. From this experiment, it
can be concluded that the behavior of an optimizer with the actual physics
model can be very different from that when used with a NN-proxy. Moreover,
most physics problems are complex and computationally expensive, and dense
sampling is not the choice, and with sparse sampling, the results can be very
erroneous. This is a matter of concern and deserves attention when using NN-
proxies to represent physics of a system.

Figure 2: Example results for the Rosenbrock function for dense (left) and
sparse (right) sampling of the function space. For each case, we show PSO and
GA results on the true function (left), results of PSO on the NN-proxy
(middle) and results of GA o the NN-proxy (right) for 5 random seeds.

Figure 3: Example results for the Rastrigin function for dense (left) and
sparse (right) sampling of the function space. For each case, we show PSO and
GA results on the true function (left), results of PSO on the NN-proxy
(middle) and results of GA o the NN-proxy (right) for 5 random seeds.

Effect of data sampling. Data sampling plays an important role in bulding a
good NN-proxy. This is clearly reflected from the results shown in Figures 2
and 3 for 2D functions. We further study this and present a quantitative
analysis in Table 1. Apart from the uniform dense and sparse samplings, we
also include Gaussian sampling around the global minimum. This choice is meant
to weakly reflect the scenario when a priori information is available on the
location of the best solution. From the obtained results, we see that for all
cases of sampling, there is a significant amount of error in the prediction
when using the NN-proxy. In general, we see that sparse sampling is not very
effective and the optimizers are very sensitive to the landscape constructed
from such data. For Rosenbrock, GA seems to be unstable even for the true
function. In terms of comparing the two optimization methods, both seem to
work better than each other for different scenarios and equally good overall.
Out of the three functions, we see that the results of PSO and GA are very
good on the NN-proxy of Ackley function. The reason could be that the globally
optimal solution of this function differs significantly from the locally
optimal solutions and the NN-proxy can represent it very well.

Table 1: Mean Euclidean distance and the standard deviation between the
globally optimal point of the three mathematical benchmarks in 2D and the
respective solutions obtained for the NN-proxies using PSO and GA. Three
different sampling methods are used, and ground-truth refers to the case where
PSO and GA are used on the true functions.

Function Rosenbrock Rastrigin Ackley Optimizer PS GA PS GA PS GA Ground-Truth
$0.0\pm 0.0$ $0.91\pm 0.53$ $0.0\pm 0.0$ $0.0\pm 0.0$ 0.0 $0.0\pm 0.0$ $0.0\pm
0.0$ Dense $0.34\pm 0.21$ $0.54\pm 0.69$ $0.75\pm 0.42$ $0.16\pm 0.22$
$0.03\pm 0.0$ $0.03\pm 0.0$ Sparse $1.17\pm 0.62$ $1.37\pm 0.49$ $1.41\pm
0.04$ $1.24\pm 0.13$ $0.14\pm 0.08$ $0.14\pm 0.05$ Gaussian $1.07\pm 0.54$
$1.15\pm 0.28$ $1.07\pm 0.54$ $1.15\pm 0.28$ $0.06\pm 0.03$ $0.07\pm 0.02$

Table 2: Mean Euclidean distance and the standard deviation between the
globally optimal point of the three mathematical benchmarks in 4D and the
respective solutions obtained for the NN-proxies using PSO and GA. Here, three
different sampling methods are used t build the training set for the NN-
proxies.

Function Rosenbrock Rastrigin Ackley Optimizer PS GA PS GA PS GA Dense
$0.60\pm 0.03$ $1.59\pm 0.74$ $0.47\pm 0.0$ $0.44\pm 0.02$ $0.76\pm 0.20$
$0.90\pm 0.20$ Sparse $0.86\pm 0.06$ $1.56\pm 0.72$ $0.72\pm 0.0$ $0.70\pm
0.02$ $0.37\pm 0.03$ $0.27\pm 0.07$ Gaussian $2.51\pm 0.57$ $2.49\pm 0.85$
$0.71\pm 0.0$ $0.67\pm 0.05$ $0.26\pm 0.07$ $0.27\pm 0.06$

Table 3: Mean Euclidean distance and standard deviation between the true
optimal solutions of the three mathematical benchmarks in 10D and the
respective solutions obtained for the NN-proxies using PSO and GA.

Function Rosenbrock Rastrigin Ackley Optimizer PS GA PS GA PS GA Dense
$2.62\pm 0.46$ $2.95\pm 0.51$ $3.69\pm 0.53$ $4.45\pm 1.06$ $0.64\pm 0.01$
$0.70\pm 0.13$ Sparse $3.37\pm 0.74$ $3.71\pm 0.51$ $2.20\pm 0.29$ $2.36\pm
0.56$ $1.19\pm 0.0$ $1.05\pm 0.08$ Gaussian $5.31\pm 1.22$ $5.67\pm 1.53$
$6.25\pm 0.43$ $7.9\pm 2.6$ $0.30\pm 0.0$ $0.30\pm 0.07$

Stability in higher dimensions. For a better approximation of the actual
physics, the mathematical benchmarks need to be studied in higher dimensions.
In this regard, we study the performance of NN-proxy in 4D and 10D also, and
the related results are reported in Tables 2 and 3, respectively. With the
increased number of dimensions, the error grows as expected. An interesting
observation is that for higher dimensions, the Gaussian sampling scheme, which
is expected to better model the landscape arond the true solution, seems to
perform inferior compared to sparse sampling. The anticipated reason behind
this issue is that too many samples around the true solution limit the number
of samples that can be used to approximate the landscape far from it. This
means that the function cannot be constructed well for regions far from it.
This can lead to landscapes comprising false good solutions for the NN-proxy.
In this regard, a better choice seems to be to go with uniform sampling.

##  3 Conclusions

In this paper, we have studied how the outcome of a physics-based optimization
process can be adversely influenced when approximating the underlying physics
process with a neural network proxy. We have demonstrated through experiments
on popular mathematical benchmarks, that neural network approximations (NN-
proxies) of such functions can lead to erroneous results for different
training setups. The correctness of the approximate model depends on the
extent of sampling conducted in the parameter space, and caution needs to be
taken when constructing this landscape with neural networks. Further, the NN-
proxies are hard to train for higher dimensional functions, and through
experiments on 4D and 10D problems, we demonstrated that the error in the
optimized solution can be very high. The results reported in this paper are
our first observations on the effect of NN-proxies on the optimization
process. Clearly, we need to develop a rigorous benchmarking criterion to
evaluate any novel NN-proxy and understand its robustness and generalization
before it can be plugged in an optimization framework to replace a certain
physics process.

##  4 Broader Impact Statement

This work focuses on understanding the error associated with the use of deep
learning to replace the complex physics-based forward models used in
optimization processes. While the initial concept in the paper is demonstrated
on mathematical benchmarks, we hope that the observations and any
countermeasures will be transferable on real world optimization problems
involving physical simulations. Example of such problems includes proxy model
for eigen value analysis in a compliant mechanism problem of topology and
design optimization, among others. Overall, there are several problems where
the results of this paper and the follow up will be useful.

Further, we do not see any ethical concerns of negative societal impact of
this work.

## References

  * [1] D.H Ackley.  A Connectionist Machine for Genetic Hill climbing, volume SECS28 of The Kluwer International Series in Engineering and Computer Science. Kluwer Academic Publishers, Boston, 1987.
  * [2] Mauricio Araya-Polo, Joseph Jennings, Amir Adler, and Taylor Dahlke.  Deep-learning tomography.  The Leading Edge, 37(1):58–66, 2018.
  * [3] K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan.  A fast and elitist multiobjective genetic algorithm: Nsga-ii.  IEEE Transactions on Evolutionary Computation, 6(2):182–197, 2002\.
  * [4] John H. Holland.  Genetic algorithms.  Scientific American, 267(1):66–73, 1992.
  * [5] J. Jumper, R. Evans, and A. et al. Pritzel.  Highly accurate protein structure prediction with alphafold.  Nature, 596:583–589, 2021.
  * [6] J. Kennedy and R. Eberhart.  Particle swarm optimization.  In Proceedings of ICNN’95 - International Conference on Neural Networks, volume 4, pages 1942–1948, 1995.
  * [7] Leonard Andreevič Rastrigin.  Systems of extremal control.  Nauka, 1974.
  * [8] H. H. Rosenbrock.  An Automatic Method for Finding the Greatest or Least Value of a Function.  The Computer Journal, 3(3):175–184, 1960.
  * [9] R. Vinuesa and S. L. Brunton.  Enhancing computational fluid dynamics with machine learning.  Nature Computational Science, 2:358–366, 2022.
  * [10] P. Zwartjes, D. Gupta, and T. Gupta.  Near surface velocity estimation using surface waves and deep learning.  Proceedings, EAGE, pages 1–5, 2022.

##  Appendix A Benchmark Functions

We use popular mathematical benchmark functions from the field of global
optimization to train our NN-proxies. Our choice of functions are Rosenbrock,
Rastrigin and Ackley. Figure 4 shows the functions plots.

Rosenbrock [8] is usually evaluated on the hypercube
$x_{i}\in[-5,10],\enskip\forall\enskip i=1,…,d$. We have restricted to the
hypercube $x_{i}\in[-2.048,2.048],\enskip\forall\enskip i=1,…,d$ and use it in
our study for the sake of simplification. The mathematical formulation is
stated in Eq. 1. The global minimum for this function is at
$\mathbf{x}=[1,..,1]$ and $f(x_{i})=0$.

Rastrigin [7] is usually evaluated on the hypercube
$x_{i}\in[-5.12,5.12],\enskip\forall\enskip i=1,…,d$, and the mathematical
formulation for this function is stated in Eq. 2). The global minimum is at
$\mathbf{x}=[0,..,0]$ with $f(x_{i})=0$.

Ackley [1] is usually evaluated on the hypercube
$x_{i}\in[-32.768,32.768]\enskip\forall\enskip i=1,…,d$, it may be restricted
to a smaller domain. For the mathematical representation of this function, see
Eq. 3. The global minimum is at $x_{i}=[0,...,0]$ with $f(x_{i})=0$.

$f(x)=\sum_{i=1}^{d-1}[100(x_{i+1}-x_{i}^{2})^{2}+(x_{i}-1)^{2}]$ (1)

$f(x)=10d+\sum_{i=1}^{d}[x_{i}^{2}-10cos(2\pi x_{i})]$ (2)

$f(x)=-20exp(-0.2\sqrt{\frac{1}{d}\sum_{i=1}^{d}x_{i}^{2}})-exp(\frac{1}{d}\sum_{i=1}^{d}cos(2\pi
x_{i}))+20+e$ (3)

(a) Rosenbrock

(b) Rastrigin

(c) Ackley

Figure 4: Three benchmark functions: Rosenbrock, Rastrigin and Ackley.

##  Appendix B Optimization Methods

###  B.1 Particle Swarm Optimization

Particle Swarm Optimization was proposed by Kennedy and Eberhart [6]. It is
inspired by the idea that a school of fish or a flock of birds that moves in a
group ‘can profit from the experience of all other members’. A brief overview
of the underlying mechanism is described below.

Algorithm Details. Let there be $P$ particles, where we denote the position of
particle $i$ at iteration $t$ as $X^{i}(t)$. This position can then be stated
as $X^{i}(t)=(x^{i}(t),y^{i}(t))$ and velocity as
$V^{i}(t)=(v_{x}^{i}(t),v_{y}^{i}(t))$. At the next iteration, the position is
updated as: $X^{i}(t+1)=X^{i}(t)+V^{i}(t+1)$ and the velocity is updated as:

$V^{i}(t+1)=wV^{i}(t)+c_{1}r_{1}(pbest^{i}\textendash
X^{i}(t))+c_{2}r_{2}(gbest\textendash X^{i}(t))$ (4)

Where $r_{1}$ and $r_{2}$ are random numbers between 0 and 1. Constants $w$,
$c_{1}$ and $c_{2}$ are inertia, learning rate for individual ability and
social influence respectively. $pbest^{i}$ is the position that gives the best
cost explored by particle $i$ and $gbest$ is the best explored solution by all
swarms.

###  B.2 Genetic Algorithm

The genetic algorithm [4] is a search heuristic that is inspired by Darwin’s
theory of natural evolution in which the fittest of individuals are the ones
who survive. This fitness is measured by a fitness function. At every
iteration, individuals with high fitness have more chance to be selected for
reproduction to produce off springs for the next generation. The algorithm
terminates if the population has converged i.e, offspring are not different
from the previous generation or number of generations has been reached. For
details related to implementation, see for example [3].

##  Appendix C Experiments: Additional Details

In this section we provide details about our experimental settings including
the proxy models’ architecture and training details.

We employ three different sampling mechanisms. For the dense sampling
strategy, we use 10,000 samples for all cases to build the NN-proxies and
these samples are equally spaced in the hyperspace. For sparse case, we use
25% of the data that we use in the dense case. For the case of Gaussian
sampling, we use the same number of points as the dense case, but the
distribution follows a multidimensional Gaussian distribution.

In Dense case we use 10K points for training our models. In Sparse case we use
25% of the data. In Gaussian case, we use 10K centred around the functions’
global minimum.

For Rosenbrock proxy model in 2D case, we employ a 3 layer multilayer
perceptron (MLP) with hidden units of $\\{15,50,15\\}$ and RelU activation
function. We train for 100 epochs with learning rate of 0.001 and Adam
optimizer. In 4D case, we use a 4 layer MLP with hidden units of
$\\{15,50,15,10\\}$ and RelU activation. We train with the same number of
epochs and similar learning rate and optimizer as 2D case. In 10D case, we use
the same model setting as 4D case and train for 500 epochs with learning rate
of 0.001 and Adam optimizer. All experiments were run on a NVIDIA RTX 6000
graphics card.

For Rastrigin and Ackley proxy models in 2D, 4D and 10D cases, we use 6 layer
MLP with hidden units of $\\{20,50,120,70,20,10\\}$. We use Relu activation
and train for 500 epochs. Adam optimizer is used with learning rate of 0.001.


$\iota_{2}=\iota_{F^{1}_{1}}\iota_{F^{2}_{2}}-\iota_{F^{2}_{1}}\iota_{F^{1}_{2}}$
on $\sigma(\psi^{1},\psi^{2})$; a cancellation analogous to the one in (5.172)
and (5.2) occurs. Really, one has

(5.174)
$\iota_{2}\sigma(\psi^{1},\psi^{2})=\left(\iota_{F^{1}_{1}}\iota_{F^{2}_{2}}-\iota_{F^{2}_{1}}\iota_{F^{1}_{2}}\right)\eta(\psi^{1})-\frac{1}{4}\left(\iota_{F^{1}_{1}}\iota_{F^{2}_{2}}-\iota_{F^{2}_{1}}\iota_{F^{1}_{2}}\right)\wedge\sum_{\alpha,\beta}\psi^{2}_{\alpha}\psi^{2}_{\beta}\iota_{F^{1}_{\alpha}}\iota_{F^{1}_{\beta}}\eta(\psi^{1})+\ldots=0\
.$

The vanishing of the first term is trivial, because of the contractions along
the vectors $F^{2}_{\bullet}$, the higher terms vanish because of the minus
sign in $\iota_{2}$, as for (5.2). For example, the second term of (5.174)
reads

$-\frac{1}{4}\left(\iota_{F^{1}_{1}}\iota_{F^{2}_{2}}-\iota_{F^{2}_{1}}\iota_{F^{1}_{2}}\right)\left[\left(\psi^{2}_{1}\right)^{2}\iota^{2}_{F^{1}_{1}}+2\psi^{2}_{1}\psi^{2}_{2}\iota_{F^{1}_{1}}\iota_{F^{1}_{2}}+\left(\psi^{2}_{2}\right)^{2}\iota^{2}_{F^{1}_{2}}\right]\eta(\psi^{1})=$
(5.175)
$=-\frac{1}{2}\left[\psi^{2}_{1}\iota_{F^{1}_{1}}^{2}\iota_{F^{1}_{2}}+\psi^{2}_{2}\iota^{2}_{F^{1}_{2}}\iota_{F^{1}_{1}}-\psi^{2}_{1}\iota^{2}_{F^{1}_{1}}\iota_{F^{1}_{2}}-\psi^{2}_{2}\iota_{F^{1}_{1}}\iota_{F^{1}_{2}}^{2}\right]\eta(\psi^{1})=0\
,$

and analogously for all the higher terms. This confirms that actually there is
a single cohomology generator for pseudoforms in the coset
$\mathfrak{osp}(2|2)/\left(\mathfrak{so}(2)\oplus\mathfrak{sp}(2)\right)$.

The results of this section are non-trivial: we have shown that the explicit
distributional realisation of pseudoforms and integral forms is a powerful
tool to calculate cohomology representatives and that the explicit results are
consistent with those obtained without referring to any realisation through
spectral sequences. Nonetheless, we have been able to construct a pseudoform
that is invariant w.r.t. $\mathfrak{so}(2)\oplus\mathfrak{sp}(2)$; the
explicit realisation, supported by the abstract counterpart, could serve as a
starting point for the introduction of pseudoforms in more general contexts,
e.g., on supermanifolds: one could think to introduce pseudoforms as integral
forms of sub-supermanifolds respecting some isometries (see also [Man]). The
general setup for these definitions is once again suggested by the algebraic
setting analysed in this paper: relative (co)homology. This will be the
subject of future investigations.

## Acknowledgements

This work has been partially supported by Università del Piemonte Orientale
research funds. We thank S. Cacciatori, R. Catenacci and S. Noja for many
useful discussions. C.A.C. is grateful to B. Jurčo for many suggestions and
discussions.

## References

  * [AETW] A. Achucarro, J. M. Evans, P. K. Townsend and D. L. Wiltshire, _Super p-Branes_ , Phys. Lett. B 198 (1987), 441-446
  * [BS] M. A. Baranov and A. S. Schwarz, _Multiloop Contribution to String Theory_ , JETP Lett. 42 (1985), 419-421
  * [Bel] A. Belopolsky, _De Rham cohomology of the supermanifolds and superstring BRST cohomology_ , Phys. Lett. B 403 (1997), 47-50
  * [Bel2] A. Belopolsky, _New geometrical approach to superstrings_ , [arXiv:hep-th/9703183 [hep-th]]
  * [Bel3] A. Belopolsky, _Picture changing operators in supergeometry and superstring theory_ , [arXiv:hep-th/9706033 [hep-th]]
  * [BL] I. N. Bernstein, D. A. Leites _Integral forms and the Stokes formula on supermanifolds_ Funkt.Anal.Pril. 11 (1977) 55
  * [BKN] B. D. Boe, J. R. Kujawa, D. K. Nakano, _Cohomology and support varieties for Lie superalgebras_ , Trans. Amer. Math. Soc. 362 (2010), no. 12, 6551–6590
  * [CNR] S. L. Cacciatori, S. Noja and R. Re, _The Unifying Double Complex on Supermanifolds_ , Doc. Math. 27, 489-518 (2022)
  * [CDF] L. Castellani, R. D’Auria and P. Fre, _Supergravity and superstrings: A Geometric perspective. Vol. 1,2,3_ , Singapore: World Scientific (1991) 1375-2162
  * [CGN] R. Catenacci, P. A. Grassi and S. Noja, _Superstring Field Theory, Superforms and Supergeometry_ , J. Geom. Phys. 148 (2020), 103559
  * [CGN2] R. Catenacci, P. A. Grassi and S. Noja, _$A_{\infty}$ -Algebra from Supermanifolds_, Annales Henri Poincare 20 (2019) no.12, 4163-4195
  * [CCGN] R. Catenacci, C. A. Cremonini, P. A. Grassi and S. Noja, _On forms, cohomology and BV Laplacians in odd symplectic geometry_ , Lett Math Phys 111, 44 (2021)
  * [CCGN2] R. Catenacci, C. A. Cremonini, P. A. Grassi and S. Noja, _Cohomology of Lie Superalgebras: Forms, Integral Forms and Coset Superspaces_ , [arXiv:2012.05246 [hep-th]]
  * [CE] C. Chevalley, S. Eilenberg, _Cohomology Theory of Lie Groups and Lie Algebras_ , Trans. Am. Math. Soc. 63 (1948) 85
  * [CG] C. A. Cremonini and P. A. Grassi, _Pictures from Super Chern-Simons Theory_ , JHEP 03 (2020), 043
  * [CG2] C. A. Cremonini and P. A. Grassi, _Super Chern-Simons theory: Batalin-Vilkovisky formalism and $A_{\infty}$ algebras_, Phys. Rev. D 102 (2020) no.2, 025009
  * [CG3] C. A. Cremonini, P. A. Grassi et alii _In preparation_
  * [DL] M. J. Duff and J. X. Lu, _Type II p-branes: The Brane scan revisited_ , Nucl. Phys. B 390 (1993), 276-290
  * [Duf] M. J. Duff, _The conformal brane-scan: an update_ , [arXiv:2112.13784 [hep-th]].
  * [FiSaSc] D. Fiorenza, H. Sati and U. Schreiber, _Super Lie n-algebra extensions, higher WZW models, and super p-branes with tensor multiplet fields_ , Int. J. Geom. Meth. Mod. Phys. 12 (2014), 1550018
  * [FrSoSc] L. Frappat, P. Sorba, A. Sciarrino, _Dictionary on Lie Algebras and Superalgebras_ , Academic Press (2000)
  * [FMS] D. Friedan, E. J. Martinec and S. H. Shenker, _Conformal Invariance, Supersymmetry and String Theory_ , Nucl. Phys. B 271 (1986), 93-165
  * [Fu] D. Fuks, _Cohomology of Infinite-Dimensional Lie Algebras_ , (1986) Springer, New York, NY
  * [FH] W. Fulton and J. Harris, _Representation Theory: A First Course_ , New York: Springer (1991)
  * [GSW] M. B. Green, J. H. Schwarz and E. Witten, _SUPERSTRING THEORY. VOL. 1: INTRODUCTION_ , Cambridge Monographs on Mathematical Physics, 978-0-521-35752-4 (1988)
  * [GHV] W. Greub, S. Halperin, R. Vanstone, _Connections, Curvature, and Cohomology: Lie groups, principal bundles, and characteristic classes_ , Vol.3. Academic Press (1973)
  * [GLS] P. Grozman, D. A. Leites and I. Shchepochkina, _Lie superalgebras of string theories_ , Acta Math. Vietnamica 26 (2005), 27-63
  * [HS] G. Hochschild, J. P. Serre, _Cohomology of Lie Algebras_ , Ann. Math. 57 (1953) 591
  * [HLP] J. Hughes, J. Liu and J. Polchinski, _Supermembranes_ , Phys. Lett. B 180 (1986), 370-374
  * [Kac] V. G. Kac, _Lie Superalgebras_ , Adv. Math. 26 (1977), 8-96
  * [Kos] J. L. Koszul, _Homologie et Cohomologie des Algebres de Lie_ , Bull. Soc. Math. France, 78 (1950) pp. 65-127
  * [LLS] A. Lebedev, D. A. Leites and I. Shereshevskii, _Lie superalgebra structures in cohomology spaces of Lie algebras with coefficients in the adjoint representation_ , Am. Math. Soc. Transl. 213 (2005), 157-172
  * [LZ] G. I. Lehrer and R. B. Zhang, _The First Fundamental Theorem of Invariant Theory for the Orthosymplectic Supergroup_ , Commun. Math. Phys. 349 (2017) no.2, 661-702
  * [LZ2] G. I. Lehrer and R. B. Zhang, _The second fundamental theorem of invariant theory for the orthosymplectic supergroup_ , Nagoya Mathematical Journal, 242, (2021) 52-76
  * [Lei] D. A. Leites, _Representations of Lie Superalgebras_ , Theor. Math. Phys. 52 (1982), 764-766
  * [Man] Y. I. Manin, _Gauge Field Theory and Complex Geometry_ , Berlin, Germany: Springer (1988) 295 p. (Grundlehren der Mathematischen Wissenschaften, 289)
  * [NR] S. Noja, R. Re, _A note on super Koszul complex and the Berezinian_ , Annali di Matematica 201, 403–421 (2022)
  * [OP] O. V. Ogievetskii, I. B. Penkov, _Serre Duality for Projective Supermanifolds_ , Funct. Anal. its Appl. 18 68-70(1984)
  * [Pen] I. B. Penkov, _D-Modules on Supermanifolds_ , Invent. Math. 71, 501-512, (1983)
  * [ScZ] M. Scheunert and R. B. Zhang, _Cohomology of Lie superalgebras and of their generalizations_ , J. Math. Phys. 39 (1998), 5024-5061
  * [SuZ] Y. Su and R. B. Zhang, _Cohomology of Lie superalgebras $sl_{m|n}$ and $osp_{2|2n}$_, Proceedings of the London Mathematical Society, Volume 94, Issue 1, January 2007, Pages 91–136
  * [SuZ2] Y. Su and R. B. Zhang, _Mixed cohomology of Lie superalgebras_ , J. Algebra 549 (2020), 1-29
  * [Sul] D. Sullivan, _Infinitesimal computations in topology_ , Publications Mathématiques de l’IHÉS, Tome 47 (1977), pp. 269-331
  * [VV] E. P. Verlinde and H. L. Verlinde, _Multiloop Calculations in Covariant Superstring Theory_ , Phys. Lett. B 192 (1987), 95-102
  * [Wit] E. Witten, _Notes on Supermanifolds and Integration_ , Pure Appl. Math. Quart. 15 (2019) no.1, 3-56
  * [Wit2] E. Witten, _Perturbative superstring theory revisited_ , Pure Appl. Math. Quart. 15 (2019) no.1, 213-516