UNIFORM ESTIMATES FOR SOLUTIONS OF NONLINEAR FOCUSING DAMPED WAVE EQUATIONS

Thomas Perrin

Abstract. For a damped wave (or Klein-Gordon) equation on a bounded domain,
with a focusing power-like nonlinearity satisfying some growth conditions, we
prove that a global solution is bounded in the energy space, uniformly in
time. Our result applies in particular to the case of a cubic equation on a
bounded domain of dimension $3$.

## Introduction

#### Setting and main result.

Let $\Omega$ be a compact Riemannian manifold of dimension $d\geq 3$, with or
without boundary. Let $\beta\in\mathbb{R}$ be such that the Poincaré
inequality

$\int_{\Omega}\left(\left|\nabla
u\right|^{2}+\beta|u|^{2}\right)\mathrm{d}x\gtrsim\int_{\Omega}|u|^{2}\mathrm{d}x,\quad
u\in H_{0}^{1}(\Omega),$ (1)

is satisfied. For real-valued initial data $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, $\gamma\in
L^{\infty}(\Omega,\mathbb{R}_{+})$, and
$f\in\mathscr{C}^{1}(\mathbb{R},\mathbb{R})$, consider the equation

$\left\\{\begin{array}[]{rcccl}\square u+\gamma\partial_{t}u+\beta
u&=&f(u)&&\text{in }\mathbb{R}\times\Omega,\\\
(u(0),\partial_{t}u(0))&=&\left(u^{0},u^{1}\right)&&\text{in }\Omega,\\\
u&=&0&&\text{on }\mathbb{R}\times\partial\Omega.\end{array}\right.$

Set $F(s)=\int_{0}^{s}f(\tau)\mathrm{d}\tau$. We assume that $f(0)=0$, and
that there exist $1<p<+\infty$, $1<p_{0}<+\infty$, $C>0$, $\varepsilon>0$ such
that

$\left|f(s_{1})-f(s_{2})\right|\leq
C|s_{1}-s_{2}|\left(1+|s_{1}|^{p-1}+|s_{2}|^{p-1}\right),\quad
s_{1},s_{2}\in\mathbb{R},$ (1) $sf(s)\geq(2+\varepsilon)F(s),\quad
s\in\mathbb{R},$ (2)

and

$\left\|u\right\|_{L^{2}}^{p_{0}+1}\leq C\int_{\Omega}F(u)\mathrm{d}x,\quad
u\in H_{0}^{1}(\Omega).$ (3)

The fact that $\Omega$ is bounded will only be used to ensure the existence of
a nonlinearity satisfying (1), (2) and (3). We will always assume that
$\Omega$ and $f$ are such that the Cauchy problem (Setting and main result.)
can be solved locally (see Theorem 2 for a precise statement).

For $\left(u^{0},u^{1}\right)\in H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, write
$\left\|u^{0}\right\|_{H_{0}^{1}}^{2}=\int_{\Omega}\left(\left|\nabla
u^{0}\right|^{2}+\beta|u^{0}|^{2}\right)\mathrm{d}x$ and set

$E\left(u^{0},u^{1}\right)=\frac{1}{2}\left\|u^{0}\right\|_{H_{0}^{1}}^{2}+\frac{1}{2}\left\|u^{1}\right\|_{L^{2}}^{2}-\int_{\Omega}F(u^{0})\mathrm{d}x.$

For $\left(u^{0},u^{1}\right)\in H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, if
the solution $u$ of (Setting and main result.) exists on some interval
$[0,T]$, then one has the energy equality

$E\left(u(t_{1}),\partial_{t}u(t_{1})\right)=E\left(u(t_{0}),\partial_{t}u(t_{0})\right)-\int_{t_{0}}^{t_{1}}\int_{\Omega}\gamma\left|\partial_{t}u\right|^{2}\mathrm{d}t\mathrm{d}x,\quad
0\leq t_{0}\leq t_{1}\leq T.$ (4)

Our main result is the following.

######  Theorem 1.

Assume that $f$ satisfies (1), with $p<\frac{d+2}{d-2}$, (2), and (3).

  1. (i)

There exists $C=C(f,\gamma)>0$ such that for all $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, if the solution $u$ of (Setting and
main result.) exists on $\mathbb{R}_{+}$, then

$E\left(u^{0},u^{1}\right)\geq E\left(u(t),\partial_{t}u(t)\right)\geq-C,\quad
t\geq 0.$

  2. (ii)

There exist $c_{0}=c_{0}(f,\gamma)>0$, $c_{1}=c_{1}(f,\gamma)>0$, and
$c_{2}=c_{2}(f)>0$ such that for all $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, if the solution $u$ of (Setting and
main result.) exists on $\mathbb{R}_{+}$, then

$\left\|u(t)\right\|_{L^{2}}^{2}\leq\left\|u^{0}\right\|_{L^{2}}^{2}e^{-c_{2}t}+\left(c_{0}+c_{1}\left|E\left(u^{0},u^{1}\right)\right|\right)\left(1-e^{-c_{2}t}\right),\quad
t\geq 0.$ (5)

  3. (iii)

Assume, in addition, that $f$ satisfies (1) with $p\leq\frac{d}{d-2}$. Then,
there exists $c=c(f,\gamma)>0$ such that the function $\alpha:s\mapsto
c\exp(cs)$ satisfies the following property. For all
$\left(u^{0},u^{1}\right)\in H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, if the
solution $u$ of (Setting and main result.) exists on $\mathbb{R}_{+}$, then

$\left\|\left(u(t),\partial_{t}u(t)\right)\right\|_{H_{0}^{1}\times
L^{2}}\leq\alpha\left(\left\|u^{0}\right\|_{H_{0}^{1}}^{2}+\left\|u^{1}\right\|_{L^{2}}^{2}\right),\quad
t\geq 0.$ (6)

In addition, there exists $T\geq 0$, which depends on $f$, $\gamma$ and
$\left\|u^{0}\right\|_{H_{0}^{1}}^{2}+\left\|u^{1}\right\|_{L^{2}}^{2}$, such
that

$\left\|\left(u(t),\partial_{t}u(t)\right)\right\|_{H_{0}^{1}\times
L^{2}}\leq\alpha\left(\left|E\left(u^{0},u^{1}\right)\right|\right),\quad
t\geq T.$ (7)

All constants in this article may depend on $\Omega$. Note that when we
specify that a constant depends on $f$, it depends, in fact, only on the
constants appearing in (1), (2), and (3).

#### Examples.

A typical example of $f$ satisfying (1) for some $p>1$, (2), and (3) for some
$p_{0}>1$, is given by

$f(s)=\lambda_{1}s|s|^{\alpha_{1}-1}+\cdots+\lambda_{n}s|s|^{\alpha_{n}-1},\quad
s\in\mathbb{R},$

with $n\in\mathbb{N}^{\ast}$, $\lambda_{i}>0$ for all $i\in\llbracket
1,n\rrbracket$, and $p_{0}\leq\alpha_{1}\leq\cdots\leq\alpha_{n}\leq p$.
Assumption (3) is satisfied by such $f$ since $\Omega$ is bounded. As
previously mentioned, the fact that $\Omega$ is bounded will not be used
elsewhere in this article. Theorem 1-_(i)_ and Theorem 1-_(ii)_ can be applied
to such $f$ if $p<\frac{d+2}{d-2}$, and Theorem 1-_(iii)_ can also be applied
if $p\leq\frac{d}{d-2}$. For example, Theorem 1 applies to the focusing cubic
wave equation on a compact manifold of dimension $3$.

#### Connection with the existing literature.

The main objective of this article is to extend the results of [4] to the case
of a damped equation. In [4], Cazenave assumes that $f$ satisfies (1) with
$p=\frac{d}{d-2}$ if $d\geq 3$. Note that this assumption is not satisfied by
all subcritical nonlinearity. The usual critical exponent $\frac{d+2}{d-2}$ is
the one that gives $F(u)\in L^{1}(\Omega)$ by the Sobolev embedding
$H^{1}\hookrightarrow L^{\frac{2d}{d-2}}$, implying that the energy is finite,
whereas (1) yields $f(u)\in L^{2}(\Omega)$. However, some of the results of
[4] are valid for $f$ satisfying (2) with $p<\frac{d+2}{d-2}$ if $d\geq 3$ :
more precisely, a careful reading of [4] shows that if $\gamma=0$, then
Theorem 1 holds true without assumption (3), and with $\Omega$ bounded or
unbounded. Unfortunately, we have not managed to find a proof that is valid
for any damping $\gamma$, without the additional assumption (3). The author of
[4] also considers the cases $d=1,2$, with an additional assumption on $f$
when $d=2$. For simplicity, we have not included the cases $d=1,2$, although
the methods of this article and of [4] can be adapted to those cases.

In [5], Feireisl considers the case $\Omega=\mathbb{R}^{d}$, with a damping
equal to a positive constant, and a nonlinearity satisfying the restrictive
condition $1<p<1+\min\left(\frac{2}{d-2},\frac{4}{d}\right)$. Lemma 4.2 of [5]
states that a solution with a bounded energy is global (for positive times),
and is bounded in $H^{1}(\mathbb{R}^{d})\times L^{2}(\mathbb{R}^{d})$.
However, we have not found in [5] a proof of the fact that a global solution
has a bounded energy.

Some related work has been done by Burq, Raugel and Schlag. More precisely, in
[3], the authors prove that a global solution is bounded in the energy space,
in the case $\Omega=\mathbb{R}^{d}$, with a damping equal to a positive
constant, for a general focusing subcritical nonlinearity, but with a radial
assumption. The proof of the boundedness of [3] uses dynamical systems
arguments. In [2], they consider the case of a damping which is constant in
space, and tends to zero as $t$ tends to infinity, still with
$\Omega=\mathbb{R}^{d}$ and a radial assumption. They adapt the proof of the
$L^{2}$-bound of Cazenave, using the fact that their damping tends to zero as
$t$ tends to infinity.

#### Outline of the article.

In Section 1, we recall the local Cauchy theory for (Setting and main
result.), and the expression of the derivative of the square of the
$L^{2}$-norm of a solution. In Section 2, 3 and 4, we prove respectively the
energy bound (Theorem 1-_(i)_), the $L^{2}$-bound (Theorem 1-_(ii)_), and the
$H^{1}$-bound (Theorem 1-_(iii)_). An appendix contains proofs of some
elementary lemmas.

#### Acknowledgements.

I warmly thank Thomas Duyckaerts and Jérôme Le Rousseau for their constant
support and guidance.

##  1 Preliminaries

We recall the local Cauchy theory for (Setting and main result.).

######  Theorem 2.

Consider $f$ satisfying (1) for some $p<\frac{d+2}{d-2}$. For any real-valued
initial data $\left(u^{0},u^{1}\right)\in H_{0}^{1}(\Omega)\times
L^{2}(\Omega)$, there exist a maximal time of existence $T\in(0,+\infty]$ and
a unique solution $u$ of (Setting and main result.) in
$\mathscr{C}^{0}([0,T),H_{0}^{1}(\Omega))\cap\mathscr{C}^{1}([0,T),L^{2}(\Omega))$.
If $T<+\infty$, then

$\left\|\left(u(t),\partial_{t}u(t)\right)\right\|_{H_{0}^{1}\times
L^{2}}\underset{t\rightarrow T^{-}}{\longrightarrow}+\infty.$

The proof of Theorem 2 relies on Strichartz estimates and a fixed point
argument (see for example [6]). If $\Omega$ is a manifold without boundary,
then Theorem 2 follows from the Strichartz estimates proved in [8]. In the
case of a manifold with boundary, we _assume_ that $f$ is such that Theorem 2
holds true : we refer to [1] for Strichartz estimates for a large range of
exponents. Ivanovici’s counterexamples in [7] show that Strichartz estimates
are not true for the full range of exponents in the case of a manifold with
boundary.

######  Lemma 3.

Consider $u$ a solution of (Setting and main result.) on some interval
$[0,T]$, and set $M(t)=\left\|u(t)\right\|_{L^{2}}^{2}$, for $t\in[0,T]$. Then
$M\in\mathscr{C}^{2}([0,T],\mathbb{R})$, and for $t\in[0,T]$, one has
$M^{\prime}(t)=2\int_{\Omega}u(t)\partial_{t}u(t)\mathrm{d}x$, and

$M^{\prime\prime}(t)=2\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}-2\left\|u(t)\right\|_{H_{0}^{1}}^{2}+2\int_{\Omega}u(t)f(u(t))\mathrm{d}x-2\int_{\Omega}\gamma
u(t)\partial_{t}u(t)\mathrm{d}x.$ (8)

The proof is standard and we omit it. See for example [4] or [9].

##  2 The energy bound

Here, we prove Theorem 1-_(i)_. Consider $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$ such that the solution $u$ of (Setting
and main result.) exists on $\mathbb{R}_{+}$. By (3), there exists
$C_{0}=C_{0}(f)>0$ such that

$\left|2\int_{\Omega}\gamma u(t)\partial_{t}u(t)\mathrm{d}x\right|\leq
C_{0}\left\|\gamma\right\|^{2}_{L^{\infty}}\left(\int_{\Omega}F(u(t))\mathrm{d}x\right)^{\frac{2}{p_{0}+1}}+\frac{\varepsilon}{4}\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2},\quad
t\geq 0,$ (9)

where $\varepsilon$ is given by (2). Using (2) and (9) in (8), one finds

$\displaystyle M^{\prime\prime}(t)\geq$
$\displaystyle\left(2-\frac{\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}-2\left\|u(t)\right\|_{H_{0}^{1}}^{2}$
$\displaystyle+2(2+\varepsilon)\int_{\Omega}F(u(t))\mathrm{d}x-C_{0}\left\|\gamma\right\|^{2}_{L^{\infty}}\left(\int_{\Omega}F(u(t))\mathrm{d}x\right)^{\frac{2}{p_{0}+1}},\quad
t\geq 0.$ (10)

Write
$I_{0}=\left(\frac{C_{0}\left\|\gamma\right\|^{2}_{L^{\infty}}}{\varepsilon}\right)^{\frac{p_{0}+1}{p_{0}-1}}\geq
0$, so that

$I\geq I_{0}\quad\Longrightarrow\quad
2(2+\varepsilon)I-C_{0}\left\|\gamma\right\|^{2}_{L^{\infty}}I^{\frac{2}{p_{0}+1}}\geq(4+\varepsilon)I.$

Assume by contradiction that there exists $T_{0}\geq 0$ such that
$E\left(u(T_{0}),\partial_{t}u(T_{0})\right)<-I_{0}$. Then for $t\geq T_{0}$,
one has
$\int_{\Omega}F(u(t))\mathrm{d}x\geq-E\left(u(t),\partial_{t}u(t)\right)>I_{0}$,
yielding

$M^{\prime\prime}(t)\geq\left(2-\frac{\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}-2\left\|u(t)\right\|_{H_{0}^{1}}^{2}+(4+\varepsilon)\int_{\Omega}F(u(t))\mathrm{d}x,\quad
t\geq T_{0}.$

By definition of the energy, this gives

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\frac{\varepsilon}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}-(4+\varepsilon)E\left(u(t),\partial_{t}u(t)\right),\quad
t\geq T_{0}.$ (11)

One has $E\left(u(t),\partial_{t}u(t)\right)\leq 0$ for $t\geq T_{0}$,
implying

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2},\quad
t\geq T_{0}.$

In particular, one obtains

$\left(M^{\prime}(t)\right)^{2}\leq\left(1+\frac{\varepsilon}{16}\right)^{-1}M(t)M^{\prime\prime}(t),\quad
t\geq T_{0}.$

We use the following lemma.

######  Lemma 4.

Consider $T\geq 0$, $M\in\mathscr{C}^{2}([T,+\infty),\mathbb{R}_{+})$ and
$0<\delta<1$ such that $\left(M^{\prime}(t)\right)^{2}\leq\delta
M(t)M^{\prime\prime}(t)$ for $t\geq T$. Then $M$ is non-increasing. In
particular, one has $M(t)\leq M(T)$ for $t\geq T$.

By (11), one has
$M^{\prime\prime}(t)\geq-(4+\varepsilon)E\left(u(T_{0}),\partial_{t}u(T_{0})\right)>0$
for $t\geq T_{0}$, a contradiction with Lemma 4. Hence, one has
$E\left(u(t),\partial_{t}u(t)\right)\geq-I_{0}$ for all $t\geq 0$.

######  Remark 5.

If $\gamma=0$, then $I_{0}=0$ : we recover the fact that a global solution of
the undamped equation has a non-negative energy.

##  3 The $L^{2}(\Omega)$-estimate

Here, we prove Theorem 1-_(ii)_. Consider $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$ such that the solution $u$ of (Setting
and main result.) exists on $\mathbb{R}_{+}$. We split the proof into 2 steps.

#### Step 1 : an estimate on $M^{\prime\prime}$.

Set $r_{0}=\frac{2}{p_{0}+1}\in(0,1)$. Using (10) and the definition of the
energy, one finds

$\displaystyle M^{\prime\prime}(t)\geq$
$\displaystyle\left(4+\frac{3\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\varepsilon\left\|u(t)\right\|_{H_{0}^{1}}^{2}-2(2+\varepsilon)E\left(u(t),\partial_{t}u(t)\right)$
$\displaystyle-
C_{0}\left\|\gamma\right\|^{2}_{L^{\infty}}\left(\frac{1}{2}\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\frac{1}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}-E\left(u(t),\partial_{t}u(t)\right)\right)^{r_{0}},\quad
t\geq 0,$

with $C_{0}=C_{0}(f)$. There exists $C_{1}=C_{1}(\gamma,f)>0$ such that

$\displaystyle M^{\prime\prime}(t)\geq$
$\displaystyle\left(4+\frac{3\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\varepsilon\left\|u(t)\right\|_{H_{0}^{1}}^{2}-2(2+\varepsilon)E\left(u(t),\partial_{t}u(t)\right)$
$\displaystyle-\frac{C_{1}}{2}\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2r_{0}}-C_{1}\left|\frac{1}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}-E\left(u(t),\partial_{t}u(t)\right)\right|^{r_{0}},\quad
t\geq 0.$

By Theorem 1-_(i)_ , one has $\sup_{t\geq
0}\left|E\left(u(t),\partial_{t}u(t)\right)\right|<+\infty$. Hence, there
exists

$C_{2}=C_{2}\left(\gamma,f,\sup_{t\geq
0}\left|E\left(u(t),\partial_{t}u(t)\right)\right|\right)>0$

such that

$\varepsilon
X^{2}-2(2+\varepsilon)E\left(u(t),\partial_{t}u(t)\right)-C_{1}\left|\frac{1}{2}X^{2}-E\left(u(t),\partial_{t}u(t)\right)\right|^{r_{0}}\geq\frac{\varepsilon}{2}X^{2}-C_{2},\quad
t\geq 0,\quad X\geq 0,$

yielding

$M^{\prime\prime}(t)\geq\left(4+\frac{3\varepsilon}{4}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}-\frac{C_{1}}{2}\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2r_{0}}+\frac{\varepsilon}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}-C_{2},\quad
t\geq 0.$

Note that $C_{2}$ can be chosen of the form

$C_{2}=\tilde{C}_{0}+\tilde{C}_{1}\sup_{t\geq
0}\left|E\left(u(t),\partial_{t}u(t)\right)\right|,$

with $\tilde{C}_{0}=\tilde{C}_{0}\left(\gamma,f\right)$ and
$\tilde{C}_{1}=\tilde{C}_{1}\left(\gamma,f\right)$. Finally, there exists
$C_{3}=C_{3}\left(\gamma,f\right)>0$ such that

$\left(4+\frac{3\varepsilon}{4}\right)X^{2}-\frac{C_{1}}{2}X^{2r_{0}}\geq\left(4+\frac{\varepsilon}{2}\right)X^{2}-C_{3},\quad
X\geq 0,$

and this gives

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{2}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\frac{\varepsilon}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}-C_{2}-C_{3},\quad
t\geq 0.$

Summarizing up, and using the Poincaré inequality (1), one obtains

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{2}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+C_{4}\left\|u(t)\right\|_{L^{2}}^{2}-C_{5},\quad
t\geq 0,$ (12)

for some $C_{4}=C_{4}(f)>0$ and $C_{5}>0$ of the form

$C_{5}=\tilde{C}_{2}+\tilde{C}_{3}\sup_{t\geq
0}\left|E\left(u(t),\partial_{t}u(t)\right)\right|,$ (13)

with $\tilde{C}_{2}=\tilde{C}_{2}(f,\gamma)>0$ and
$\tilde{C}_{3}=\tilde{C}_{3}(f,\gamma)>0$.

#### Step 2 : proof of (5).

The proof is based on the following elementary lemma.

######  Lemma 6.

Consider $M_{0}\in\mathscr{C}^{2}(\mathbb{R}_{+},\mathbb{R})$ such that there
exists $C>0$ such that $M_{0}^{\prime\prime}(t)\geq CM_{0}(t)$, for all $t\geq
0$. Then either $M_{0}(t)\rightarrow+\infty$ as $t\rightarrow+\infty$, or
$M_{0}(t)\leq M_{0}(0)e^{-\sqrt{C}t}$, for all $t\geq 0$.

Set $M_{0}(t)=C_{4}M(t)-C_{5}$, $t\geq 0$. By (12), one has
$M_{0}^{\prime\prime}(t)\geq C_{4}M_{0}(t)$, for all $t\geq 0$. Assume that
$M_{0}(t)\rightarrow+\infty$ as $t\rightarrow+\infty$. Then, there exists
$T\geq 0$ such that $M_{0}(t)\geq 0$, for $t\geq T$, implying

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{2}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2},\quad
t\geq T,$

by (12). This gives

$\left(M^{\prime}(t)\right)^{2}\leq\left(1+\frac{\varepsilon}{8}\right)^{-1}M(t)M^{\prime\prime}(t),\quad
t\geq T.$

Applying Lemma 4, one finds that $M$ is bounded, a contradiction with
$M_{0}(t)\rightarrow+\infty$. Hence, by Lemma 6, one obtains $M_{0}(t)\leq
M_{0}(0)e^{-\sqrt{C_{4}}t}$ for $t\geq 0$, that is,

$M(t)\leq\frac{C_{5}}{C_{4}}+\left(M(0)-\frac{C_{5}}{C_{4}}\right)e^{-\sqrt{C_{4}}t},\quad
t\geq 0.$ (14)

Note that using Theorem 1-_(i)_ and the fact that the energy is non-
increasing, one has

$\sup_{t\geq
0}\left|E\left(u(t),\partial_{t}u(t)\right)\right|\leq\max\left(\left|E\left(u^{0},u^{1}\right)\right|,C_{6}\right)\leq\left|E\left(u^{0},u^{1}\right)\right|+C_{6},$
(15)

for some $C_{6}=C_{6}(f,\gamma)>0$. By (13), this gives

$\frac{C_{5}}{C_{4}}\leq c_{0}+c_{1}\left|E\left(u^{0},u^{1}\right)\right|,$

for some $c_{0}=c_{0}(f,\gamma)>0$ and $c_{1}=c_{1}(f,\gamma)>0$. Using this
in (14), one obtains (5).

##  4 The $H^{1}(\Omega)$-estimates

Here, we prove Theorem 1-_(iii)_. We first need an estimate on $M^{\prime}$.
Note that it does not require the restrictive exponent $\frac{d}{d-2}$.

######  Lemma 7.

Assume that $d\geq 3$, and that $f$ satisfies (1), with $p<\frac{d+2}{d-2}$,
(2), and (3). Then, there exist $c_{0}=c_{0}(f,\gamma)>0$,
$c_{1}=c_{1}(f,\gamma)>0$, and $c_{2}=c_{2}(f)>0$ such that for all
$\left(u^{0},u^{1}\right)\in H_{0}^{1}(\Omega)\times L^{2}(\Omega)$, if the
solution $u$ of (Setting and main result.) exists on $\mathbb{R}_{+}$, then

$M^{\prime}(0)e^{-c_{2}t}-\left(c_{0}+c_{1}\left|E\left(u^{0},u^{1}\right)\right|\right)\left(1-e^{-c_{2}t}\right)\leq
M^{\prime}(t)\leq c_{0}+c_{1}\left|E\left(u^{0},u^{1}\right)\right|,\quad
t\geq 0.$ (16)

###### Proof.

Using the proof of Theorem 1-_(ii)_ (see (12), (13) and (15)), one has

$M^{\prime\prime}(t)\geq\left(4+\frac{\varepsilon}{2}\right)\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+C_{0}\left\|u(t)\right\|_{L^{2}}^{2}-C_{1},\quad
t\geq 0,$ (17)

for some $C_{0}=C_{0}(f)>0$ and $C_{1}>0$ of the form

$\tilde{C}_{0}+\tilde{C}_{1}\left|E\left(u^{0},u^{1}\right)\right|,$ (18)

with $\tilde{C}_{0}=\tilde{C}_{0}(f,\gamma)>0$ and
$\tilde{C}_{1}=\tilde{C}_{1}(f,\gamma)>0$. For $\eta>0$, using (17), one finds

$\displaystyle\left|M^{\prime}(t)\right|$
$\displaystyle\leq\eta\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}+\frac{1}{\eta}\left\|u(t)\right\|_{L^{2}}^{2}$
$\displaystyle\leq\eta\left(4+\frac{\varepsilon}{2}\right)^{-1}M^{\prime\prime}(t)+\left(\frac{1}{\eta}-C_{0}\eta\left(4+\frac{\varepsilon}{2}\right)^{-1}\right)\left\|u(t)\right\|_{L^{2}}^{2}+C_{1}\eta\left(4+\frac{\varepsilon}{2}\right)^{-1},\quad
t\geq 0.$

Choosing $\eta$ such that
$\frac{1}{\eta}=C_{0}\eta\left(4+\frac{\varepsilon}{2}\right)^{-1}$, one
obtains

$\left|M^{\prime}(t)\right|\leq C_{2}M^{\prime\prime}(t)+C_{3},\quad t\geq 0,$
(19)

with $C_{2}=C_{2}(f)>0$ and $C_{3}>0$ of the form (18). Integrating (19), one
finds

$\left(M^{\prime}(t_{0})-C_{3}\right)e^{-\frac{t_{0}}{C_{2}}}\leq\left(M^{\prime}(t_{1})-C_{3}\right)e^{-\frac{t_{1}}{C_{2}}},\quad
0\leq t_{0}\leq t_{1},$ (20)

and

$M^{\prime}(0)+C_{3}\leq\left(M^{\prime}(t_{1})+C_{3}\right)e^{\frac{t_{1}}{C_{2}}},\quad
t_{1}\geq 0.$ (21)

If there exists $t_{0}\geq 0$ such that $M^{\prime}(t_{0})-C_{3}>0$, then (20)
implies that $M^{\prime}(t)\rightarrow+\infty$ as $t\rightarrow+\infty$, a
contradiction with Theorem 1-_(ii)_. Hence, one has $M^{\prime}(t)\leq C_{3}$,
for all $t\geq 0$. This gives the upper bound of (16). The lower bound of (16)
is given by (21). ∎

Now, we prove Theorem 1-_(iii)_. Consider $\left(u^{0},u^{1}\right)\in
H_{0}^{1}(\Omega)\times L^{2}(\Omega)$ such that the solution $u$ of (Setting
and main result.) exists on $\mathbb{R}_{+}$. Set

$E_{\mathscr{L}}(t)=\frac{1}{2}\left\|u(t)\right\|_{H_{0}^{1}}^{2}+\frac{1}{2}\left\|\partial_{t}u(t)\right\|_{L^{2}}^{2}=E\left(u(t),\partial_{t}u(t)\right)+\int_{\Omega}F(u(t))\mathrm{d}x,\quad
t\geq 0.$

By (17) and (18), one has

$E_{\mathscr{L}}(t)\leq C_{0}M^{\prime\prime}(t)+C_{1},\quad t\geq 0,$ (22)

for some $C_{0}=C_{0}(f)>0$, and $C_{1}>0$ of the form
$\tilde{C}_{0}+\tilde{C}_{1}\left|E\left(u^{0},u^{1}\right)\right|$, with
$\tilde{C}_{0}=\tilde{C}_{0}(f,\gamma)>0$ and
$\tilde{C}_{1}=\tilde{C}_{1}(f,\gamma)>0$.

######  Remark 8.

In particular, without the restrictive assumption of Theorem 1-_(iii)_ , Lemma
7 and (22) imply

$\int_{0}^{T}E_{\mathscr{L}}(t)\mathrm{d}t\lesssim 1+T,\quad T\geq 0,$

with a constant depending on $f$, $\gamma$, and $\left(u^{0},u^{1}\right)$.
This is an average estimate of $E_{\mathscr{L}}$. It won’t be used later on :
to derive a non-average estimate, it seems necessary to estimate the
derivative of $E_{\mathscr{L}}$, which is undefined without the additional
assumption of Theorem 1-_(iii)_.

For the rest of the proof, we write $C$ and $C^{\prime}$ for some positive
constants, that may change from line to line, and which depend on $\Omega$,
$f$ and $\gamma$. We assume that $f$ satisfies (1) with $p\leq\frac{d}{d-2}$.
Together with the Sobolev embedding $H^{1}(\Omega)\hookrightarrow
L^{\frac{2d}{d-2}}$, it gives $f(u(t))\in L^{2}(\Omega)$, implying

$\int_{\Omega}\left|\partial_{t}u(t)f(u(t))\right|\mathrm{d}x\leq
C\left\|\partial_{t}u(t)\right\|_{L^{2}}\left(\left\|u(t)\right\|_{L^{2}}+\left\|u(t)\right\|_{H_{0}^{1}}^{p}\right),\quad
t\geq 0.$ (23)

In particular, the derivative of $E_{\mathscr{L}}$ is well-defined, and using
also (4), one has

$E_{\mathscr{L}}^{\prime}(t)=-\int_{\Omega}\gamma\left|\partial_{t}u\right|^{2}\mathrm{d}x+\int_{\Omega}\partial_{t}u(t)f(u(t))\mathrm{d}x\leq\int_{\Omega}\partial_{t}u(t)f(u(t))\mathrm{d}x,\quad
t\geq 0.$

Together with (1) and (23), one finds

$E_{\mathscr{L}}^{\prime}(t)\leq
C\left\|\partial_{t}u(t)\right\|_{L^{2}}\left\|u(t)\right\|_{H_{0}^{1}}\left(1+\left\|u(t)\right\|_{H_{0}^{1}}^{p-1}\right),\quad
t\geq 0,$

Using a second time the fact that $p\leq\frac{d}{d-2}$, and $d\geq 3$, one has
$p-1\leq 2$, implying

$E_{\mathscr{L}}^{\prime}(t)\leq
CE_{\mathscr{L}}(t)\left(1+E_{\mathscr{L}}(t)\right),\quad t\geq 0.$

Integrating, this gives

$E_{\mathscr{L}}(t_{1})\leq
E_{\mathscr{L}}(t_{0})\exp\left(C\int_{t_{0}}^{t_{1}}\left(1+E_{\mathscr{L}}(s)\right)\mathrm{d}s\right),\quad
t_{1}\geq t_{0}\geq 0.$ (24)

#### Estimate for $t\leq 1$.

Choosing $t_{0}=0$ in (24) and using (22), one finds

$E_{\mathscr{L}}(t_{1})\leq
E_{\mathscr{L}}(0)\exp\left(\int_{0}^{1}\left(C+C^{\prime}M^{\prime\prime}(s)\right)\mathrm{d}s\right),\quad
1\geq t_{1}\geq 0.$

Note that $\left|E\left(u^{0},u^{1}\right)\right|\lesssim E_{\mathscr{L}}(0)$,
with a constant depending on $f$, and that $\left|M^{\prime}(0)\right|\lesssim
E_{\mathscr{L}}(0)$. Hence, using Lemma 7, one obtains

$E_{\mathscr{L}}(t_{1})\leq
E_{\mathscr{L}}(0)\exp\left(C+C^{\prime}E_{\mathscr{L}}(0)\right),\quad 1\geq
t_{1}\geq 0,$

We choose the function $\alpha$ of Theorem 1-_(iii)_ so that $\alpha(s)\geq
2s\exp\left(C+2C^{\prime}s\right)$, for all $s\in[0,1]$.

#### Estimate for $t\geq 1$.

Integrating (24), and using (22), one finds

$E_{\mathscr{L}}(t_{1})\leq\left(\int_{t_{1}-1}^{t_{1}}\left(C+C^{\prime}M^{\prime\prime}(t)\right)\mathrm{d}t\right)\exp\left(\int_{t_{1}-1}^{t_{1}}\left(C+C^{\prime}M^{\prime\prime}(s)\right)\mathrm{d}s\right),\quad
t_{1}\geq 1.$ (25)

Lemma 7 gives

$\int_{t_{1}-1}^{t_{1}}M^{\prime\prime}(t)\mathrm{d}t\leq
2\left(c_{0}+c_{1}\left|E\left(u^{0},u^{1}\right)\right|\right)+\left|M^{\prime}(0)\right|e^{c_{2}(1-t_{1})},\quad
t_{1}\geq 1,$ (26)

where $c_{0}$, $c_{1}$ and $c_{2}$ are the constants of Lemma 7.

Using (25) and (26), we can complete the proof of Theorem 1-_(iii)_. We first
prove (6). By (26), one has

$\int_{t_{1}-1}^{t_{1}}M^{\prime\prime}(t)\mathrm{d}t\leq
C+C^{\prime}E_{\mathscr{L}}(0),\quad t_{1}\geq 1.$

Hence, (25) gives

$E_{\mathscr{L}}(t_{1})\leq\left(C+C^{\prime}E_{\mathscr{L}}(0)\right)\exp\left(C+C^{\prime}E_{\mathscr{L}}(0)\right),\quad
t_{1}\geq 1,$

We choose the function $\alpha$ of Theorem 1-_(iii)_ so that
$\alpha(s)\geq\left(C+2C^{\prime}s\right)\exp\left(C+2C^{\prime}s\right)$, for
all $s\geq 1$. This completes the proof of (6).

Finally, we prove (7). There exists $T>0$, which depends on $f$, $\gamma$, and
$E_{\mathscr{L}}(0)$, such that

$\left|M^{\prime}(0)\right|e^{c_{2}(1-t_{1})}\leq c_{0},\quad t_{1}\geq T.$

Hence, (25) and (26) imply

$E_{\mathscr{L}}(t_{1})\leq\left(C+C^{\prime}\left|E\left(u^{0},u^{1}\right)\right|\right)\exp\left(C+C^{\prime}\left|E\left(u^{0},u^{1}\right)\right|\right),\quad
t_{1}\geq T.$

Up to increasing the constant $c$ defining the function $\alpha$, this
completes the proof of (7).

##  Appendix A Proof of the elementary lemmas

###  A.1 Proof of Lemma 4

Assume by contradiction that there exists $T_{0}\geq T$ such that
$M^{\prime}(T_{0})>0$. In particular, one has $M(T_{0})>0$. We claim that
$M(t)>0$ for all $t\geq T_{0}$. Indeed, assume by contradiction that there
exists $T_{1}>T_{0}$ such that $M(T_{1})=0$ and $M(t)>0$ for
$t\in[T_{0},T_{1})$. Then, one has $M^{\prime\prime}(t)\geq 0$ for
$t\in[T_{0},T_{1})$, implying that $M^{\prime}(t)\geq M^{\prime}(T_{0})>0$ for
$t\in[T_{0},T_{1})$. In particular, this gives $0=M(T_{1})\geq M(T_{0})>0$, a
contradiction.

As $M(t)>0$ for all $t\geq T_{0}$, one has $M^{\prime\prime}(t)>0$ for $t\geq
T_{0}$, implying $M^{\prime}(t)\geq M^{\prime}(T_{0})$, for $t\geq T_{0}$. In
particular, $M(t)\rightarrow+\infty$ as $t\rightarrow+\infty$. Integrating two
times the inequality

$\frac{M^{\prime}(t)}{M(t)}\leq\delta\frac{M^{\prime\prime}(t)}{M^{\prime}(t)},\quad
t\geq T_{0},$

one obtains

$\left(t-T_{0}\right)\frac{M^{\prime}(T_{0})}{M(T_{0})^{\frac{1}{\delta}}}\leq\left(\frac{1}{\delta}-1\right)^{-1}\left(\frac{1}{M(T_{0})^{\frac{1}{\delta}-1}}-\frac{1}{M(t)^{\frac{1}{\delta}-1}}\right),\quad
t\geq T_{0}.$

For $t$ sufficiently large, this gives a contradiction.

###  A.2 Proof of Lemma 6

Set $m(t)=M_{0}(t)e^{\sqrt{C}t}$, $t\geq 0$. One has $m^{\prime\prime}\geq
2\sqrt{C}m^{\prime}$, implying

$m(t_{1})\geq
m(t_{0})+\frac{m^{\prime}(t_{0})}{2\sqrt{C}}\left(e^{2\sqrt{C}(t_{1}-t_{0})}-1\right),\quad
t_{1}\geq t_{0}\geq 0.$

This gives

$M_{0}(t_{1})\geq
m(t_{0})e^{-\sqrt{C}t_{1}}+\frac{m^{\prime}(t_{0})}{2\sqrt{C}}\left(e^{\sqrt{C}(t_{1}-2t_{0})}-e^{-\sqrt{C}t_{1}}\right),\quad
t_{1}\geq t_{0}\geq 0.$

If there exists $t_{0}\geq 0$ such that $m^{\prime}(t_{0})>0$, then
$M_{0}(t)\rightarrow+\infty$. Conversely, if $m^{\prime}(t)\leq 0$ for all
$t\geq 0$, then one has $m(t)\leq m(0)$ for $t\geq 0$, that is, $M_{0}(t)\leq
M_{0}(0)e^{-\sqrt{C}t}$, for $t\geq 0$.

## References

  * [1] Matthew D. Blair, Hart F. Smith and Christopher D. Sogge  “Strichartz estimates for the wave equation on manifolds with boundary”  In _Ann. Inst. Henri Poincaré, Anal. Non Linéaire_ 26.5, 2009, pp. 1817–1829  DOI: 10.1016/j.anihpc.2008.12.004
  * [2] Nicolas Burq, Genevieve Raugel and Wilhelm Schlag  “Long time dynamics for weakly damped nonlinear Klein-Gordon equations”  arXiv, 2018  URL: https://arxiv.org/abs/1801.06735
  * [3] Nicolas Burq, Geneviève Raugel and Wilhelm Schlag  “Long time dynamics for damped Klein-Gordon equations”  In _Annales scientifiques de l’École normale supérieure_ 50.6 Societe Mathematique de France, 2017, pp. 1447–1498  DOI: 10.24033/asens.2349
  * [4] Thierry Cazenave  “Uniform estimates for solutions of nonlinear Klein-Gordon equations”  In _Journal of Functional Analysis_ 60, 1985, pp. 36–55  URL: https://api.semanticscholar.org/CorpusID:121770305
  * [5] Eduard Feireisl  “Finite Energy Travelling Waves for Nonlinear Damped Wave Equations”  In _Quarterly of Applied Mathematics_ 56.1 Brown University, 1998, pp. 55–70  URL: http://www.jstor.org/stable/43638210
  * [6] J. Ginibre and G. Velo  “The global Cauchy problem for the nonlinear Klein-Gordon equation. II”  In _Ann. Inst. Henri Poincaré, Anal. Non Linéaire_ 6.1, 1989, pp. 15–35  DOI: 10.1016/S0294-1449(16)30329-8
  * [7] Oana Ivanovici  “Counterexamples to the Strichartz inequalities for the wave equation in general domains with boundary”  In _Journal of the European Mathematical Society_ 014.5 European Mathematical Society Publishing House, 2012, pp. 1357–1388
  * [8] Markus Keel and Terence Tao  “Endpoint Strichartz Estimates”  In _American Journal of Mathematics_ 120.5 Johns Hopkins University Press, 1998, pp. 955–980
  * [9] L.. Payne and D.. Sattinger  “Saddle points and instability of nonlinear hyperbolic equations”  In _Israel Journal of Mathematics_ 22, 1975, pp. 273–303  DOI: 10.1007/BF02761595

Perrin Thomas<EMAIL_ADDRESS>
Laboratoire Analyse Géométrie et Application, Institut Galilée - UMR 7539,
CNRS/Université Sorbonne Paris Nord, 99 avenue J.B. Clément, 93430
Villetaneuse, France


# Charmonium production in high multiplicity pp collisions and the structure
of the proton

R. Terra1 and F. S. Navarra1 1Instituto de Física, Universidade de São Paulo,
Rua do Matão, 1371, CEP 05508-090, São Paulo, SP, Brazil

###### Abstract

In this work we study charmonium production in high multiplicity proton-proton
collisions. We investigate the role of the spatial distribution of partons in
the protons and assume that the proton has a Y shape. In this configuration
quarks are more at the surface and gluons in the inner part of the proton.
Going from peripheral to more central and then to ultra-central proton-proton
collisions, we go from quark-quark collisions to gluon-gluon collisions. Since
gluons are much more abundant, the cross sections grow. In the case of charm
production this growth is enhanced by the fact that, $\sigma(g+g\to
c+\bar{c})>>\sigma(q+\bar{q}\to c+\bar{c})$. These effects can explain the
growth seen in the data.

##  I Introduction

Surprising new features of proton proton collisions were revealed when the LHC
collaborations became able to trigger on very high multiplicity events alice-
psi7 ; alice-cc ; star18 ; alice19 ; alice20 ; alice-psi13 ; alice22 . In
these events, the data presented evidence of collective behavior, which could
be interpreted in terms of a hydrodynamical expansion of the system. However,
to the best of our knowledge, so far no particle production model,
hydrodynamic or non-hydrodynamic, matches all the features of the high-
multiplicity pp data.

High multiplicity events come from ultracentral collisions, with very small
impact parameter. In this regime we expect to observe several effects which
can be responsible for the anomalous features of particle production, such as
double parton scattering and parton saturation, caused by the larger
saturation scale at lower impact parameters. Apart from these changes in the
dynamics of the collisions, it is also possible that geometric effects
associated with the spatial distribution of strongly interacting matter play a
significant role.

Successfull models of the static proton are based on lattice QCD simulations,
which show that quarks are bound by gluonic strings. This leads to the “Y”
picture of the proton, where quarks are at the extremities tied by an Y-like
gluon string, called gluon junction suga15  or baryon junction. This structure
leads to a spatial configuration where the gluons are mostly in the center and
the quark closer to the proton surface.

At first sight the static geometric configuration of the proton could be
blurred or even completely washed out by a boost with the consequent parton
evolution and branching. However, there are indications that the geometric
organization of matter persists even at LHC energies. In manti23  it was shown
that exclusive vector meson production is sensitive to the geometric
deformation of the target. The authors conclude that the fluctuations in the
nuclear geometry originating from the deformed structure are not washed out by
the JIMWLK evolution and hence the deformations previously inferred from low-
energy experiments will be visible in high-energy collisions. This suggests
the the Y shape of the proton could survive the quantum evolution and manifest
itself in high energy proton proton collisions.

The existence of the Y-shape configuration of the proton (along with its gluon
junction) might shift the leading baryon distribution to smaller rapidities,
providing an additional mechanism of baryon stopping. In these events the
baryon number would be carried by the baryon junction, whereas in the
traditional approach it is carried by the valence quarks (for a discussion see
marb89 ). This interesting idea was proposed in khar96 , implemented in the
Monte-Carlo event generator HIJING/B top04  and also in analytical models such
as in bopp06 , being successfull in explaining the data on forward baryon
production. The search for new effects of the baryon junction is in progress
bran22 .

Having in mind what was said above, we would expect to see more manifestations
of the proton Y shape. This was explored in deb20 . One of the surprising
aspects of high multiplicity proton proton collisions is the observation of
the ridge effect and of elliptic flow. In Ref. deb20  a simple model based on
the proton spatial configuration was developed to explain these phenomena. The
authors improved the existing partonic Glauber model for proton-proton
collisions loi16 , including anisotropic and inhomogeneous density profiles
for the proton. They obtained a very good description of the $v_{2}$ measured
by the ALICE collaboration in pp collisions at $\sqrt{s}=13$ TeV alice19 .

Another interesting observation made in high multiplicity events refers to
charm production. In alice-cc  the ALICE collaboration measured the $D^{0}$
yield as a function of the central rapidity density and found out an
unexpectedly strong growth with the charged particle multiplicity. A similar
trend was observed in the case of $J/\psi$ production. There are already some
possible explanations of this growth of the charm yield given in Refs. cgc1 ,
boris20 , siddi20  and tripa22 .

In this work we will try to understand the data on charmonium production in
high multiplicity events using the same geometrical model proposed in deb20
with the parameters fixed in that work. Charmonium production will be computed
with the Color Evaporation Model (CEM). Leaving aside the quantitative
aspects, the idea is fairly simple and it is as follows. The proton-proton
collisions can be described by Y-Y collisions. Since we are not interested in
the anisotropy aspects of the collision, we can take the average configuration
of the Y over different orientations, i.e. we “rotate” the Y, obtaining a
circular configuration with an inner gluonic shell and an outer quark shell.
In other words, we have a core (gluon) - corona (quark) model of the proton.
Going from peripheral to more central and then to ultra-central, we go from
quark-quark collisions to gluon-gluon collisions. Since gluons are much more
abundant, and since $\sigma(gg\to c\bar{c})\,>>\,\sigma(q\bar{q}\to c\bar{c})$
the cross sections grow strongly. These effects combined should explain the
growth seen in the data.

In the next section we make some remarks concerning the proton structure, in
Section III we review the version of the Glauber model adapted for proton-
proton collisions and calculate the basic quantities, which are $N_{part}$ and
$N_{coll}$. In Section IV we review the main formulas of the Color Evaporation
Model used to study charm production. In Section V we show the results and
compare them with data. Finally, some concluding remarks are presented in the
last section.

##  II Remarks on the proton structure

In textbooks thomson  (Chap. 7) we learn that in low energy elastic $e-p$
scattering we can determine the charge distribution of the proton. This is
done by introducing electric, $G_{E}(Q^{2})$, and magnetic, $G_{M}(Q^{2})$,
form factors and fitting the resulting (Rosenbluth) formula to the
experimental data. Then, in a very specific limit, when $Q^{2}<<m^{2}_{p}$, we
find that $Q^{2}\simeq\mathbf{q}^{2}$ and hence $G_{E}(Q^{2})\simeq
G_{E}(\mathbf{q}^{2})$. In this limit the electric form factor can be
interpreted as the Fourier transform of the charge distribution:

$G_{E}(\mathbf{{q^{2}}})=\int\rho(\mathbf{r})e^{i\mathbf{q}\cdot\mathbf{r}}d^{3}\mathbf{r}$
(1)

Hence, measuring $G_{E}(Q^{2})$ we find the $\rho$ shown in Fig. 1a. Moving
away from the very low $Q^{2}$ limit, we do not have a well defined
prescription to extract charge densities from the measured form factors (see,
however, lorce  for progress in this direction). From electron-proton high
energy scattering in the high $Q^{2}$ limit, we know thomson  (Chap. 8) that
the incoming probe (photon or gluon) will indentify pointlike particles
(partons) inside the proton, but we do not know how these partons are
distributed in the transverse plane.

In the parton model description of deep inelastic scattering, the parton
distribution functions depend on $Q^{2}$ and this dependence comes from the
solution of the DGLAP evolution equations. In gla11 ; gla12 ; gla17  a
different kind of evolution was proposed: the Renormalization Group Procedure
for Effective Particles (REGPEP). In this approach, increasing the resolution
scale we go from a configuration where the quarks are unresolved to a
configuration with three effective quarks (quarks + antiquarks + gluons)
disposed in the vertices of an equilateral triangle with a star-like juntion
between them. Hence at high $Q^{2}$ the charge density tends to be moved from
the origin, as shown qualitatively in Fig. 1b).

|  |
---|---|---
(a) |  | (b)
Figure 1: a) Charge distribution of the proton measured in elastic electron
proton scattering at low energies and $Q^{2}\simeq 0$. b) Charge distribution
suggested by the Y shape model of the proton.

The proton structure derived from the REGPEP approach was used to contruct a
model successfully applied to describe several aspects of proton-proton
collisions with high multiplicities kubi14 ; kubi15 ; gla16  and also inspired
the model presented in deb20  and used here. In the next section we will
discuss this model in more detail.

##  III The Glauber model for proton-proton collisions

In this section we briefly describe the model proposed in kubi15 ; kubi14  and
successfully applied in deb20 . It combines the parton spatial distribution in
the proton advanced in kubi15 ; kubi14  with the geometry machinery provided
by the Glauber model glau07 . The quarks are anisotropically distributed, at
the edges of the Y shape junction and the gluons are isotropically distributed
arround the center kubi15 ; kubi14 .

Originally, the Glauber model was designed to represent the geometry of heavy-
ions, including a dynamical component given by the nucleon cross sections. In
loi16  it was adapted to proton-proton collisions. The effective number of
partonic (subnucleonic) degrees of freedom was called $N_{c}$. The analysis of
data performed in loi16  indicated that this number is $N_{c}\simeq 3-10$.
Following loi16  and deb20  we will work with the effective number of partonic
degrees of freedom (here we follow deb20  and call it $N_{g}$) and reinterpret
$N_{part}$ and $N_{coll}$ as “number of partons that participate in the
collision” and “number of binary collisions between partons” respectively.

In order to represent the proton internal structure, the matter distribution
is given by kubi14 :

$\rho_{p}(\mathbf{r};\mathbf{r}_{1},\mathbf{r}_{2},\mathbf{r}_{3})=\sum_{i=1}^{3}\rho_{q}(\mathbf{r}-\mathbf{r}_{i})+\rho_{g}\bigg{(}\mathbf{r}-\sum_{i=1}^{3}\frac{\mathbf{r}_{i}}{3}\bigg{)}$
(2)

In this distribution, there are three effective quarks, called $i=1,2,3$, with
the gaussian distribution:

$\rho_{q}(r)=(1-\kappa)\frac{N_{g}}{3}\frac{e^{-r^{2}/2r_{q}^{2}}}{(2\pi)^{3/2}r_{q}^{3}}$
(3)

and a gluon gaussian distribution:

$\rho_{g}(r)=\kappa N_{g}\frac{e^{-r^{2}/2r_{g}^{2}}}{(2\pi)^{3/2}r_{g}^{3}}$
(4)

centered at the average coordinate of the quarks. In these distributions,
$N_{g}$ is the total number of partons in the proton, $\kappa$ is the fraction
of $N_{g}$ that corresponds to the gluon body, $r_{q}$ is the radius of the
effective quark and $r_{g}$ the radius of the gluon body. These parameters are
fixed and we take them from previous works with this model. In Fig. 2 we show
the contour plot of these distributions. Fig. 3 shows a projection of these
distributions and it represents what we would see moving (up and to the right)
from the center of the proton to the peak of the quark density (in Fig. 2a)
plotted together with the gluon density from Fig. 2b. It is an illustration of
the “core-corona” aspect of the model and helps to understand the results.

|  |
---|---|---
(a) |  | (b)
Figure 2: a) Quark and b) gluon distributions in the transverse plane. Figure
3: Quark and gluon distributions. Projection of Fig. 2 along a diagonal
direction, starting from the center of the proton and passing through one
maximum of the quark density.

With the quark and gluon densities we can calculate the inputs to be used in
the Glauber model. The proton thickness, $T_{p}(x,y)$, is given by:

$T_{p}(x,y)=\int\rho_{p}(x,y,z)dz=T_{p}^{q}(x,y)+T_{p}^{g}(x,y)$ (5)

where $T_{p}^{q}(x,y)$ and $T_{p}^{g}(x,y)$ come from the quark and gluon
terms of Eq.(2) and are given by:

$T_{p}^{q}(x,y)=\frac{N_{g}(1-\kappa)}{6\pi
r_{q}^{2}}\sum_{i=1}^{3}e^{-\frac{(x-x_{i})^{2}+(y-y_{i})^{2}}{2r_{q}^{2}}}$
(6)

and

$T_{p}^{g}(x,y)=\frac{N_{g}\kappa}{2\pi
r_{g}^{2}}e^{-\frac{(x-\sum_{i=1}^{3}x_{i}/3)^{2}+(y-\sum_{i=1}^{3}y_{i}/3)^{2}}{2rg^{2}}}$
(7)

The overlap function can be calculated as:

$T_{pp\prime}(b)=\int T_{p}(x-b/2,y)T_{p\prime}(x+b/2,y)dxdy$ (8)

The above equation can be split into four terms. Because of the gaussian forms
involved, all the integrals can be done analytically. The quark-quark term is
given by:

$T_{pp\prime}^{qq}(b)=\int
T_{p}^{q}(x-b/2,y)T_{p\prime}^{q}(x+b/2,y)dxdy=\frac{N_{g}^{2}(1-\kappa)}{36\pi
r_{q}^{2}}\sum_{i=1}^{3}\sum_{j=1}^{3}\exp\bigg{(}-\frac{\big{(}b+(x_{i}-x_{j})\big{)}^{2}+\big{(}y_{i}-y_{j}\big{)}^{2}}{4rq^{2}}\bigg{)}$
(9)

with the quarks $i=1,2,3$ from the proton $p$, and $j=1,2,3$ from the proton
$p\prime$. The gluon-gluon term can be calculated as:

$\displaystyle T_{pp\prime}^{gg}(b)$ $\displaystyle=$ $\displaystyle\int
T_{p}^{g}(x-b/2,y)T_{p\prime}^{g}(x+b/2,y)dxdy$ (10) $\displaystyle=$
$\displaystyle\frac{N_{g}^{2}\kappa^{2}}{4\pi
r_{g}^{2}}\exp\bigg{(}-\frac{\big{(}b+\sum_{i=1}^{3}x_{i}/3-\sum_{j=1}^{3}x_{j}/3\big{)}^{2}+\big{(}\sum_{i=1}^{3}y_{i}/3-\sum_{j=1}^{3}y_{j}/3\big{)}^{2}}{4rg^{2}}\bigg{)}$

The gluon-quark term, which accounts for the interactions between the gluons
from $p$ with quarks from $p\prime$, is given by:

$\displaystyle T_{pp\prime}^{gq}(b)$ $\displaystyle=$ $\displaystyle\int
T_{p}^{g}(x-b/2,y)T_{p\prime}^{q}(x+b/2,y)dxdy$ (11) $\displaystyle=$
$\displaystyle\frac{N_{g}^{2}\kappa(1-\kappa)}{6\pi(r_{q}^{2}+r_{g}^{2})}\sum_{j=1}^{3}\exp\bigg{(}-\frac{\big{(}b+\sum_{i=1}^{3}x_{i}/3-x_{j}\big{)}^{2}+\big{(}\sum_{i=1}^{3}y_{i}/3-y_{j}\big{)}^{2}}{2(r_{q}^{2}+r_{g}^{2})}\bigg{)}$

The last term refers to the interaction between the quarks from $p$ and the
gluons from $p\prime$. It is given by:

$\displaystyle T_{pp\prime}^{qg}(b)$ $\displaystyle=$ $\displaystyle\int
T_{p}^{q}(x-b/2,y)T_{p\prime}^{g}(x+b/2,y)dxdy$ (12) $\displaystyle=$
$\displaystyle\frac{N_{g}^{2}\kappa(1-\kappa)}{6\pi(r_{q}^{2}+r_{g}^{2})}\sum_{i=1}^{3}\exp\bigg{(}-\frac{\big{(}b+x_{i}-\sum_{j=1}^{3}x_{j}/3\big{)}^{2}+\big{(}y_{i}-\sum_{j=1}^{3}y_{j}/3\big{)}^{2}}{2(r_{q}^{2}+r_{g}^{2})}\bigg{)}$

The total thickness function of the model is obtained by inserting Eqs. (9),
(10), (11) and (12) into Eq.(8). It reads:

$T_{pp\prime}(b)=T_{pp\prime}^{qq}(b)+T_{pp\prime}^{gg}(b)+T_{pp\prime}^{gq}(b)+T_{pp\prime}^{qg}(b)$
(13)

The number of binary collisions is then given by:

$N_{coll}=T_{pp\prime}(b)\sigma^{pp}$ (14)

where $\sigma^{pp}$ is the parton-parton cross section, which is a parameter
often used in the literature. The obtained values in the wounded quark model
of Ref. sigpp1 , in the hydrodynmical analysis of Ref. sigpp2  and in the
parton transport model of Ref. sigpp3  are all in the range $\sigma^{pp}=5-10$
mb. The number of participants is calculated according to

$N_{part}(b)=N_{coll}^{x}(b)$ (15)

where $x=0.75$ glau07 . In the above equations, $N_{part}$ and $N_{coll}$ are
calculated for a given quark configuration, i.e. a given choice of the quark
coordinates in the projectile proton $\mathbf{r}_{i}=(x_{i},y_{i})$ and in the
target proton $\mathbf{r^{\prime}}_{i}=({x^{\prime}}_{i},{y^{\prime}}_{i})$.
In order to simulate a proton-proton collision we must now take the average
over these configurations. This is best done writting the quark position in
polar coordinates:

$\mathbf{r}_{i}=\frac{d}{2}\big{(}\cos(\phi_{i}+\alpha),\sin(\phi_{i}+\alpha)\big{)}\,\,\,\,\,\,\,\,\,\,\,\,\,\mbox{and}\,\,\,\,\,\,\,\,\,\,\,\,\,\mathbf{r^{\prime}}_{i}=\frac{d}{2}\big{(}\cos(\phi_{i}+\beta),\sin(\phi_{i}+\beta)\big{)}$
(16)

where $\phi_{1}=\pi/3$, $\phi_{2}=-\pi/3$ and $\phi_{3}=-\pi$. This choice of
positions generates protons in which the quarks are in the vertices of
equilateral triangles in the $x-y$ plane, rotated by the angles $\alpha$ and
$\beta$ ($\alpha\text{, }\beta\in[0,2\pi]$) around the $z$ axis. The angles
are randomly chosen. For $d=1.3\text{ fm}$, one example of initial condition
is shown in Fig. 4a.

|  |
---|---|---
(a) |  | (b)
Figure 4: a) Example of initial spatial configuration of the two colliding
protons. b) Average of $N_{part}$ and $N_{coll}$ over different quark spatial
configurations for each impact parameter.

Using the parameters of Table 1 and fixing the impact parameter we choose the
angles $\alpha$ and $\beta$. With them we calculate the coordinates
$\mathbf{r}_{i}$ and then the thickness function, which we substitute into
Eqs. (14) and (15) obtaining $N_{coll}$ and $N_{part}$. We repeat this
procedure for 10000 choices and take the average. After that, we move to the
next impact parameter and repeat the steps. In the end we obtain $N_{coll}$
and $N_{part}$ as a function of the impact parameter. The results for $pp$
collisions at $\sqrt{s}=7\text{ TeV}$ are shown in Fig. 4b. Here the impact
parameter is such that $b\in[0,2.2]\text{ fm}$, with steps $\Delta b=0.1\text{
fm}$. The same procedure was repeated for $\sqrt{s}=13\text{ TeV}$. The
resulting curves are quite similar to those shown in Fig. 4b but with higher
maximum values. For conciseness we do not show them here.

Table 1: Input used in the parton densities and in the Glauber model. See text for definitions. $N_{g}$ deb20  | $10$
---|---
$\kappa$ deb20  | $0.5$
$r_{q}\text{ (fm)}$deb20  | $0.25$
$r_{g}\text{ (fm)}$ kubi14  | $0.5$
$\sigma^{pp}\text{ (mb)}$ ($\sqrt{s}=7$ TeV) sigpp1 ; sigpp2  | $4.3$
$\sigma^{pp}\text{ (mb)}$ ($\sqrt{s}=13$ TeV) sigpp1 ; sigpp2  | $7.6$
$d\text{ (fm)}$ deb20  | 1.3
$x$ glau07  | 0.75

##  IV Charmonium production

Charm production can be described by perturbative QCD and there are currently
several calculations which reproduce the data with very good accuracy. For our
purposes a leading order calculation is sufficient. Therefore we shall employ
the Color Evaporation Model vogt18 ; vogt07  (CEM). This model provides a way
to calculate the production of charm pairs through the processes
$gg\rightarrow c\bar{c}$ and $q\bar{q}\rightarrow c\bar{c}$. The cross section
for $J/\psi$ production is then given simply by:

$\displaystyle\sigma^{CEM}$ $\displaystyle=$
$\displaystyle\mathbf{F}K\sum_{i,j}\int_{(2m_{c})^{2}}^{(\Lambda)^{2}}dm^{2}\int
dx_{1}dx_{2}f_{i}(x_{1},\mu^{2})f_{j}(x_{2},\mu^{2})\sigma_{ij}(m^{2})\delta(m^{2}-x_{1}x_{2}s)$
(17) $\displaystyle=$ $\displaystyle\sigma_{gg}^{CEM}+\sigma_{q\bar{q}}^{CEM}$

where $f(x,\mu^{2})$ are the parton distribution functions at the
renormalization scale $\mu$ and $\Lambda$ is a cut-off. In the case of open
charm production $\Lambda=\sqrt{s}$ and for $J/\psi$ production
$\Lambda=2\,m_{D}$. The parameter $F$ is equal to the percentage of the
$c\bar{c}$ states with $2\,m_{c}<m<2\,m_{D}$ which becomes a $J/\psi$. We will
assume that $m_{c}=1200\text{ MeV}$ and $m_{D}=1800\text{ MeV}$. The symbol
$\sigma_{ii}(m^{2})$ represents the elementary $gg\rightarrow c\bar{c}$ and
$q\bar{q}\rightarrow c\bar{c}$ cross sections. At leading order, they are
given by:

$\sigma_{gg}(m^{2})=\frac{\pi\alpha_{s}^{2}(m^{2})}{3m^{2}}\bigg{\\{}\big{(}1+\frac{4m_{c}^{2}}{m^{2}}+\frac{m_{c}^{4}}{m^{4}}\big{)}\ln\bigg{(}\frac{1+\lambda}{1-\lambda}\bigg{)}-\frac{1}{4}\big{(}7+\frac{31m_{c}^{2}}{m^{2}}\big{)}\lambda\bigg{\\}}$
(18)

and

$\sigma_{q\bar{q}}(m^{2})=\frac{8\pi\alpha_{s}^{2}(m^{2})}{27m^{2}}\bigg{(}1+\frac{2m_{c}^{2}}{m^{2}}\bigg{)}\lambda$
(19)

where $\alpha_{s}$ (the strong coupling constant) and $\lambda$ are given by:

$\alpha_{s}(\mu^{2})=\frac{12\pi}{(33-2N_{f})\ln\big{(}\frac{\mu^{2}}{\Lambda_{QCD}^{2}}\big{)}}\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\lambda=\bigg{(}1-\frac{(2m_{c})^{2}}{m^{2}}\bigg{)}^{1/2}$
(20)

where $N_{f}$ is the number of flavors and $\Lambda_{QCD}=200\text{ MeV}$. The
parameter $K$ is introduced to account for higher order corrections. We can
fix it setting $F=1$, $\Lambda=\sqrt{s}$ and adjusting the cross section Eq.
(17) to the experimental data on open charm production phenix02 ; phenix06 ;
alice12 . The result is shown in Fig. 5. We obtain a good fit of these data
with $K=2$. Imposing that $\lambda$ is real leads to the kinematical
constraint $x_{1}x_{2}s\geq 4m_{c}^{2}$. Since the smallest value of $x_{1}$
occurs for $x_{2}=1$, the parton momentum fraction must be such that
$(2m_{c})^{2}/s\leq x_{1}\leq 1$ and $(2m_{c})^{2}/(sx_{1})\leq x_{2}\leq 1$.
The parton distribution functions, used with $Q^{2}=2.4^{2}\text{ GeV}^{2}$,
are taken from the CTEQ5 set PDFs . Finally, the number of $c\bar{c}$ pairs is
given by:

$N_{c\bar{c}}(b)=T_{pp\prime}^{gg}(b)\,\sigma_{gg}^{CEM}+T_{pp\prime}^{qq}(b)\,\sigma_{qq}^{CEM}$
(21)

The $T_{pp\prime}^{qg}(b)$ and $T_{pp\prime}^{gq}(b)$ terms of Eq.(13) are not
included because the quark-gluon interactions do not produce $c\bar{c}$ pairs.

Figure 5: Charm production cross sections. Contribution of $gg\rightarrow
c\bar{c}$ (dashed line), $q\bar{q}\rightarrow c\bar{c}$ (dotted line) and
total cross-section (solid line). The experimental data are from phenix02 ;
phenix06 ; alice12 .

##  V Results

We now address the experimental data of alice-psi7 ; alice-psi13 ; alice-cc .
Following the experimental papers, we present the yields in terms of the
normalized pseudo-rapidity and rapidity densities:

$\frac{({dN}/{d\eta}(\eta=0))}{\langle{dN}/{d\eta}\rangle}\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\frac{({dN}/{dy}(y=0))}{\langle{dN}/{dy}\rangle}$
(22)

The charged particle pseudorapidity density at $\eta=0$ is calculated as in
kn01 :

$\frac{dN}{d\eta}(\eta=0)=n_{pp}(s)\big{\\{}(1-f)\frac{N_{part}}{2}+fN_{coll}\big{\\}}$
(23)

where $N_{part}$ and $N_{coll}$ are given by (15) and (14) respectively and
$n_{pp}(s)$ is given by kn01 :

$n_{pp}(s)=2.5-0.25\log[s]+0.023(\log[s])^{2}$ (24)

The factor $f$ is the fraction of hard processes in the collision and can be
estimated through the ratio: $f=\sigma_{minijet}/\sigma_{inel}$ as in mini1 ;
mini2 . A minijet is defined as the result of a parton-parton collision with
$p_{T}>p_{T0}$ with $p_{T0}$ being of the order of a few GeV (with a possible
dependence on $\sqrt{s}$). The average density $\langle dN/d\eta\rangle$ is a
number given by each experimental group and we show them in Table II, together
with the $F$ parameter.

Table 2: Input used for $J/\psi$ production. See text for definitions. |  $\sqrt{s}=7$ TeV |  $\sqrt{s}=13$ TeV
---|---|---
f | 19 % mini1 ; mini2  | 16 % mini1 ; mini2
F | 1.5 % vogt18 ; vogt07  | 4.3 % vogt18 ; vogt07
$\langle dN/d\eta\rangle$ |  $6.1$ alice-dn-13  |  $6.4$ alice-dn-13
$\langle dN_{J/\psi}/dy\rangle$ |  $8.2\times 10^{-5}$ alice-psi7  |  $7.9\times 10^{-5}$ alice-psi13

The $J/\psi$ rapidity density, $dN_{J/\psi}/dy$, is obtained from the
differential form of Eq.(21):

$\frac{dN_{c\bar{c}}}{dy}(b)=T_{pp\prime}^{gg}(b)\,\frac{d\sigma_{gg}^{CEM}}{dy}+T_{pp\prime}^{qq}(b)\,\frac{d\sigma_{qq}^{CEM}}{dy}$
(25)

We can evaluate this expression starting from Eq.(17) and applying the
following change of variables:

$x_{1}=\frac{p_{T}}{\sqrt{s}}e^{y}\,\,\,\,\,\,\,\,\,\,\,\,x_{2}=\frac{p_{T}}{\sqrt{s}}e^{-y}\,\,\,\,\,\,\,\,\,\,\,\,x_{1}x_{2}s=p_{T}^{2}\,\,\,\,\,\,\,\,\,\,\,\,dx_{1}dx_{2}\rightarrow\frac{2p_{T}}{s}dydp_{T}$
(26)

After changing the variables from $(x_{1},x_{2})$ to $(y,p_{T})$, we integrate
over $p_{T}$, differentiate with respect to $y$ and take $y=0$. The results
for the $J/\psi$ yields for $\sqrt{s}=7\text{ TeV}$ and $\sqrt{s}=13\text{
TeV}$ are shown in Fig. 6. As it can be seen we obtain a very good description
of data. One might argue that we have too many input parameters and this
reduces the predictive power of the model. However, all these input numbers
are strongly constrained by other independent studies of other observables. As
we notice in the tables, all the numbers are consistent with the equivalent
numbers found elesewhere. We took the precaution of nowhere using exotic
values for any of these numbers. In view of the results, we conclude that the
$J/\psi$ dependence on the charged multiplicity is compatible with the
geometrical picture of the proton used here. The analysis can be extended to
$D$ production, to beauty production and to different cuts in $p_{T}$ and
rapidity. We are already working in these topics.

|  |
---|---|---
(a) |  | (b)
Figure 6:  a) $J/\psi$ relative yield in $pp$ collisions at
$\sqrt{s}_{NN}=7\text{ TeV}$. The experimental data are from alice-psi7 . b)
$J/\psi$ relative yield in $pp$ collisions at $\sqrt{s}_{NN}=13\text{ TeV}$.
The experimental data are from alice-psi13 .

##  VI Concluding remarks

We have developed the idea that proton-proton collisions can be described by
Y-Y collisions. Averaging over orientations this yields a circular
configuration with an inner gluonic shell and an outer quark shell. This is a
core (gluon) - corona (quark) model of the proton. Going from peripheral to
more central and then to ultra-central, we go from quark-quark collisions to
gluon-gluon collisions. Since gluons are much more abundant, and since
$\sigma(gg\to c\bar{c})\,>>\,\sigma(q\bar{q}\to c\bar{c})$ the cross sections
grow strongly. These effects combined explain the growth seen in the data.

This behavior was qualitatively expected and we have used a combination of
models to implement this idea quantitatively: a model for the parton spatial
densities deb20 ; the Glauber model for proton-proton collisions loi16  and
the color evaporation model for charmonium production vogt18 ; vogt07 . These
models had been previously tested in other contexts and were shown to give
good results. All the parameters used were strongly constrained by the
analysis of other data and here we had little room to change them. We obtained
a good description of charmonium production in proton-proton events with high
multiplicities. This encourages us to further extend this model to open charm
production and to proton-nucleus collisions.

###### Acknowledgements.

We are grateful to K. Werner for instructive discussions. This work was
financed by the Brazilian funding agencies CNPq, FAPESP and by the INCT-FNA.

## References

  * (1) B. Abelev et al. [ALICE], Phys. Lett. B 712, 165 (2012).
  * (2) J. Adam et al. [ALICE Collaboration], JHEP 09, 148 (2015).
  * (3) J. Adam et al. [STAR], Phys. Lett. B 786, 87 (2018).
  * (4) K. Gajodsov, [ALICE Collaboration], Nucl. Phys. A 982, 487 (2019).
  * (5) S. Acharya et al. [ALICE], Eur. Phys. J. C 80, 167 (2020).
  * (6) S. Acharya et al. [ALICE], Phys. Lett. B 810, 135758 (2020).
  * (7) S. Acharya et al. [ALICE], JHEP 06, 015 (2022).
  * (8) N. Sakumichi and H. Suganuma, Phys. Rev. D 92, 034511 (2015).
  * (9) H. Mantysaari, B. Schenke, C. Shen and W. Zhao, arXiv:2303.04866.
  * (10) G. N. Fowler, F. S. Navarra, M. Plumer, A. Vourdas, R. M. Weiner and G. Wilk, Phys. Rev. C 40, 1219 (1989).
  * (11) D. Kharzeev, Phys. Lett. B 378, 238 (1996),
  * (12) V. Topor Pop, M. Gyulassy, J. Barrette, C. Gale, X. N. Wang and N. Xu, Phys. Rev. C 70, 064906 (2004).
  * (13) F. W. Bopp and M. Shabelski, Eur. Phys. J. A 28, 237 (2006).
  * (14) J. D. Brandenburg, N. Lewis, P. Tribedy and Z. Xu, [arXiv:2205.05685 [hep-ph]].
  * (15) S. Deb, G. Sarwar, D. Thakur, P. Subramani, R. Sahoo and J. e. Alam, Phys. Rev. D 101, 014004 (2020).
  * (16) C. Loizides, Phys. Rev. C 94, 024914 (2016).
  * (17) Y.-Q. Ma, P. Tribedy, R. Venugopalan and K. Watanabe Phys. Rev. D 98, 074025 (2018).
  * (18) B. Z. Kopeliovich, H. J. Pirner, I. K. Potashnikova, K. Reygers and I. Schmidt, Phys. Rev. D 101, 054023 (2020).
  * (19) E. Levin, I. Schmidt and M. Siddikov, Eur. Phys. J. C 80, 560 (2020).
  * (20) T. Tripathy et al., arXiv:2209.08784.
  * (21) M. Thomson, “Modern Particle Physics”, Cambridge University Press, (2013).
  * (22) C. Lorcé, Phys. Rev. Lett. 125, 232002 (2020).
  * (23) S. Glazek, Few-Body Syst. 52, 367, (2012).
  * (24) S. D. Glazek, Phys. Rev. D 85, 125018 (2012).
  * (25) S. D. Głazek and A. P. Trawiński, Few Body Syst. 58, 49 (2017).
  * (26) P. Kubiczek, “Geometrical model of azimuthal correlations in high-multiplicity proton-proton collisions”, Thesis, University of Warsaw, (2014).
  * (27) P. Kubiczek and S. D. Głazek, Lith. J. Phys. 55, 155 (2015).
  * (28) S. D. Glazek and P. Kubiczek, Few Body Syst. 57, 425 (2016).
  * (29) M.L. Miller et al., Annu. Rev. Nucl. Part. Sci. 57, 205 (2007).
  * (30) P. Bożek, W. Broniowski, M. Rybczyński, Phys. Rev. C 94, 014902 (2016).
  * (31) H. Drescher, A. Dumitru, C. Gombeaud, J. Ollitraut, Phys. Rev. C 76, 024905 (2007).
  * (32) Z. Xu and C. Greiner, Nucl. Phys. A 774, 787 (2006)
  * (33) V. Cheung and R. Vogt, Phys. Rev. D 98, 114029 (2018);
  * (34) R. Vogt, “Ultrarelativistic heavy-ion collisions”, Elsevier, (2007).
  * (35) K. Adcox et al. [PHENIX Collaboration], Phys. Rev. Lett. 88, 192303 (2002).
  * (36) A. Adare et al. [PHENIX Collaboration], Phys. Rev. Lett. 97, 252002 (2006).
  * (37) B. Abelev et al. [ALICE Collaboration], JHEP 1207, 191 (2012).
  * (38) H.L. Lai et al. [CTEQ Collaboration], Eur. Phys. J. C 12, 375 (2000).
  * (39) D. Kharzeev, M. Nardi, Phys. Lett. B 507, 121 (2001).
  * (40) P. Kotko, A. M. Stasto and M. Strikman, Phys. Rev. D 95, 054009 (2017).
  * (41) I. Sarcevic, S. D. Ellis and P. Carruthers, Phys. Rev. D 40, 1446 (1989).
  * (42) J. Adam et al. [ALICE Collaboration], Phys. Lett. B 753, 319 (2016).