#  An extension of Schur’s irreducibility result

ANKITA JINDAL111Stat-Math Unit, ISI Bangalore Centre, Bangalore, India. Email
<EMAIL_ADDRESS>AND SUDESH KAUR KHANDUJA222Corresponding author
333Indian Institute of Science Education and Research Mohali, SAS Nagar, India
& Department of Mathematics, Panjab University, Chandigarh, India. Email :
<EMAIL_ADDRESS>
(Dedicated to Professor André Leroy on his $65^{th}$ birthday)

Abstract

Let $n\geq 2$ be an integer. Let $\phi(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is irreducible modulo all
primes less than or equal to $n$. Let $a_{0}(x),a_{1}(x),\dots,a_{n-1}(x)$
belonging to $\mbox{$\mathbb{Z}$}[x]$ be polynomials each having degree less
than $\deg\phi(x)$ and $a_{n}$ be an integer. Assume that $a_{n}$ and the
content of $a_{0}(x)$ are coprime with $n!$. In the present paper, we prove
that the polynomial
$\sum\limits_{i=0}^{n-1}a_{i}(x)\frac{\phi(x)^{i}}{i!}+a_{n}\frac{\phi(x)^{n}}{n!}$
is irreducible over the field $\mathbb{Q}$ of rational numbers. This
generalizes a well known result of Schur which states that the polynomial
$\sum\limits_{i=0}^{n}a_{i}\frac{x^{i}}{i!}$ is irreducible over $\mathbb{Q}$
for all $n\geq 1$ when each $a_{i}\in\mbox{$\mathbb{Z}$}$ and
$|a_{0}|=|a_{n}|=1$. The present paper also extends a result of Filaseta
thereby leading to a generalization of the classical Schönemann Irreducibility
Criterion.

Keywords : Irreducible polynomials, Truncated exponential series.

2020 Mathematics Subject Classification : 11C08, 11R04.

##  1 Introduction

A long established theorem of Schur [17] states that the polynomial
$\sum\limits_{i=0}^{n}a_{i}\frac{x^{i}}{i!}$ is irreducible over the field
$\mathbb{Q}$ of rational numbers for all $n\geq 1$ when each
$a_{i}\in\mbox{$\mathbb{Z}$}$ and $|a_{0}|=|a_{n}|=1$. An alternate proof of
this result when $|a_{i}|=1$ for $i\in\\{0,1,\dots,n\\}$ was given by Coleman
using Newton polygons with respect to primes (cf. [3]). In the present paper,
we extend the result of Schur as well as of Coleman in different directions
using $\phi$-Newton polygons defined in this section and prove the following.

######  Theorem 1.1.

Let $n\geq 2$ be an integer. Let $\phi(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is irreducible modulo all
primes less than or equal to $n$. Let $a_{0}(x),a_{1}(x),\dots,a_{n-1}(x)$
belonging to $\mbox{$\mathbb{Z}$}[x]$ be polynomials each having degree less
than $\deg\phi(x)$ and $a_{n}$ be an integer. Assume that $a_{n}$ and the
content444The content of a polynomial with coefficients in the ring
$\mathbb{Z}$ of integers is the greatest common divisor of its coefficients.
of $a_{0}(x)$ are coprime with $n!$. Then the polynomial

$\sum\limits_{i=0}^{n-1}a_{i}(x)\frac{\phi(x)^{i}}{i!}+a_{n}\frac{\phi(x)^{n}}{n!}$

is irreducible over $\mathbb{Q}$.

It may be pointed out that in the above theorem the assumptions “$a_{n}$ and
the content of $a_{0}(x)$ are coprime with $n!$" cannot be dispensed with. For
example, consider the polynomial $\phi(x)=x^{2}+x+1$ which is irreducible
modulo 2. The polynomial
$F(x)=\frac{\phi(x)^{2}}{2!}+\phi(x)-4=\frac{1}{2}(\phi(x)+4)(\phi(x)-2)$ is
reducible over $\mbox{$\mathbb{Q}$}.$ Similarly the polynomial
$G(x)=12\frac{\phi(x)^{2}}{2!}-\phi(x)-1=(3\phi(x)+1)(2\phi(x)-1)$ is
reducible over $\mbox{$\mathbb{Q}$}.$

We also give below an example to show that Theorem 1.1 may not hold if $a_{n}$
is replaced by a (monic) polynomial $a_{n}(x)$ with integer coefficients
having degree less than $\deg\phi(x)$. Consider $\phi(x)=x^{2}+x+1$,
$a_{2}(x)=x+1$, $a_{1}(x)=x+2$ and $a_{0}(x)=4x+3$. Then

$a_{2}(x)\frac{\phi(x)^{2}}{2!}+a_{1}(x)\phi(x)+a_{0}(x)$

has $-1$ as a root.

Assuming a weaker condition on $\phi(x)$ in Theorem 1.1, we also prove the
following result using the Factorisation Theorem of Ore (cf. [13, Theorem 1.1]
and [14]).

######  Theorem 1.2.

Let $n\geq 2$ be an integer. Let $\phi(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is irreducible modulo all
primes dividing $n$. Let $a_{0}(x),a_{1}(x),\dots,a_{n-1}(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be polynomials each having degree less than
$\deg\phi$. Assume that the content of the polynomial
$\prod\limits_{i=0}^{n-1}a_{i}(x)$ is coprime with $n$. Then the polynomial

$\sum\limits_{i=0}^{n-1}a_{i}(x)\frac{\phi(x)^{i}}{i!}+\frac{\phi(x)^{n}}{n!}$

is irreducible over $\mathbb{Q}$.

We wish to point out that analogues of Theorems 1.1 and 1.2 do not hold for
$n=1$ because if $\phi(x)\in\mbox{$\mathbb{Z}$}[x]$ is a monic polynomial of
degree $m\geq 2$, then the polynomial $\phi(x)-(\phi(x)-x^{m})$ is reducible
over $\mathbb{Q}$.

In the proof of Theorem 1.1, a key role is played by Theorem 1.3 which extends
a well known result of Filaseta [6, Lemma 2]. For stating this theorem, we
first introduce the notion of $\phi$-Newton polygon with respect to a prime
number modulo which $\phi(x)\in\mbox{$\mathbb{Z}$}[x]$ is irreducible.

For a prime $p$, $v_{p}$ will denote the $p$-adic valuation of $\mathbb{Q}$
defined for any non-zero integer $c$ to be the highest power of $p$ dividing
$c$. We shall denote by $v_{p}^{x}$ the Gaussian valuation extending $v_{p}$
defined on the polynomial ring $\mbox{$\mathbb{Z}$}[x]$ by

$\displaystyle
v_{p}^{x}(\sum\limits_{i}c_{i}x^{i})=\min_{i}\\{v_{p}(c_{i})\\},\
c_{i}\in\mbox{$\mathbb{Z}$}.$

If $\phi(x)$ is a fixed monic polynomial with coefficients in $\mathbb{Z}$,
then any $f(x)\in\mbox{$\mathbb{Z}$}[x]$ can be uniquely written as a finite
sum $\sum\limits_{i}f_{i}(x)\phi(x)^{i}$ with $\deg f_{i}(x)<\deg\phi(x)$ for
each $i$; this expansion will be referred to as the $\phi$-expansion of
$f(x)$.

Let $\phi(x)\in\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is
irreducible modulo a given prime $p$. Let $f(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be a polynomial having $\phi$-expansion
$\sum\limits_{i=0}^{n}f_{i}(x)\phi(x)^{i}$ with $f_{0}(x)f_{n}(x)\neq 0$. Let
$P_{i}$ stand for the point in the plane having coordinates
$(i,v_{p}^{x}(f_{n-i}(x)))$ when $f_{n-i}(x)\neq 0$, $0\leq i\leq n$. Let
$\mu_{ij}$ denote the slope of the line joining the points $P_{i}$ and $P_{j}$
if $f_{n-i}(x)f_{n-j}(x)\neq 0$. Let $i_{1}$ be the largest index $0<i_{1}\leq
n$ such that

$\displaystyle\mu_{0i_{1}}=\min\\{\mu_{0j}\ |\ 0<j\leq n,\ f_{n-j}(x)\neq
0\\}.$

If $i_{1}<n$, let $i_{2}$ be the largest index $i_{1}<i_{2}\leq n$ satisfying

$\displaystyle\mu_{i_{1}i_{2}}=\min\\{\mu_{i_{1}j}\ |\ i_{1}<j\leq n,\
f_{n-j}(x)\neq 0\\}$

and so on. The $\phi$-Newton polygon of $f(x)$ with respect to $p$ is the
polygonal path having segments
$P_{0}P_{i_{1}},P_{i_{1}}P_{i_{2}},\dots,P_{i_{k-1}}P_{i_{k}}$ with $i_{k}=n$.
These segments are called the edges of the $\phi$-Newton polygon of $f(x)$ and
their slopes from left to right form a strictly increasing sequence. The
$\phi$-Newton polygon with respect to $p$ minus the horizontal part (if any)
is called its principal part.

With the above notation, we prove the following theorem, which was proved by
Filaseta [6, Lemma 2] in 1995 in the particular case when $\phi(x)=x$; the
proof given here is similar to the one given by him.

######  Theorem 1.3.

Let $n,k$ and $\ell$ be integers with $0\leq\ell<k\leq\frac{n}{2}$ and $p$ be
a prime. Let $\phi(x)\in\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is
irreducible modulo $p$. Let $f(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$ be a
monic polynomial not divisible by $\phi(x)$ having $\phi$-expansion
$\sum\limits_{i=0}^{n}f_{i}(x)\phi(x)^{i}$ with $f_{n}(x)\neq 0$. Assume that
$v_{p}^{x}(f_{i}(x))>0$ for $0\leq i\leq n-\ell-1$ and the right-most edge of
the $\phi$-Newton polygon of $f(x)$ with respect to $p$ has slope less than
$\frac{1}{k}$. Let $a_{0}(x),a_{1}(x),\dots,a_{n}(x)$ be polynomials over
$\mathbb{Z}$ satisfying the following conditions.

  * (i)

$\deg a_{i}(x)<\deg\phi(x)-\deg f_{i}(x)$ for $0\leq i\leq n$,

  * (ii)

$v_{p}^{x}(a_{0}(x))=0$, i.e., the content of $a_{0}(x)$ is not divisible by
$p$,

  * (iii)

the leading coefficient of $a_{n}(x)$ is not divisible by $p$.

Then the polynomial $\sum\limits_{i=0}^{n}a_{i}(x)f_{i}(x)\phi(x)^{i}$ does
not have a factor in $\mbox{$\mathbb{Z}$}[x]$ with degree lying in the
interval $[(\ell+1)\deg\phi,(k+1)\deg\phi).$

The classical Schönemann Irreducibility Criterion (cf. [16], [15, Chapter 3,
Theorem D]) stated below will be quickly deduced from the above theorem.

######  Corollary 1.4 (Classical Schönemann Irreducibility Criterion).

Let $\phi(x)\in\mbox{$\mathbb{Z}$}[x]$ be a monic polynomial which is
irreducible modulo a given prime $p$. Let $f(x)$ belonging to
$\mbox{$\mathbb{Z}$}[x]$ be of the form $f(x)=\phi(x)^{n}+pM(x)$ where
$M(x)\in\mbox{$\mathbb{Z}$}[x]$ has degree less than $n\deg\phi$. If $\phi(x)$
is coprime with $M(x)$ modulo $p$, then $f(x)$ is irreducible over
$\mathbb{Q}$.

##  2 Proof of Theorem 1.2

The following simple result is well known (see [2, p. 122], [3]). Its proof is
omitted.

######  Lemma 2.1.

Let $p$ be a prime and $m$ be a positive integer. If
$m=a_{0}+a_{1}p+\cdots+a_{r}p^{r}$ where $\ 0\leq a_{i}<p$ for each $i$, then

$\displaystyle v_{p}(m!)=\frac{m-(a_{0}+a_{1}+\cdots+a_{r})}{p-1}.$

###### Proof of Theorem 1.2.

Set $a_{n}(x)=1$ and consider the monic polynomial
$F(x)\in\mbox{$\mathbb{Z}$}[x]$ defined by

$F(x)=\sum\limits_{i=0}^{n}\frac{n!}{i!}a_{i}(x)\phi(x)^{i}.$

We shall prove that $F(x)$ is irreducible over $\mathbb{Q}$. Let $p$ be a
prime dividing $n$. By hypothesis, the content of $a_{i}(x)$ is coprime with
$p$ for $0\leq i\leq n-1$, i.e., $v_{p}^{x}(a_{i}(x))=0$ for $0\leq i\leq
n-1$. Recall that $a_{n}(x)=1$. Therefore the $\phi$-Newton polygon of $F(x)$
with respect to $p$ is the polygonal path formed by the lower edges along the
convex hull of points of the set $S$ defined by

$S=\\{(i,v_{p}(n!/(n-i)!))\ |\ 0\leq i\leq n\\}.$

Write

$n=c_{1}p^{m_{1}}+c_{2}p^{m_{2}}+\cdots+c_{s}p^{m_{s}},$

where $0<m_{1}<m_{2}<\cdots<m_{s}$, $0<c_{i}<p$ for each $i$. Set $z_{0}=0$
and

$z_{i}=c_{1}p^{m_{1}}+\cdots+c_{i}p^{m_{i}}$

for $1\leq i\leq s$. As in [3], making use of Lemma 2.1, it can be easily
shown that the polygonal path along the lower convex hull of points of the set
$S$ consists of $s$ edges. The $i^{th}$ edge from left to right is the segment
joining the points

$(z_{i-1},v_{p}(n!/(n-z_{i-1})!)),\ (z_{i},v_{p}(n!/(n-z_{i})!)).$

Again using Lemma 2.1, we see that the slope $\lambda_{i}$ of the $i^{th}$
edge of the $\phi$-Newton polygon of $F(x)$ is given by

$\displaystyle\lambda_{i}$
$\displaystyle=\frac{-v_{p}((n-z_{i})!)+v_{p}((n-z_{i-1})!)}{z_{i}-z_{i-1}}$
$\displaystyle=\frac{z_{i}+(c_{i+1}+\cdots+c_{s})-z_{i-1}-(c_{i}+\cdots+c_{s})}{(z_{i}-z_{i-1})(p-1)}.$

So

$\displaystyle\lambda_{i}=\frac{c_{i}p^{m_{i}}-c_{i}}{c_{i}p^{m_{i}}(p-1)}=\frac{p^{m_{i}}-1}{p^{m_{i}}(p-1)}.$
(1)

Since $0<m_{1}<\cdots<m_{s}$, it is immediate from (1) that $p^{m_{1}}$
divides the denominator of slope of each edge of the $\phi$-Newton polygon of
$F(x)$ with respect to $p$. Therefore in view of Ore’s Factorisation Theorem
[13, Theorem 1.1(ii)], $p^{m_{1}}\deg\phi$ divides the degree of each factor
of $F(x)$ over the field $\mbox{$\mathbb{Q}$}_{p}$ of $p$-adic numbers.
Keeping in mind that $p^{m_{1}}$ is the exact power of an arbitrary prime $p$
dividing $n$, the assertion proved above implies that $n\deg\phi$ divides the
degree of each irreducible factor of $F(x)$ over $\mathbb{Q}$. Hence $F(x)$
having degree $n\deg\phi$ is irreducible over $\mathbb{Q}$. ∎

######  Example 2.2.

Consider $\phi(x)=x^{3}-x^{2}+1$ or $x^{4}-x-1$. It can be easily checked that
in both situations $\phi(x)$ is irreducible modulo 2 as well as modulo 3. Let
$a_{0}(x),a_{1}(x),\dots,a_{5}(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$ be
polynomials with degree less than $\deg\phi$ and each having content coprime
with $6$. By Theorem 1.2, the polynomial

$\displaystyle
F(x)=\phi(x)^{6}+\frac{6!}{5!}a_{5}(x)\phi(x)^{5}+\frac{6!}{4!}a_{4}(x)\phi(x)^{4}+\frac{6!}{3!}$
$\displaystyle a_{3}(x)\phi(x)^{3}+\frac{6!}{2!}a_{2}(x)\phi(x)^{2}$
$\displaystyle+6!a_{1}(x)\phi(x)+6!a_{0}(x)$

is irreducible over $\mathbb{Q}$.

##  3 Proof of Theorem 1.3

We shall use the following lemma whose first assertion is proved in [13, Lemma
2.4].

Lemma 3.A. Let $\phi(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$ be a monic
polynomial which is irreducible modulo a given prime $p$. Let $g(x),h(x)$
belonging to $\mbox{$\mathbb{Z}$}[x]$ be polynomials not divisible by
$\phi(x)$ having leading coefficients coprime with $p$. Then the following
hold.

  * (i)

The principal part of the $\phi$-Newton polygon of $g(x)h(x)$ with respect to
$p$ is obtained by constructing a polygonal path beginning with a point on the
non-negative side of $x$-axis and using translates of edges of the principal
part in the $\phi$-Newton polygons of $g(x),h(x)$ in the increasing order of
slopes.

  * (ii)

The length of the edge of the $\phi$-Newton polygon of $g(x)h(x)$ with respect
to $p$ having slope zero is either the sum of the lengths of the edges of the
$\phi$-Newton polygons of $g(x)$, $h(x)$ with respect to $p$ having slope zero
or it exceeds this sum by one.

###### Proof.

We only prove assertion (ii). For any $f(x)\in\mbox{$\mathbb{Z}$}[x]$, let
$\bar{f}(x)$ denote the polynomial obtained by replacing each coefficient of
$f(x)$ modulo $p$. Let $g(x)=\sum\limits_{i=0}^{t}g_{i}(x)\phi(x)^{i}$ and
$h(x)=\sum\limits_{i=0}^{u}h_{i}(x)\phi(x)^{i}$ be the $\phi$-expansions of
$g(x),h(x)$ respectively with $g_{t}(x)h_{u}(x)\neq 0$. Let $r\geq 0$, $s\geq
0$ denote respectively the length of the edges having slope zero of the
$\phi$-Newton polygons of $g(x),h(x)$ with respect to $p$. Then

$\displaystyle v_{p}^{x}(g_{t-r}(x))=0,\ v_{p}^{x}(g_{t-i}(x))>0\textrm{ for
}i\in\\{r+1,\dots,t\\}$

and

$\displaystyle v_{p}^{x}(h_{u-s}(x))=0,\ v_{p}^{x}(h_{u-i}(x))>0\textrm{ for
}i\in\\{s+1,\dots,u\\}.$

So the highest power of $\bar{\phi}(x)$ dividing $\bar{g}(x)$, $\bar{h}(x)$ is
$t-r$ and $u-s$ respectively. Consequently the highest power of
$\bar{\phi}(x)$ dividing $\bar{g}(x)\bar{h}(x)$ is $t+u-r-s$. Note that if we
write the $\phi$-expansion of $g(x)h(x)$ as
$\sum\limits_{i=0}^{k}d_{i}(x)\phi(x)^{i}$ with $d_{k}(x)\neq 0$, $\deg
d_{i}(x)<\deg\phi(x)$, then either $k=t+u$ or $k=t+u+1$. Therefore the length
of the edge having slope zero in the $\phi$-Newton polygon of $g(x)h(x)$ is
either $r+s$ or $r+s+1$; this proves assertion (ii) of the lemma. ∎

###### Proof of Theorem 1.3.

Write

$\displaystyle F(x)=\sum\limits_{i=0}^{n}a_{i}(x)f_{i}(x)\phi(x)^{i}.$

We first consider the case when $a_{i}(x)=1$ for all $i\in\\{0,1,\dots,n\\}$
so that $F(x)=f(x)$. Recall that $f(x)\in\mbox{$\mathbb{Z}$}[x]$ is monic.
Suppose to the contrary that $f(x)$ has a factor over $\mathbb{Z}$ with degree
lying in the interval $[(\ell+1)\deg\phi,(k+1)\deg\phi)$. Then there exist
monic polynomials $g(x),h(x)\in\mbox{$\mathbb{Z}$}[x]$ with $f(x)=g(x)h(x)$
and

$\displaystyle(\ell+1)\deg\phi\leq\deg g<(k+1)\deg\phi.$ (2)

Let $\sum\limits_{i=0}^{r}g_{i}(x)\phi(x)^{i}$ be the $\phi$-expansion of
$g(x)$ with $g_{r}(x)\neq 0$. Note that $r$ is the sum of lengths of
horizontal projections of all edges of the $\phi$-Newton polygon of $g(x)$
with respect to $p$. We have

$\displaystyle r\deg\phi\leq\deg g=\deg g_{r}+r\deg\phi\leq(r+1)\deg\phi-1.$
(3)

The above inequality together with (2) implies that

$\displaystyle r+1\geq\frac{\deg
g+1}{\deg\phi}\geq\frac{(\ell+1)\deg\phi+1}{\deg\phi}=\ell+1+\frac{1}{\deg\phi},$

which shows that

$\displaystyle r\geq\ell+\frac{1}{\deg\phi}>\ell.$ (4)

We now consider the $\phi$-Newton polygon of $f(x)$ (with respect to $p$). By
hypothesis, each edge of this Newton polygon has slope in the interval
$[0,\frac{1}{k})$. For now, we consider an edge of the $\phi$-Newton polygon
of $f(x)$ which has positive slope. Let $(a,b)$ and $(c,d)$ be two points with
integer entries on such an edge. Then slope of the line joining these points
is the slope of that edge, so that

$\displaystyle\frac{1}{|c-a|}\leq\frac{|d-b|}{|c-a|}<\frac{1}{k}.$

Hence $|c-a|>k$. Therefore any two such points on an edge with positive slope
of the $\phi$-Newton polygon of $f(x)$ have their $x$-coordinates separated by
a distance strictly greater than $k$.

Note that in view of (3) and (2), we have

$\displaystyle r\leq\frac{\deg g}{\deg\phi}<k+1.$

The above inequality implies that $r$, which is the sum of lengths of
horizontal projections of all edges of the $\phi$-Newton polygon of $g(x)$, is
less than $k+1$. So by what has been proved in the above paragraph, the
translates of the edges of the $\phi$-Newton polygon of $g(x)$ cannot be found
within those edges of the $\phi$-Newton polygon of $f(x)$ which have positive
slope. Therefore by Lemma 3.A, the left-most edge of the $\phi$-Newton polygon
of $f(x)$ must have slope zero and its length is greater than or equal to $r$.
In view of the hypothesis

$\displaystyle v_{p}^{x}(f_{n-i}(x))>0\textrm{ for }i\in\\{\ell+1,\dots,n\\},$

we see that the length of the edge of the $\phi$-Newton polygon of $f(x)$
having slope zero is less than or equal to $\ell$. On recalling that $\ell<r$
by virtue of (4), we obtain a contradiction in view of Lemma 3.A and hence the
theorem is proved when all $a_{i}(x)$ are 1.

Next we consider the general case when

$\displaystyle F(x)=\sum\limits_{i=0}^{n}a_{i}(x)f_{i}(x)\phi(x)^{i},$

where the polynomials $a_{i}(x)$ satisfy the conditions $(i),(ii)$ and
$(iii)$. It is clear that the left-most endpoint and the right-most endpoint
of the $\phi$-Newton polygon of $F(x)$ are the same as the left-most endpoint
and the right-most endpoint of the $\phi$-Newton polygon of $f(x)$. Recall
that

$\displaystyle v_{p}^{x}(a_{i}(x)f_{i}(x))\geq v_{p}^{x}(f_{i}(x))$
$\displaystyle>0\quad\quad\textrm{ for }i\in\\{0,1,\dots,n-\ell-1\\},$

and

$\displaystyle v_{p}^{x}(a_{i}(x)f_{i}(x))\geq v_{p}^{x}(f_{i}(x))$
$\displaystyle\geq 0\quad\quad\textrm{ for }i\in\\{0,1,\dots,n\\}.$

Therefore all the edges of the $\phi$-Newton polygon of $F(x)$ lie above or on
the line containing the right-most edge of the $\phi$-Newton polygon of
$f(x)$. Since the right-most endpoints of the $\phi$-Newton polygons of $F(x)$
and $f(x)$ are the same, we deduce that the slope of the right-most edge of
the $\phi$-Newton polygon of $F(x)$ is less than or equal to the slope of the
right-most edge of the $\phi$-Newton polygon of $f(x)$ and hence the former is
also less than $\frac{1}{k}$. Recall that the leading coefficient of
$a_{n}(x)$ is coprime with $p$ by hypothesis. Thus by appealing to the first
part of the proof and applying Lemma 3.A, Theorem 1.3 follows. ∎

######  Remark 3.1.

It may be pointed out that the proof of Theorem 1.3 carries over verbatim when
the polynomials $f(x)$, $a_{0}(x)$, $\dots$, $a_{n}(x)$ have coefficients in
the valuation ring $R_{v}$ of a discrete valuation $v$ having value group
$\mathbb{Z}$. Note that Lemma 3.A used in the above proof also remains valid
for polynomials with coefficients in a discrete valuation ring $R_{v}$ having
leading coefficient a unit in $R_{v}$. Indeed assertion $(i)$ of Lemma 3.A is
already proved for such polynomials in [13, Lemma 2.4], whereas the proof of
assertion $(ii)$ given here carries over verbatim in the case of such
polynomials.

###### Proof of Corollary 1.4.

Let $f(x)=\sum\limits_{i=0}^{n}f_{i}(x)\phi(x)^{i}$ be the $\phi$-expansion of
$f(x)$. In view of the hypothesis, we have $f_{n}(x)=1$,
$v_{p}^{x}(f_{i}(x))>0$ for $0\leq i\leq n-1$ and $v_{p}^{x}(f_{0}(x))=1$. So
the $\phi$-Newton polygon of $f(x)$ with respect to $p$ consists of a single
edge joining the points $(0,0)$ and $(n,1)$ which has slope $\frac{1}{n}$.
Choose $\ell=0$ and $k$ to be the largest integer not exceeding $\frac{n}{2}$.
Applying Theorem 1.3, we conclude that $f(x)$ does not have a factor in
$\mbox{$\mathbb{Z}$}[x]$ with degree in the interval
$[\deg\phi,(k+1)\deg\phi)$.

Suppose to the contrary that $f(x)$ is reducible over $\mathbb{Q}$. Then
$f(x)=g(x)h(x)$ where $g(x)$, $h(x)$ are monic polynomials belonging to
$\mbox{$\mathbb{Z}$}[x]$ of positive degrees. Therefore on passing to
$\mbox{$\mathbb{Z}$}/p\mbox{$\mathbb{Z}$}$, we have
$\bar{f}(x)=\bar{g}(x)\bar{h}(x)$. By hypothesis,
$\bar{f}(x)=\bar{\phi}(x)^{n}$. Since $\bar{\phi}(x)$ is irreducible over
$\mbox{$\mathbb{Z}$}/p\mbox{$\mathbb{Z}$}$, it follows that
$\bar{g}(x)=\bar{\phi}(x)^{d}$ and $\bar{h}(x)=\bar{\phi}(x)^{e}$ for some
positive integers $d$ and $e$. In particular, both $g(x)$ and $h(x)$ have
degrees greater than or equal to the degree of $\phi(x)$. Since one of $g(x)$
or $h(x)$ has degree strictly less than $(k+1)\deg\phi$, we obtain a
contradiction in view of what has been proved above. ∎

##  4 Proof of Theorem 1.1

The following result to be used in the sequel was proved independently by
Sylvester [19] and Schur [18]. This was reproved by Erdős in 1934 (cf. [5]).

Theorem 4.A. Let $m$ and $k$ be positive integers with $m\geq k$. Then there
is a prime $p\geq k+1$ which divides one of the integers $m+1,m+2,\dots,m+k$.

###### Proof of Theorem 1.1.

Consider the polynomial
$F(x)=\sum\limits_{i=0}^{n-1}\frac{n!}{i!}a_{i}(x)\phi(x)^{i}+a_{n}\phi(x)^{n}$
with coefficients in $\mathbb{Z}$. In view of the hypothesis, $a_{0}(x)$ is
non-zero. So $F(x)$ is not divisible by $\phi(x)$. Since $\deg
a_{i}(x)<\deg\phi(x)$ for each $i$, it follows that the leading coefficient of
$F(x)$ is $a_{n}$, which is coprime with $n!$ by hypothesis. Thus the content
$c$ (say) of $F(x)$ is coprime with $n!$. It suffices to prove that
$\frac{F(x)}{c}$ is irreducible over $\mathbb{Z}$. We divide the proof of the
theorem in two steps.

Step I. In this step, we show that $\frac{F(x)}{c}$ does not have a non-
constant factor over $\mathbb{Z}$ with degree less than $\deg\phi(x)$. Suppose
to the contrary that $\frac{F(x)}{c}=g(x)h(x)$ where $g(x),h(x)$ are non-
constant polynomials having coefficients in $\mathbb{Z}$ with $\deg
g(x)<\deg\phi(x)$. Fix a prime $q$ dividing $n$. Then as pointed out above,
the leading coefficient of $F(x)$ and hence those of $g(x)$ and $h(x)$ are
coprime with $q$. Therefore on passing to
$\mbox{$\mathbb{Z}$}/q\mbox{$\mathbb{Z}$}$, we see that the degree of
$\bar{g}(x)$ is same as that of $g(x)$. Thus $\deg\bar{g}(x)$ is positive and
less than $\deg\phi(x)$. This is impossible because $\bar{g}(x)$ is a divisor
of $\frac{\bar{F}(x)}{\bar{c}}=\frac{\bar{a_{n}}}{\bar{c}}\bar{\phi}(x)^{n}$
and $\bar{\phi}(x)$ is irreducible over
$\mbox{$\mathbb{Z}$}/q\mbox{$\mathbb{Z}$}$. This contradiction proves the
assertion of Step 1.

Step II. In this step, we show that $F(x)$ does not have a factor in
$\mbox{$\mathbb{Z}$}[x]$ with degree lying in the interval
$[k\deg\phi,(k+1)\deg\phi)$ for any integer $k$, $1\leq k\leq\frac{n}{2}$.
This along with Step I will prove the irreducibility of $F(x)$ over
$\mathbb{Q}$.

Suppose to the contrary that $F(x)$ has a factor in $\mbox{$\mathbb{Z}$}[x]$
with degree lying in the interval $[k\deg\phi,(k+1)\deg\phi)$ for some
positive integer $k\leq\frac{n}{2}$. By Theorem 4.A, there exists a prime
$p\geq k+1$ dividing one of the members of the set

$\\{n-k+1,\dots,n-1,n\\}.$

Choose $\ell\in\\{0,1,\dots,k-1\\}$ such that $p|(n-\ell)$; such an integer
$\ell$ is unique.

Define a monic polynomial $f(x)\in\mbox{$\mathbb{Z}$}[x]$ by
$f(x)=\sum\limits_{i=0}^{n}c_{i}\phi(x)^{i}$ where $c_{i}=\frac{n!}{i!}$. Note
that $p|c_{i}$ for $i\in\\{0,1,\dots,n-\ell-1\\}$. Let $\lambda$ denote the
slope of the right-most edge of the $\phi$-Newton polygon of $f(x)$ with
respect to $p$, where the prime $p$ is as above. We shall show that

$\displaystyle\lambda<\frac{1}{k}.$ (5)

Note that the prime $p$, being a divisor of $n!$, does not divide the content
of $a_{0}(x)$ nor does it divide $a_{n}$ by virtue of the hypothesis.
Therefore (5) together with Theorem 1.3 immediately yields that $F(x)$ does
not have a factor in $\mbox{$\mathbb{Z}$}[x]$ with degree lying in the
interval $[(\ell+1)\deg\phi,(k+1)\deg\phi)$. Since $\ell+1\leq k$, this will
imply that $F(x)$ does not have a factor in $\mbox{$\mathbb{Z}$}[x]$ with
degree in the interval $[k\deg\phi,(k+1)\deg\phi)$, which leads to a
contradiction. This contradiction will complete the proof of the theorem.

It only remains to verify (5). By definition

$\displaystyle\lambda=\max\limits_{1\leq i\leq
n}\left\\{\frac{v_{p}(c_{0})-v_{p}(c_{i})}{i}\right\\}.$

Note that $v_{p}(c_{0})-v_{p}(c_{i})=v_{p}(n!)-v_{p}(n!/i!)=v_{p}(i!)$. It is
immediate from Lemma 2.1 that

$\displaystyle\frac{v_{p}(i!)}{i}<\frac{1}{p-1};$

consequently $\lambda<\frac{1}{p-1}$ which in view of the fact $p\geq k+1$
proves (5). ∎

######  Example 4.1.

Consider $\phi(x)=x^{3}-x^{2}+1$ or $x^{4}-x-1$. It can be easily checked that
$\phi(x)$ is irreducible modulo 2 as well as modulo 3. Let $a_{0}(x)$,
$a_{1}(x)$, $a_{2}(x)$, $a_{3}(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$ be
polynomials each having degree less than $\deg\phi$ and $a_{4}$ be an integer.
Assume that $a_{4}$ and the content of $a_{0}(x)$ are coprime with $6$. By
Theorem 1.1 the polynomial

$\displaystyle
F(x)=a_{4}\phi(x)^{4}+\frac{4!}{3!}a_{3}(x)\phi(x)^{3}+\frac{4!}{2!}a_{2}(x)\phi(x)^{2}+4!a_{1}(x)\phi(x)+4!a_{0}(x)$

is irreducible over $\mathbb{Q}$.

######  Example 4.2.

Consider $\phi(x)=x^{4}-x-1$. It can be easily checked that $\phi(x)$ is
irreducible modulo 2, 3 and 5. Let $a_{i}(x)\in\mbox{$\mathbb{Z}$}[x]$ be
polynomials each having degree less than 4 for $0\leq i\leq 5$ and $a_{6}$ be
an integer. Assume that $a_{6}$ and the content of $a_{0}(x)$ are coprime with
$30$. By Theorem 1.1, the polynomial

$\displaystyle
F(x)=a_{6}\phi(x)^{6}+\frac{6!}{5!}a_{5}(x)\phi(x)^{5}+\frac{6!}{4!}a_{4}(x)\phi(x)^{4}+$
$\displaystyle\frac{6!}{3!}a_{3}(x)\phi(x)^{3}+\frac{6!}{2!}a_{2}(x)\phi(x)^{2}$
$\displaystyle\hskip 14.22636pt+6!a_{1}(x)\phi(x)+6!a_{0}(x)$

is irreducible over $\mathbb{Q}$.

It may be pointed out that irreducibility of the polynomials $F(x)$ in
Examples 2.2 and 4.2 does not seem to follow from any known irreducibility
criterion (cf. [1], [7], [8], [9], [10], [11] and [12]) including Generalized
Schönemann Irreducibility Criterion555This criterion for polynomials over
$\mathbb{Z}$ asserts that if $\phi(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$
is a monic polynomial which is irreducible modulo a fixed prime $p$ and if the
$\phi$-expansion of a polynomial $f(x)$ belonging to $\mbox{$\mathbb{Z}$}[x]$
given by $\sum\limits_{i=0}^{n}f_{i}(x)\phi(x)^{i}$ satisfies (i)
$f_{0}(x)\neq 0$, $f_{n}(x)=1$, (ii)
$\frac{v_{p}^{x}(f_{i}(x))}{n-i}\geq\frac{v_{p}^{x}(f_{0}(x))}{n}>0$ for
$0\leq i\leq n-1$, (iii) $v_{p}^{x}(f_{0}(x)),n$ are coprime, then $f(x)$ is
irreducible over $\mathbb{Q}$. .

Acknowledgements
The first author is grateful to the National Board for Higher Mathematics,
Department of Atomic Energy, India for postdoctoral fellowship. This work was
mainly done during the visit of the first author to IISER Mohali in February
2023. She is thankful to IISER Mohali for the facilities provided to her. The
second author is thankful to the Indian National Science Academy for Honorary
Scientistship.

## References

  * [1] R. Brown, Roots of generalized Schönemann polynomials in henselian extension fields, Indian J. Pure Appl. Math., 39 (2008) 403-410.
  * [2] D. M. Burton, Elementary Number Theory, McGraw-Hill, Seventh edition, 2017.
  * [3] R. F. Coleman, On the Galois groups of the exponential Taylor polynomials, L’Enseignement Math., 33 (1987) 183-189.
  * [4] G. Dumas, Sur quelques cas d’irréductibilité des polynômes á coefficients rationnels, Journal de Math. Pure et Appl., 2 (1906) 191-258.
  * [5] P. Erdős, A theorem of Sylvester and Schur, J. London Math. Soc., 9 (1934) 282-288.
  * [6] M. Filaseta, The irreducibility of all but finitely many Bessel polynomials, Acta Math., 174 (1995) 383-397.
  * [7] A. Jakhar, On the factors of a polynomial, Bull. London. Math. Soc., 52 (2020) 158-160.
  * [8] A. Jakhar, On the irreducible factors of a polynomial, Proc. Amer. Math. Soc., 148 (2020) 1429-1437.
  * [9] A. Jakhar, On the irreducible factors of a polynomial II, J. Algebra, 556 (2020) 649–655.
  * [10] A. Jakhar, A simple generalization of the Schönemann-Eisenstein irreducibility criterion, Arch. Math., 117 (2021) 375–378.
  * [11] A. Jakhar and N. Sangwan, On a mild generalization of the Schönemann irreducibility criterion, Comm. Algebra, 45 (2017) 1757-1759.
  * [12] B. Jhorar and S. K. Khanduja, A Generalization of the Eisenstein-Dumas-Schönemann Irreducibility Criterion, Proc. Edinburgh Math. Soc., 60 (2017) 937-945.
  * [13] S. K. Khanduja and S. Kumar, On prolongations of valuations via Newton polygons and liftings of polynomials, J. Pure Appl. Algebra, 216 (2012) 2648-2656.
  * [14] $\emptyset$. Ore, Newtonsche Polygone in der Theorie der algebraischen Körper, Mathematische Annalen, 99 (1928) 84-117
  * [15] P. Ribenboim, The Theory of Classical Valuations, Springer-Verlag, New York, 1999.
  * [16] T. Schönemann, Von denjenigen moduln, welche potenzen von primzahlen sind, J. Reine Angew. Math., 32 (1846) 93-105.
  * [17] I. Schur, Einige Sätze über Primzahlen mit Anwendungen auf Irreduzibilit$\ddot{a}$tsfragen, I, Sitzungsber. Preuss. Akad. Wiss. Berlin Phys.-Math. Kl., 14 (1929) 125-136.
  * [18] I. Schur, Einige Sätze über Primzahlen mit Anwendungen auf Irreduzibilitätsfragen,, Sitzungsber. Preuss. Akad. Wiss. Berlin Phys.-Math. Kl., 23 (1929) 1-24.
  * [19] J. J. Sylvester, On arithmetical series, Messenger of Math., 21 (1892) 1-19, 87-120.


# BEAM-BEAM EFFECTS UNDER THE INFLUENCE OF EXTERNAL NOISE

K. Ohmi  KEK  Oho  Tsukuba  305-0801  Japan

###### Abstract

Fast external noise, which gives fluctuation into the beam orbit, is discussed
in connection with beam–beam effects. Phase noise from crab cavities and
detection devices (position monitor) and kicker noise from the bunch by bunch
feedback system are the sources. Beam-beam collisions with fast orbit
fluctuations with turn by turn or multi-turn correlations, cause emittance
growth and luminosity degradation. We discuss the tolerance of the noise
amplitude for LHC and HL-LHC.

##  1 INTRODUCTION

Beam-beam effects under external noise are studied with the weak-strong model
in this paper. The strong beam is regarded as a target with a Gaussian charge
distribution. In the model, an external noise is introduced into the
transverse position of the strong beam at the collision point.

We first discuss an orbit (transverse position) shift of the strong beam given
as:

$\Delta x_{i+1}=(1-1/\tau)\Delta x_{i}+\delta x\cdot\hat{r}$ (1)

where $\Delta x_{i}$ is the orbit shift at the ith turn. $\tau$, $\delta x$
and $\hat{r}$ are the damping times, a constant characterizing the random
fluctuation amplitude and a Gaussian random number with unit standard
deviation. This is known as the Ornstein-Uhlenbeck process. This type of noise
is referred to as first type later.

All particles in the weak beam experience the fluctuation of the strong beam,
thus a transverse collective motion is induced. The collective motion results
in an emittance growth due to filamentation caused by the nonlinear beam–beam
force.

The stable amplitude of the fluctuation of the strong beam is given by:

$\Delta x^{2}=\langle\Delta x^{2}_{n\rightarrow\infty}\rangle=\frac{\tau\delta
x^{2}}{2}.$ (2)

The correlation function between ith and in-th turns is expressed by the
damping time as:

$\langle\Delta x_{\ell}\Delta x_{\ell+n}\rangle=\Delta x^{2}e^{-|n|/\tau}.$
(3)

The damping time is regarded as the correlation time of the fluctuation. For
white noise, which corresponds to $\tau=1$, the correlation function is
expressed as:

$\langle\Delta x_{\ell}\Delta x_{\ell+n}\rangle=\Delta x^{2}\delta_{n0}$ (4)

where $\delta_{n0}$ is the Kronecker delta.

The beam oscillates with the betatron frequency. We consider a second type of
noise as:

$\Delta x_{i+1}=(1-1/\tau)(\Delta x_{i}\cos\mu_{o}+\Delta
p_{i}\sin\mu_{o})+\delta x\hat{r}$ $\Delta p_{i+1}=(1-1/\tau)(-\Delta
x_{i}\sin\mu_{o}+\Delta p_{i}\cos\mu_{o})+\delta x\hat{r}$ (5)

where $x$ and $p$ are the coordinate and canonical momentum normalized by the
beta function, so that $J=(x^{2}+p^{2})/2$. $\mu_{o}=2\pi\nu_{o}$ is the
betatron tune multiplied by $2\pi$. In collision the offset causes an
emittance growth. The stable dipole oscillation amplitude is expressed by the
same equation as Eq. (2). The correlation function contains the betatron tune
as:

$\langle\Delta x_{\ell}\Delta x_{\ell+n}\rangle=\Delta x^{2}e^{-|n|/\tau}\cos
n\mu_{o}.$ (6)

We discuss the effect of noise for the cases of the LHC and High Luminosity-
LHC. The parameters are listed in Table 1. The phenomena depend on the
beam–beam parameter, the noise amplitude normalized by the beam size and the
Piwinski angle.

Table 1: Parameters for LHC (50 ns bunch spacing) and HL-LHC (25 ns bunch spacing). | LHC | HL-LHC(25ns) | HL-LHC(50ns)
---|---|---|---
Circumference (m) | 26 658
Energy (TeV) | 7
Tunes $Q_{x},Q_{y},Q_{s}$ | 64.31, 59.32, 0.0019
Normalized Emittance ($\mu$m) | 2.0 | 2.5 | 3.0
$\beta^{*}$ (m) | 0.55 | 0.15 | 0.15
Bunch length (m) | 0.0755
Bunch population ($10^{11}$) | $1.65$ | $2.2$ | $3.5$
Number of bunches | 1380 | 2808 | 1404
Beam-beam parameter/IP | 0.0034 | 0.005-0.011 | 0.005-0.014

##  2 Emittance growth due to the external noise

The emittance growth under an external noise and the non-linear force of the
beam–beam interaction is discussed in [1, 2, 3]. Previous work is reviewed in
this section.

The beam–beam potential for a bunch population $N_{p}$ and the transverse size
$\sigma_{r}$ is expressed as:

$U(x)=\frac{N_{p}r_{p}}{\gamma_{p}}\int_{0}^{\infty}\frac{1-e^{-x^{2}/(2\sigma_{r}^{2}+q)}}{2\sigma_{r}^{2}+q}dq$
(7)

where $r_{p}$ and $\gamma_{p}$ are the classical radius of the proton and the
relativistic factor of the (weak) beam, respectively. The potential is
expanded as a Fourier series:

$U(x)=\frac{N_{p}r_{P}}{\gamma_{p}}\sum^{\infty}_{k=0}U_{k}(a)\cos 2k\psi$ (8)

where

$U_{k}(a)=\int_{0}^{a}\left[\delta_{0k}-(2-\delta_{0k})(-1)^{k}e^{-w}I_{k}(w)\right]\frac{dw}{w},$
(9)

and $a=\beta^{*}J/2\sigma_{r}^{2}=J/2\varepsilon$. The change of $J$ per
revolution is given by the derivative of the beam–beam potential with respect
to $\psi$ as:

$\Delta J=-\frac{\partial
U}{\partial\psi}=\frac{Nr_{p}}{\gamma}\sum_{k=0}^{\infty}2kU_{k}\sin 2k\psi.$
(10)

This change, which indicates a stable sinusoidal modulation of the betatron
amplitude, does not induce emittance growth.

We consider the case in which the strong beam has a small offset ($\Delta x$).
The beam–beam potential with the offset is expanded for $\Delta x$:

$U(x+\Delta x)=U(x)+U^{\prime}(x)\Delta x.$ (11)

Here $\Delta x$ is a random variable fluctuating described by Eq. (1) or (5).

The potential with the offset is expanded as a Fourier series:

$\displaystyle U^{\prime}(J,\psi)=\frac{\partial U}{\partial J}\frac{\partial
J}{\partial x}+\frac{\partial U}{\partial\psi}\frac{\partial\psi}{\partial x}$
(13)
$\displaystyle=\frac{N_{p}r_{p}}{2\gamma\sigma_{r}}\sum_{k=0}^{\infty}G_{k}(a)\cos(2k+1)\psi.$

The Fourier coefficients as a function of $a$ are expressed as:

$G_{k}(a)=\sqrt{a}\left[U^{\prime}_{k+1}+U^{\prime}_{k}\right]+\frac{1}{\sqrt{a}}\left[(k+1)U_{k+1}-kU_{k}\right].$
(14)

where $U^{\prime}_{k}$ is the derivative with respect to $a$.

The diffusion of $J^{2}$ after N revolutions is given by:

$\langle\Delta
J^{2}(N)\rangle=\sum_{\ell=1}^{N}\sum_{n=-\ell+1}^{N-\ell}\frac{\partial
U^{\prime}(\ell)}{\partial\psi}\frac{\partial
U^{\prime}(\ell+n)}{\partial\psi}\langle\Delta x_{\ell}\Delta
x_{\ell+n}\rangle$ (15)

For turn-by-turn white noise, the correlation function is replaced by the
Kronecker delta, $\delta_{n0}$. The diffusion of $J$ is expressed by:

$\langle\Delta J^{2}\rangle=\frac{\langle\Delta
J^{2}(N)\rangle}{N}\approx\frac{N_{p}^{2}r_{p}^{2}}{8\gamma^{2}\sigma_{r}^{2}}\sum_{k=0}^{\infty}(2k+1)^{2}G_{k}(a)^{2}.$
(16)

The diffusion of $J$ per revolution is given for the fluctuation in Eq. (1)
by:

$\displaystyle\langle\Delta
J^{2}\rangle\approx\frac{N_{p}^{2}r_{p}^{2}}{8\gamma^{2}\sigma_{r}^{2}}$ (17)
$\displaystyle\sum_{n=-\infty}^{\infty}\sum_{k=0}^{\infty}(2k+1)^{2}G_{k}^{2}\cos[(2k+1)n\mu_{o}]e^{-|n|/\tau}$
$\displaystyle\approx\frac{N_{p}^{2}r_{p}^{2}}{8\gamma^{2}\sigma_{r}^{2}}\sum_{k=0}^{\infty}\frac{(2k+1)^{2}G_{k}(a)^{2}\sinh
1/\tau}{\cosh 1/\tau-\cos(2k+1)\mu_{o}}.$

The diffusion of $J$ for the second type of noise (Eq. (5)) is given using the
correlation of Eq. (6):

$\displaystyle\langle\Delta
J^{2}\rangle\approx\frac{N_{p}^{2}r_{p}^{2}}{16\gamma^{2}\sigma_{r}^{2}}\sum_{k=0}^{\infty}(2k+1)^{2}G_{k}(a)^{2}\sinh
1/\tau$ $\displaystyle\left[\frac{1}{\cosh
1/\tau-\cos(2k\mu-\delta\mu)}+\right.$ $\displaystyle\left.\frac{1}{\cosh
1/\tau-\cos(2(k+1)\mu+\delta\mu)}\right]$

where $\delta\mu$ is the tune difference between the weak and strong beam
oscillations ($\delta\mu=\mu-\mu_{o}$).

Figure 1 shows the diffusion rate of $J$ as a function of $J$. The diffusion
rate is proportional to the square of the fluctuation amplitude $\Delta x$ and
the square of the beam–beam parameter $\propto N_{p}$. The rate is normalized
by the factor $C=(N_{p}r_{p}\Delta x/\gamma\sigma_{r})^{2}/8$ in the figure.

Figure 1:  Diffusion rate given by Eq. (17). The rate is normalized by
$C=(N_{p}r_{p}\Delta x/\gamma\sigma_{r})^{2}/8$.

The emittance growth is evaluated from the diffusion rate, when the rate
$\langle\Delta J^{2}\rangle$ is proportional to $J$.

$\frac{\Delta\varepsilon}{\varepsilon}=\frac{\langle\Delta
J^{2}\rangle}{2\varepsilon J}=\frac{1}{4\varepsilon^{2}}\frac{d\langle\Delta
J^{2}\rangle}{da}.$ (19)

Figure 1 shows that the rate is proportional to $J$ for small
$J/2\varepsilon<2$. The slope of $\langle\Delta J^{2}\rangle$ for turn-by-turn
noise ($\tau~{}=~{}1$) is:

$\frac{\langle\Delta
J^{2}\rangle}{a}=\frac{N_{p}^{2}r_{p}^{2}}{8\gamma^{2}}\frac{\Delta
x^{2}}{\sigma_{r}^{2}}\times 4.4.$ (20)

The luminosity degradation rate per collision is estimated by the emittance
growth rate as:

$\Delta L/L=\left(\xi\frac{\Delta x}{\sigma_{r}}\right)^{2}\times 21.7.$ (21)

For two IPs, the formula is corrected by a factor two, i.e. $21.7\rightarrow
10.8$ and $\xi\rightarrow\xi_{tot}$, The tolerance for the noise amplitude is
given for a luminosity life time $\Delta L/L=10^{-9}$:

$\xi_{tot}\frac{\Delta x}{\sigma_{r}}=9.8\times 10^{-6}.$ (22)

We now discuss the second type of noise given by Eq. (5). Figure 2 shows the
diffusion rates. Figures 2 (a) and (b) are given for the beam-orbit
oscillation with the same tune ($\delta\mu=0$) and a difference of
$\delta\mu=\xi=0.01$, respectively. A strong enhancement of the diffusion is
seen at small amplitudes at a large correlation time in shown in Fig. 2 (a).
This behavior mainly comes from a contribution at $k=0$.

$\langle\Delta
J^{2}\rangle\approx\frac{N_{p}^{2}r_{p}^{2}}{16\gamma^{2}\sigma_{r}^{2}}G_{0}(a)^{2}\tau$
(23)

The strong beam modulation with the same tune gives an external force
oscillation to the weak beam particles. For colliding beams, the assumption,
that beam-orbit oscillations have the same tune, is not obvious. The diffusion
rate for $\delta\mu=\xi$ in Fig. 2 (b) may be better to represent the
beam–beam system. The diffusion rate, which is saturated at
$J/2\varepsilon=1$, is similar as that of $\tau=1$ on the whole. Therefore we
study the diffusion rate for $\tau=1$ in simulations.

Figure 2:  Diffusion rate given by Eq. (2). The rate is normalized by
$C=(N_{p}r_{p}\Delta x/\gamma\sigma_{r})^{2}/8$.

It may be better that the noise effects are studied in the framework of a
strong-strong model, especially for the second type of noise. The noise
induces either coherent $\sigma$ or $\pi$ modes or a continuous frequency
spectrum. The $\sigma$ mode does not contribute the emittance growth.
Emittance growth based on the strong-strong model had been discussed in [3].
The author discussed that 18% of the dipole motion induced by offset collision
into the mode with continuous frequency spectrum. The emittance increases by
smearing the dipole motion. The growth rate is expressed by:

$\frac{\delta\varepsilon}{\varepsilon}\approx\frac{K}{\left(1+\frac{1}{2\pi\tau|\xi|}\right)^{2}}\frac{\delta
x^{2}}{\sigma_{x}^{2}}=\frac{2K}{\tau\left(1+\frac{1}{2\pi\tau|\xi|}\right)^{2}}\frac{\Delta
x^{2}}{\sigma_{x}^{2}}$ (24)

where $K=0.089$ is a form factor for the emittance change induced by a dipole
amplitude, and the damping rate $1/\tau$ of the coherent motion. The emittance
growth rate is independent of the beam–beam tune shift, when $1/\tau\ll
2\pi|\xi|$, while the rate is proportional to the square of the beam–beam tune
shift, when $1/\tau\gg 2\pi|\xi|$.

Figure 3 shows the emittance growth given by Eq. (24) and by a strong-strong
beam–beam simulation [4], where the beam–beam tune shift is $\xi=0.0034/IP$.
The results agree fairly well. The strong-strong simulation suffers from
numerical noise related to the statistics of macro-particles. One million
macro-particles are used in the simulation, thus 0.1% of the offset noise is
induced by the statistics.

Figure 3:  Emittance growth given by Eq. (24) and by a strong-strong beam–beam
simulation [4].

##  3 Simulation of external noise

###  3.1 Study based on LHC

The analytical theory is based on the near solvable system far from
resonances. There is no such limitation in beam–beam simulations, while
simulations take considerable computing time to evaluate a slow emittance
growth. Simulations considering external noise are straightforward: a
modulation is applied to the strong beam with Eq. (1) or (5). Effects of
resonances, longitudinal motion and a crossing angle are taken into account in
simulations.

We only discuss weak-strong simulations taking into account external noise.
The weak beam is represented by 131072 macro-particles. The particles are
tracked one million turns interacting a strong beam located at two interaction
points. The luminosity is calculated turn-by-turn, and averaged every 100
turns. Luminosity degradation is evaluated by fitting its evolution.

Figure 4 shows the luminosity degradation for collisions without a crossing
angle. The degradation is plotted as a function of the fluctuation amplitude
for three total beam–beam parameters, $\xi_{tot}=0.02$, 0.04 and 0.05. Three
lines given by the analytical formula Eq. (22) are shown in the figure. The
simulation results agree with the formula fairly well.

Figure 4:  Diffusion rate given by a weak-strong simulation using Eq. (22).

The luminosity degradation for collision with a crossing angle ($\phi_{c}=290$
$\mu$rad) is shown in Fig. 5. The Piwinski angle is
$\phi_{c}\sigma_{z}/2\sigma_{r}=0.89$.

Figure 5:  Diffusion rate for crossing collision given by weak-strong
simulation.

Figure 6 shows luminosity degradation as a function of the beam–beam parameter
under offset noise. The tune shift is reduced to 70% for the crossing
collision. The luminosity degradation for noise is independent of the crossing
angle. At high beam–beam parameters $>0.05$, the luminosity degradation due to
the crossing angle is dominant.

Figure 6:  Luminosity degradation as a function of the beam–beam parameter
with offset noise.

There was no qualitative change for collision without a crossing angle. For
$\xi_{tot}=0.035$, a degradation due to the crossing angle is seen, but a
significant cross-talk is not observed. The degradation of the luminosity due
to the fluctuation depends on $\xi_{tot}$, but hardly on the presence of the
crossing angle.

###  3.2 High Luminosity LHC (HL-LHC)

For the HL-LHC, a higher luminosity is the target and obtained by increasing
the bunch population and squeezing to smaller beta function, while the pile up
of collision events sets an upper limit for the luminosity at
$L/coll=2.6\times 10^{31}$ cm-2s-1. The luminosity at $\beta=0.15$ m is
expected to be $L/coll=8.6$ or $18\times 10^{31}$ cm-2s-1 for the bunch
population of 2.2 or 3.5$\times 10^{11}$, respectively. Therefore luminosity
levelling keeping the luminosity constant at $L/coll=2.6\times 10^{31}$
cm-2s-1 is proposed. The levelling can be done by controlling the crab cavity
voltage or the beta function at the interaction point (IP). Leveling with the
beta function, the total beam–beam parameter (2IP) is $0.011\times 2=0.022$
(25ns) or $0.014\times 2=0.028$ (50ns) at the early stage of the collision,
where the beta function is 0.49 m or 1.02 m. The results given in the previous
subsection are applied for the parameters:

$\frac{\Delta x}{\sigma_{r}}=4.5\times 10^{-4}~{}\mbox{ or }~{}3.5\times
10^{-4}.$ (25)

for 25 ns and 50 ns, respectively.

Using a levelling with crab cavities, the crab voltage increases to keep the
luminosity constant while the beam current decreases. At the early stage of
collision, the crab voltage is low and two beams collide with a large Piwinski
angle, where $\phi_{c}\sigma_{z}/2\sigma_{r}=3.14$ or 2.87 for 25 ns or 50 ns,
respectively. We study the effects of noise for collisions with a large
Piwinski angle.

Figure 7 shows the luminosity degradation rate as a function of the offset
amplitude. The simulation is performed for two IPs. The tune shift is 0.0015
or 0.0050 in the crossing or orthogonal plane for the design bunch population
of $N_{p}=2.2\times 10^{11}$ (25ns). The tune shift is 0.0065 in both planes,
due to the combination of the horizontal and vertical crossing. The
fluctuation amplitude 0.2% is a tolerable limit for $\Delta L/L_{0}=10^{-9}$
as shown in the figure. The simple formula Eq. (22) is satisfied for the HL-
LHC, $0.0065\times 0.002=1.3\times 10^{-5}$, with 30% difference from the
formula.

Figure 7:  Luminosity degradation as a function of noise amplitude.

Figure 8 shows the luminosity degradation as a function of the correlation
time. The luminosity degradation, which scales as $1/\tau$, is consistent with
Eq. (17).

Figure 8: Luminosity degradation as a function of the correlation time.

The effect of noise on its amplitude and the dependence on the correlation
time is similar for collisions with a large Piwinski angle and those for
ordinary collisions without a crossing angle. The luminosity degradation
depends on the beam–beam parameter and the noise amplitude, but with little
dependence on the Piwinski angle.

###  3.3 Tolerance for crab cavity phase noise in HL-LHC

Crab cavities are used to compensate the crossing angle ($\phi_{c}=590$
$\mu$rad) at IP. The relation of the phase noise and collision offset is given
by

$\Delta\varphi_{cc}=\frac{\omega_{cc}}{c\phi_{c}/2}\Delta x,$ (26)

where $\Delta\varphi_{cc}$ and $\omega_{cc}$ are the phase fluctuation and
frequency of the crab cavity.

Using beta function levelling the beam–beam parameter is very high,
$\xi_{tot}=0.022$ or 0.028 for 25 ns or 50 ns, respectively. The tolerance of
the noise amplitude is given by Eq. (25). The corresponding phase error is
$\Delta\varphi=1.6\times 10^{-4}$ or $2.3\times 10^{-4}$ rad.

For the crab cavity levelling, the beam–beam parameter is $\xi_{tot}=0.0065$.
The tolerance of the noise amplitude is $\Delta x/\sigma_{r}=0.002$ and the
corresponding phase error is $\Delta\varphi=4\times 10^{-3}$, where the crab
angle is 10% of the crossing angle $\phi_{c}=59$ $\mu$rad ($L/coll=2.7\times
10^{31}$ cm-2s-1).

The crab cavity noise was measured at KEKB, $1.7\times~{}10^{-4}$ rad for
frequencies above 1 kHz ($\tau<10$). The value is critical for beta function
levelling, because of the large beam–beam parameter. Using four crab cavities,
the noise tolerance is twice as large, while for the crab voltage levelling,
the measured phase error is tolerable.

###  3.4 Incoherent noise due to intra-beam scattering

Emittance growth times due to intra-beam scattering (IBS) are 105 h and 63 h
for the horizontal and longitudinal planes, respectively, in the nominal LHC
[13]. The transverse emittance and bunch population in the nominal are
$5.0\times 10^{-10}$ and $1.15\times 10^{11}$, respectively. The horizontal
IBS growth rate is approximately proportional to the particle density in the
six dimensional phase space. The growth time is 40 h for $\xi_{tot}$ = 0.02 in
this paper ($\varepsilon=2.7\times 10^{-10}$ and $Np=1.63\times 10^{11}$). The
fluctuation is $\delta x/\sigma_{x}=5.5\times 10^{-5}$ for $\xi_{tot}$ = 0.05.
The luminosity degradation is determined by geometrical emittance growth
$\delta L/L_{0}=\delta x^{2}/\sigma_{x}^{2}$ for incoherent noise.

##  4 Coherent beam–beam effects under external noise

Effects of external noise in crab cavity were performed in KEKB during 2008
and 2009 [6]. Sinusoidal noise is applied to the crab cavity RF system. Near
the $\sigma$ mode tune, a strong luminosity drop of 80 % was seen when
suddenly exceeding a threshold excitation amplitude. A smaller luminosity drop
($L=0.9L_{0}$) was seen near $\pi$ mode frequency. Strong-strong simulations
reproduced these luminosity drops. A systematic study using the strong-strong
simulation showed that these characteristic phenomena for coherent nonlinear
beam–beam interactions. A similar phenomenon was observed in Ref.[7]. The
detailed analysis is published in Ref.[6].

##  5 CONCLUSIONS

Fast noise of the collision offset degrades the luminosity performance in
hadron colliders. The luminosity degradation depends on the product of the
noise amplitude and the beam–beam parameter as shown in Eq. (22), with little
dependence on the crossing angle. A tolerance for the crab cavity phase error
was obtained for the HL-LHC.

The crab cavity noise was measured at KEKB, $1.7\times 10^{-4}$ rad above 1
kHz ($\tau<10$). The value is critical for beta function levelling, because of
the high beam–beam parameter. For the crab voltage levelling, the measured
phase error is tolerable because of the small beam–beam tune shift.

Further studies related to beam–beam modes should be done using strong-strong
models.

##  6 Acknowledgement

The author would like to thank Drs. Y. Alexahin, R. Calaga, W. Herr, S. Paret,
T. Pieloni, J. Qiang, R. Tomas and F. Zimmermann for fruitful discussions.

## References

  * [1] G. Stupakov, SSCL-560 (1991).
  * [2] T. Sen and J.A. Ellison, Phys. Rev. Lett. 77, (1996) 1051.
  * [3] Y. Alexahin, Nucl. Instr. Methods in Phys. Res. A391, (1996) 73.
  * [4] K. Ohmi $et~{}al.$, “Beam-beam effect with an external noise in LHC”, Proc. Particle Accelerator Conference, Albuquerque, June, 2007, (2007) p. 1496; http://www.JACoW.org. K. Ohmi $et~{}al.$, CERN-LHC-PROJECT-REPORT-1048.
  * [5] LHC Design Report, CERN-2004-003-V-1; http://lhc.web.cern.ch/lhc/lhc-designreport.html.
  * [6] K. Ohmi $et~{}al.$,Phys. Rev. ST-AB, 14, 111003 (2011).
  * [7] T. Ieiri and K. Hirata, Proc. Particle Accelerator Conference, March, 1989, Chicago, USA.


# Mathematical analysis and numerical simulation of the
Guyer–Krumhansl heat equation

A. J. A. Ramos<EMAIL_ADDRESS>R. Kovács<EMAIL_ADDRESS>M. M. Freitas
<EMAIL_ADDRESS>D. S. Almeida Júnior<EMAIL_ADDRESS>Faculty of
Mathematics, Federal University of Pará, Raimundo Santana Street s/n,
Salinópolis–PA, 68721-000, Brazil  Department of Energy Engineering, Faculty
of Mechanical Engineering, Budapest University of Technology and Economics,
Muegyetem rkp. 3., H-1111 Budapest, Hungary  Department of Theoretical
Physics, Wigner Research Centre for Physics, Budapest, Hungary  PhD Program in
Mathematics, Federal University of Pará, Augusto Corrêa Street 01, Belém–PA,
66075-110, Brazil

###### Abstract

The Guyer-Krumhansl heat equation has numerous important practical
applications in both low-temperature and room temperature heat conduction
problems. In recent years, it turned out that the Guyer-Krumhansl model can
effectively describe the thermal behaviour of macroscale heterogeneous
materials. Thus, the Guyer-Krumhansl equation is a promising candidate to be
the next standard model in engineering. However, to support the Guyer-
Krumhansl equation’s introduction into the engineering practice, its
mathematical properties must be thoroughly investigated and understood. In the
present paper, we show the basic structure of this particular heat equation,
focusing on the differences in comparison to the Fourier heat equation
obtained when $(\tau_{q},\mu^{2})\rightarrow(0,0)$. Additionally, we prove the
well-posedness of a particular, practically significant initial and boundary
value problem. The stability of the solution is also investigated in the
discrete space using a finite difference approach.

###### keywords:

non-Fourier heat conduction , finite difference discretization , irreversible
thermodynamics

###### MSC:

[2010] 35E15 , 65M06 , 93D20

##  1 Introduction

Among the various heat conduction models, the Guyer-Krumhansl (GK) equation
has an exceptional background and possesses significant potential in
engineering practice. It is first derived on the basis of kinetic theory and
applied for low-temperature problems [1]. It served as a starting point to
find the celebrated window condition, i.e., the relation that helped the
researchers detect the second sound in solids, too [2]. The kinetic approach
is limited due to the prescribed heat conduction mechanisms with phonons [3].
Therefore, one cannot apply such a framework to more general room temperature
macroscale problems. However, there is the possibility of utilizing a
continuum approach – called internal variables – to derive the Guyer-Krumhansl
equation based on irreversible thermodynamics [4]. Note that different
continuum approaches can also be suitable such as GENERIC [5, 6] and Extended
Irreversible Thermodynamics [7, 8]. A continuum approach results in a model
with a notably wider range of validity since it allows to fit the parameters
to the observed phenomenon [9].

In this work we consider the one-dimensional version of the GK equation, i.e.,

$\displaystyle\rho cT_{t}+q_{x}$
$\displaystyle=0\quad\mbox{in}\quad(0,l)\times(0,\infty),$ (1)
$\displaystyle\tau_{q}q_{t}+q-\mu^{2}q_{xx}+kT_{x}$
$\displaystyle=0\quad\mbox{in}\quad(0,l)\times(0,\infty),$ (2)

subject to boundary conditions

$\displaystyle q(0,t)=q(l,t)=0,\quad\mbox{for all}\quad t\geq 0,$ (3)

and initial data

$\displaystyle T(x,0)=T_{0}(x),\quad q(x,0)=q_{0}(x),\quad x\in(0,l),$ (4)

wherein Eq. (1) is the internal energy balance with $\rho$ and $c$ being the
mass density and specific heat, and no source terms are considered in a rigid
conductor. Furthermore, $T$ denotes the temperature, $q$ is the heat flux, and
the indicies $t$ and $x$ denote the corresponding derivatives respect to time
and space. Eq. (2) is a constitutive equation in which $\tau_{q},\mu^{2}$ and
$k$ are positive coefficients following from the second law of thermodynamics
but no further restrictions are acting on them in the linear regime
(T-independent coefficients). Consequently, that continuum model can be
compatible with kinetic theory if it adapts the particular form for the
coefficients, if necessary [10]. Thanks to the flexibility of the
coefficients, it was possible to explain and model several room temperature
experiments successfully without any need for phonon hydrodynamics [9, 11,
12]. All these experiments were performed using a heterogeneous material such
as rocks, foams and layered samples. In that sense, the GK model is promising
and could be the next standard model in the engineerint practice in the
future. However, it would not be possible without thoroughly investigating and
understanding its mathematical properties first.

Therefore, in the following, we present an essential attribute of the
coefficients, characteristic only for non-Fourier heat equations, which is a
rarely known feature of the GK equation. Namely, the coefficients are
functionally connected, which becomes crucial for T-dependent material
parameters [13]. Then, considering only the linear problem, we prove the well-
posedness of the model (1)–(4). This is also crucial since there are popular
models in the corresponding literature (such as the famous dual-phase-lag
equation), which are not necessarily well-posed due to the missing physical
and mathematical consistency [14, 15, 16].

##  2 Temperature-dependent coefficients

In the framework of irreversible continuum thermodynamics, the constitutive
equations are derived utilizing the second law of thermodynamics.
Mathematically speaking, the constitutive equations are found as the solutions
of the entropy inequality

$\displaystyle\rho s_{t}+J_{x}=\sigma_{s}\geq 0,$ (5)

where $s$ is the entropy density, $J$ is the flux of the entropy density and
$\sigma_{s}$ stands for the entropy production. The solutions of Eq. (5) are
called Onsagerian relations. The detailed derivation procedure can be found in
[4]. In the case of the GK equation, the Onsagerian relation is

$\displaystyle\left(\frac{1}{T}+l_{1}q_{x}\right)_{x}-\rho mq_{t}-l_{2}q=0,$
(6)

with $l_{1},m\geq 0$, $l_{2}>0$, furthermore

$\displaystyle\tau_{q}=\frac{\rho m}{l_{2}},\quad
k=\frac{1}{l_{2}T^{2}},\quad\mu^{2}=\frac{l_{1}}{l_{2}}$ (7)

forming the coefficients of the GK constitutive equation, which reads

$\displaystyle\tau_{q}q_{t}+q-\mu^{2}q_{xx}+kT_{x}=0.$ (8)

Eq. (8) is valid only for constant coefficients, i.e., when $l_{1},l_{2}$ and
$m$ do not depend on the state variable T, which is also apparent from Eqs.
(6)–(7). Notably, Eq. (7) presents the functional connections between the
coefficients. In case of a temperature-dependent thermal conductivity $k(T)$,
$l_{2}$ becomes a function of $T$, too, consequently both $\tau_{q}$ and
$\mu^{2}$ become $T$-dependent. This is a non-trivial consequence of
thermodynamics and characteristic only for non-Fourier heat equations. As
$k(T)$ is usual in the engineering practice, and Eq. (7) is the direct
consequence of the second law of thermodynamics, it is an essential physical
property. That attribute is entirely missing from the dual-phase-lag approach.

Furthermore, the $T$-dependency of $\tau_{q}$ introduced by $l_{2}(T)$ is not
necessarily the same, which is observed in the experiments for $\tau_{q}(T)$.
In other words, one must adjust the $T$-dependence of $\tau_{q}$ according to
the experiments using $\rho$ or $m$. Nevertheless, it has further
consequences. If $\rho$ is $T$-dependent, then mechanics must be involved as
well [13]. If $m$ becomes a function, then the entire derivation must be
repeated since Eq. (6) is valid only for constant $m$. Either way, Eq. (6)
loses its validity, and one must be cautious. Additionally, due to $l_{2}(T)$,
$\mu^{2}$ is also $T$-dependent. If one has direct measurements about
$\mu^{2}(T)$, then $l_{1}$ can be used to adjust the proper $T$-dependent
behaviour. Following Eq. (6), $l_{1}(T)$ would lead to a new term, i.e.,

$\displaystyle\left(\frac{1}{T}\right)_{x}+\big{(}l_{1}(T)\big{)}_{x}q_{x}+l_{1}(T)q_{xx}-\rho
mq_{t}-l_{2}q=0,$ (9)

holds, where

$\displaystyle\big{(}l_{1}(T)\big{)}_{x}q_{x}=\frac{\textrm{d}l_{1}(T)}{\textrm{d}T}T_{x}q_{x},$
(10)

and

$\displaystyle\frac{1}{l_{2}}\left(-\frac{1}{T^{2}}+\frac{\textrm{d}l_{1}(T)}{\textrm{d}T}q_{x}\right)T_{x}+\frac{l_{1}(T)}{l_{2}}q_{xx}-\frac{\rho
m}{l_{2}}q_{t}-q=0.$ (11)

Hence, beyond the $T$-dependence of $\mu^{2}(T)$, $l_{1}(T)$ introduces a new
term, which could modify the thermal conductivity, leading to state dependence
in $k$, seemingly. Furthermore, the positivity of $k$ could provide further
restrictions on $l_{1}(T)$ since it appears with opposite sign in Eq. (11).
Interestingly, $k$ is not simply $T$-dependent but also begins to depend on
the spatial derivative (divergence) of $q$. It stands as an experimental
challenge to distinguish this effect from the usual $T$-dependence of thermal
conductivity, that is, by introducing $l_{2}(T)$. These have significantly
different consequences and physical implications.

Consequently, to determine the thermal conductivity of a material with complex
inner structure, it is strongly suggested to apply only steady-state
measurement methods. In steady-state, the equilibrium is not modified by the
non-Fourier terms, and that approach could provide a reliable method to
separate the $T$-dependence of these coefficients.

Overall, if the thermal conductivity depends on the temperature through
$l_{2}(T)$, then the other parameters inherit that dependence, and one must
include mechanics in the model because of $\rho(T)$. Furthermore, if $\mu^{2}$
depends on the temperature through $l_{1}(T)$, then it influences the thermal
conductivity as well but leaves the relaxation time $\tau_{q}$ unchanged.

##  3 Well-posedness

Among the numerous non-Fourier heat conduction models, the dual-phase-lag
(DPL) concept, i.e.,

$\displaystyle q(x,t+\tau_{q})=-kT_{x}(x,t+\tau_{T})$ (12)

stands out in popularity. Here, it is assumed that both quantities possess
some time delay, and various models are derived by Taylor series expansion.
Although the DPL concept is used to model various heat conduction problems,
such an approach hides multiple substantial physical and mathematical aspects.
The Eq. (12) has no immediate compatibility with the first and second laws of
thermodynamics, and the DPL model needs additional conditions. In other words,
neither the energy conservation nor the stability is guaranteed, contrary to
the Guyer-Krumhansl equation. It has been shown that Eq. (12) could lead to
ill-posed problems in certain situations [14, 15, 16]. Therefore, it is of
crucial importance to have a model with a solid mathematical and physical
basis in order to avoid such stability or ill-posedness issues.

Here, in the following, we show that well-posedness is fulfilled for the GK
equation in the most crucial situation from a physical point of view. Let us
recall that in Eqs. (3) and (4) we use adiabatic (perfect thermal insulation)
boundary conditions for inhomogeneous initial conditions. The other
possibilities, the prescribed temperature and heat convection (i.e., the usual
first and third-type) boundary conditions, would implement dissipation into
the system. These could result in stable and convergent equilibrium solutions
even for unstable and physically inconsistent models. In our case, the
adiabatic boundary condition means that the initial energy content of the
entire system remains the same and, in equilibrium, becomes homogeneous in
space. This holds for the internal energy for which we have a balance
equation, and it is a conserved quantity.

However, from a mathematical point of view, energy as a notion is used
differently, usually referring to a variational background, and here we follow
the usual convention. Therefore, let us introduce the term ‘functional energy’
of the system (1)–(4) as

$\displaystyle E(t):=\frac{\rho
c}{2}\int_{0}^{l}T^{2}dx+\frac{\tau_{q}}{2k}\int_{0}^{l}q^{2}dx,$ (13)

which decreases in time as $q$ will be uniformly zero in equilibrium, i.e.,

$\frac{\textrm{d}}{\textrm{d}t}E(t)=-\frac{1}{k}\int_{0}^{l}q^{2}\textrm{d}x-\frac{\mu^{2}}{k}\int_{0}^{l}q_{x}^{2}\textrm{d}x<0,\quad\mbox{for
all}\quad t\geq 0,$ (14)

consequently, $E$ is a non-increasing function of the time variable $t$.

In this section we will show that the system (1)–(4) is well-posed using the
semigroup technique. The following notations will be used in this paper:

$\displaystyle\big{<}u,v\big{>}:=\big{<}u,v\big{>}_{L^{2}(0,l)}\quad\mbox{and}\quad\|u\|:=\sqrt{\big{<}u,u\big{>}},\quad\mbox{for
all}\quad u,v\in L^{2}(0,l).$ (15)

Let us start by introducing the phase space by

$\displaystyle\mathcal{H}:=L^{2}(0,l)\times L^{2}(0,l),$ (16)

which is a Hilbert space endowed with the norm

$\displaystyle\|U\|^{2}_{\cal H}=\rho c\|u\|^{2}+\frac{\tau_{q}}{k}\|v\|^{2}$
(17)

and the corresponding inner product $\big{<}\cdot,\cdot\big{>}_{\cal H}$, for
all $U=(u,v)^{\top}\in\mathcal{H}$. Thus, denoting $u:=T$ and $v:=q$ we can
transform the system (1)–(4) into the first-order set of equations

$\displaystyle\begin{cases}U_{t}(t)={\cal A}U(t),\quad\mbox{for all}\quad
t>0,\\\ U(0)=U_{0},\end{cases}$ (18)

where $U_{0}=(T_{0},\,q_{0})^{\top}$ and $\cal A:D(\cal A)\subset\cal
H\rightarrow\cal H$ is the differential operator given by

$\displaystyle\mathcal{A}:=\left(\begin{array}[]{ccccc}0&-\frac{1}{\rho
c}\frac{\partial(\cdot)}{\partial x}\\\ \\\
-\frac{k}{\tau_{q}}\frac{\partial(\cdot)}{\partial
x}&-\frac{1}{\tau_{q}}I_{d}(\cdot)+\frac{\mu^{2}}{\tau_{q}}\frac{\partial^{2}(\cdot)}{\partial
x^{2}}\end{array}\right),$ (22)

where $I_{d}(\cdot)$ denotes the identity operator, furthermore

$\displaystyle D(\mathcal{A}):=\Big{\\{}U=(u,v)^{\top}\in\mathcal{H};\ \
\mu^{2}v_{x}-ku\in H^{1}(0,l),\ v\in H^{1}_{0}(0,l)\Big{\\}}.$

######  Theorem 3.1.

The system (18) is well-posed, i.e., for any $U_{0}\in\mathcal{H}$, the system
(18) has a unique weak solution $U(t)=e^{\mathcal{A}t}U_{0}\in
C\big{(}[0,\infty\big{)};\,\mathcal{H})$. Furthermore, if $U_{0}\in
D(\mathcal{A})$, $U(t)\in C^{1}\big{(}[0,\infty);\,\mathcal{H}\big{)}\cap
C\big{(}[0,\infty);\,D(\mathcal{A})\big{)}$ becomes the classic solution for
(18).

Proof:  It is not difficult to see that $\mathcal{A}$ is a dissipative
operator in the space $\mathcal{H}$. More precisely, we have

$\displaystyle\big{\langle}\mathcal{A}U,U\big{\rangle}_{\mathcal{H}}=-\frac{1}{k}\|q\|^{2}-\frac{\mu^{2}}{k}\|q_{x}\|^{2}<0.$
(23)

Next, we prove that $\lambda I-\mathcal{A}$ is surjective, for all
$\lambda>0$. Since $\mathcal{A}$ is dissipative and $D(\mathcal{A})$ is dense
in $\mathcal{H}$, it is sufficient to show that $\mathcal{A}$ is maximal, that
is, $\lambda I-\mathcal{A}$ is surjective. For this purpose, given
$F=\big{(}f_{1},f_{2}\big{)}^{\top}\in\mathcal{H}$, we seek
$U=\big{(}u,v\big{)}^{\top}\in D(\mathcal{A})$ which is the solution of

$\displaystyle(\lambda I-\mathcal{A})U=F,$ (24)

that is, the entries of $U$ satisfy the system of equations

$\displaystyle\lambda\rho c\,u+v_{x}=\rho cf_{1}\ $ $\displaystyle\in$
$\displaystyle\ L^{2}(0,l),$ (25)
$\displaystyle\lambda\tau_{q}v+v-\mu^{2}v_{xx}+ku_{x}=\tau_{q}f_{2}\ $
$\displaystyle\in$ $\displaystyle\ L^{2}(0,l).$ (26)

Solving the system (25)–(26) is equivalent to finding

$(u,v)\in L^{2}(0,l)\times H^{1}_{0}(0,l),$

such that

$\displaystyle\lambda\rho
c\,\big{<}u,\,\widetilde{u}\big{>}+\big{<}v_{x},\,\widetilde{u}\big{>}$
$\displaystyle=$ $\displaystyle\rho c\,\big{<}f_{1},\,\widetilde{u}\big{>},$
(27)
$\displaystyle\frac{\lambda\tau_{q}}{k}\big{<}v,\,\widetilde{v}\big{>}+\frac{1}{k}\big{<}v,\,\widetilde{v}\big{>}+\frac{1}{k}\big{<}\mu^{2}v_{x}-ku,\,\widetilde{v}_{x}\big{>}$
$\displaystyle=$
$\displaystyle\frac{\tau_{q}}{k}\big{<}f_{2},\,\widetilde{v}\big{>},$ (28)

for all $(\widetilde{u},\widetilde{v})\in L^{2}(0,l)\times H^{1}_{0}(0,l)$.

Now we observe that solving the system (25)–(26) is equivalent to solve the
problem

$\mathcal{B}\big{(}(u,v),\left(\widetilde{u},\widetilde{v}\right)\big{)}=\mathcal{G}\left(\widetilde{u},\widetilde{v}\right),$
(29)

where $\mathcal{B}:\big{[}L^{2}(0,l)\times
H^{1}_{0}(0,l)\big{]}^{2}\rightarrow\mathbb{R}$ is the bilinear form given by

$\displaystyle\mathcal{B}\big{(}(u,v),\big{(}\widetilde{u},\widetilde{v}\big{)}\big{)}:=\lambda\rho
c\,\big{<}u,\,\widetilde{u}\big{>}+\frac{1}{k}(1+\lambda\tau_{q})\big{<}v,\,\widetilde{v}\big{>}+\big{<}v_{x},\,\widetilde{u}\big{>}+\frac{1}{k}\big{<}\mu^{2}v_{x}-ku,\,\widetilde{v}_{x}\big{>}$
(30)

and $\mathcal{G}:L^{2}(0,l)\times H^{1}_{0}(0,l)\rightarrow\mathbb{R}$ is the
linear form given by

$\mathcal{G}\big{(}\widetilde{u},\,\widetilde{v}\big{)}:=\rho
c\,\big{<}f_{1},\,\widetilde{u}\big{>}+\frac{\tau_{q}}{k}\big{<}f_{2},\,\widetilde{v}\big{>}.$
(31)

Now, we introduce the Hilbert space $\mathcal{V}:=L^{2}(0,l)\times
H^{1}_{0}(0,l)$ equipped with the norm

$\big{\|}(u,v)\big{\|}_{\mathcal{V}}^{2}=\|u\|^{2}+\|v\|^{2}+\mu^{2}\|v_{x}\|^{2}.$
(32)

It is clear that $\mathcal{B}$ and $\mathcal{G}$ are bounded. Furthermore, we
can obtain that there exists a positive constant $C$ such that

$\displaystyle\mathcal{B}\big{(}(u,v),\,(u,v)\big{)}$ $\displaystyle=$
$\displaystyle\lambda\rho
c\|u\|^{2}+\frac{1}{k}(1+\lambda\tau_{q})\|v\|^{2}+\frac{\mu^{2}}{k}\|v_{x}\|^{2}$
(33) $\displaystyle\geq$ $\displaystyle C\|(u,v)\|_{\mathcal{V}}^{2},$

which implies that $\mathcal{B}$ is coercive.

Hence, we assert that $\mathcal{B}$ is continuous and coercive form on
$\mathcal{V}\times\mathcal{V}$ and $\mathcal{G}$ is continuous form on
$\mathcal{V}$. So applying the Lax-Milgram Theorem, we deduce that for all

$(\widetilde{u},\widetilde{v})\in L^{2}(0,l)\times H^{1}_{0}(0,l)$ (34)

the problem (29) admits a unique solution

$(u,v)\in L^{2}(0,l)\times H^{1}_{0}(0,l).$ (35)

Furthermore, it follows from Eqs. (25) and (26) that $\mu^{2}v_{x}-ku\in
H^{1}(0,l)$.

Thus the operator $\lambda I-\mathcal{A}$ is surjective for all $\lambda>0$.
Consequently, $\mathcal{A}$ is a maximal operator. Hence, the result of
Theorem 3.1 follows from Lumer–Phillips Theorem (see [18]).

##  4 Uniform stabilization

In this section we prove that the solution of the system (1)–(4) decay
exponentially uniformly with respect to the parameters $\tau_{q}$ and
$\mu^{2}$ to an equilibrium point. In the $(\tau_{q},\mu^{2})\rightarrow(0,0)$
limit, we obtain Fourier model for which the energy decays exponentially as
well in agreement with the thermodynamic requirements. Furthermore, we prove
that the equilibrium point (or state) depends directly on the space where the
initial condition $T_{0}$ is inserted.

######  Theorem 4.1.

Let $E(t)$ be the functional energy for the system (1)–(4). Suppose that the
initial data $(T_{0},q_{0})\in D(\mathcal{A})$. Then, there exist constants
$M_{0}$, $M_{1}$, $\omega>0$ uniformly bounded as
$(\tau_{q},\mu^{2})\rightarrow(0,0)$, such that $E(t)$ satisfies the following
decay estimates:

  * $(i)$

If the initial data $T_{0}\notin L_{*}^{2}(0,l):=\big{\\{}\mathcal{T}\in
L^{2}(0,l);\ \int_{0}^{l}\mathcal{T}dx=0\big{\\}}$, then

$\displaystyle E(t)\leq M_{0}E(0)e^{-\omega
t}+M_{1}\|C_{T}\|_{\infty},\quad\mbox{for all}\quad t\geq 0,$ (36)

where $C_{T}(t):=\big{(}\mu^{2}q_{x}(0,t)-kT(0,t)\big{)}\int_{0}^{l}T(s,t)ds$
and $\|C_{T}\|_{\infty}:=\sup\big{\\{}|C_{T}(t)|;\
t\in\mathbb{R}^{+}\big{\\}}$ uniformly bounded as $\mu^{2}\rightarrow 0$.

  * $(ii)$

If the initial data $T_{0}\in L_{*}^{2}(0,l)$, then

$\displaystyle E(t)\leq M_{0}E(0)e^{-\omega t},\quad\mbox{for all}\quad t\geq
0.$ (37)

From a physical point of view, Theorem 4.1 separates two cases of initial
conditions. The first item, $(i)$, describes a physically realistic situation
in which the initial total energy content of the system (integration of the
initial temperature distribution) is not zero, therefore the temperature
history must converge to the equilibrium prescribed by that energy content for
an insulated (thermally closed) system.

The second part, $(ii)$, however, is only for the mathematical completeness,
since that situation is not physically realizable. The temperature - when
calculating the energy content - must be measured in Kelvin degrees, therefore
the zero initial energy content can be achieved only by uniformly zero
temperature distribution, otherwise one implies negative absolute temperature.
In the following, both parts are proved, but we emphasize that only part $(i)$
has practical importance.

Proof:  Integrating Eq. (1) over $[x,l]\subset[0,l]$, multiplying by
$\int_{x}^{l}T(s,t)ds$ and integrating over $[0,l]$ we obtain

$\displaystyle\rho c\
\big{<}\int_{x}^{l}T_{t}(s,t)\textrm{d}s,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}-\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}=0.$
(38)

On the other hand, multiplying Eq. $(\ref{Eq2})$ by
$\int_{x}^{l}T(s,t)\textrm{d}s$ and integrating over $[0,l]$ we have

$\displaystyle\tau_{q}\big{<}q_{t},\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}+\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}+\big{(}\mu^{2}q_{x}(0,t)-kT(0,t)\big{)}\int^{l}_{0}T(s,t)\textrm{d}s-\big{<}\mu^{2}q_{x}-kT,\,T\big{>}=0.$
(39)

Adding Eqs. (38) and (39), and considering that
$\big{<}\int_{x}^{l}T_{t}(s,t)\textrm{d}s,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}=\frac{1}{2}\frac{\textrm{d}}{\textrm{d}t}\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2}$,
we achieve

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\frac{\rho
c}{2}\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2}+\underbrace{\tau_{q}\big{<}q_{t},\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}-\mu^{2}\big{<}q_{x},\,T\big{>}}_{I:=}+\big{(}\mu^{2}q_{x}(0,t)-kT(0,t)\big{)}\int^{l}_{0}T(s,t)\textrm{d}s+k\|T\|^{2}=0.$
(40)

Exploiting the identity
$\big{<}q_{t},\,\int_{x}^{l}T(s,t)ds\big{>}=\frac{\textrm{d}}{\textrm{d}t}\big{<}q,\,\int_{x}^{l}T(s,t)ds\big{>}-\big{<}q,\,\int_{x}^{l}T_{t}(s,t)\textrm{d}s\big{>}$
together with the Eq. (1), it yields

$\displaystyle I$ $\displaystyle:=$
$\displaystyle\tau_{q}\big{<}q_{t},\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}-\mu^{2}\big{<}q_{x},\,T\big{>}$
$\displaystyle=$
$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\tau_{q}\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}-\tau_{q}\big{<}q,\,\int_{x}^{l}T_{t}(s,t)\textrm{d}s\big{>}+\mu^{2}\rho
c\big{<}T_{t},\,T\big{>}$ $\displaystyle=$
$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\bigg{(}\tau_{q}\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}+\frac{\rho
c}{2}\mu^{2}\|T\|^{2}\bigg{)}-\frac{\tau_{q}}{\rho c}\|q\|^{2}.$

Substituting $I$ into $(\ref{eq.I})$ we obtain

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}F(t)=-k\|T\|^{2}+\frac{\tau_{q}}{\rho
c}\|q\|^{2}+C_{T}(t),$ (41)

in which
$C_{T}(t):=\big{(}\mu^{2}q_{x}(0,t)-kT(0,t)\big{)}\int_{0}^{l}T(s,t)\textrm{d}s$,
and

$\displaystyle F(t):=\frac{\rho
c}{2}\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2}+\frac{\rho
c}{2}\mu^{2}\|T\|^{2}+\tau_{q}\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}.$
(42)

Now, let us consider the following Lyapunov functional,

$\displaystyle\mathcal{L}(t):=\Big{(}2l^{2}+2\mu^{2}+\frac{\tau_{q}k}{\rho
c}\Big{)}E(t)+F(t),\quad\mbox{for all}\quad t\geq 0.$ (43)

Utilizing Cauchy–Schwarz and Young inequalities, we obtain

$\displaystyle\Big{|}\mathcal{L}(t)-\Big{(}2l^{2}+2\mu^{2}+\frac{\tau_{q}k}{\rho
c}\Big{)}E(t)\Big{|}$ $\displaystyle\leq$ $\displaystyle\frac{\rho
c}{2}\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2}+\frac{\rho
c}{2}\mu^{2}\|T\|^{2}+\tau_{q}\Big{|}\big{<}q,\,\int_{x}^{l}T(s,t)\textrm{d}s\big{>}\Big{|}$
$\displaystyle\leq$ $\displaystyle\frac{\rho
c}{2}\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2}+\frac{\rho
c}{2}\mu^{2}\|T\|^{2}+\frac{\tau_{q}}{4\varepsilon}\|q\|^{2}+\tau_{q}\varepsilon\big{\|}\int_{x}^{l}T(s,t)\textrm{d}s\big{\|}^{2},$

for all $\varepsilon>0$. Choosing $\varepsilon:=k/2l^{2}$ and applying
Wirtinger’s inequality (see Corollary 3, [19]), it leads to the following
inequality:

$\displaystyle\Big{|}\mathcal{L}(t)-\Big{(}2l^{2}+2\mu^{2}+\frac{\tau_{q}k}{\rho
c}\Big{)}E(t)\Big{|}$ $\displaystyle\leq$ $\displaystyle\Big{(}\frac{\rho
c}{2}l^{2}+\frac{\rho
c}{2}\mu^{2}+\frac{\tau_{q}k}{2}\Big{)}\|T\|^{2}+\frac{\tau_{q}}{2k}l^{2}\|q\|^{2}$
(44) $\displaystyle\leq$
$\displaystyle\Big{(}l^{2}+\mu^{2}+\frac{\tau_{q}k}{\rho c}\Big{)}\frac{\rho
c}{2}\|T\|^{2}+\frac{\tau_{q}}{2k}l^{2}\|q\|^{2}$ $\displaystyle\leq$
$\displaystyle\Big{(}l^{2}+\mu^{2}+\frac{\tau_{q}k}{\rho c}\Big{)}E(t).$

Consequently,

$\displaystyle\Big{(}l^{2}+\mu^{2}\Big{)}E(t)\leq\mathcal{L}(t)\leq\Big{(}3l^{2}+3\mu^{2}+2\frac{\tau_{q}k}{\rho
c}\Big{)}E(t),\quad\mbox{for all}\quad t\geq 0.$ (45)

Using Eqs. (14) and (41), we have

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\mathcal{L}(t)\leq-k\|T\|^{2}-\frac{2}{k}\big{(}l^{2}+\mu^{2}\big{)}\|q\|^{2}+|C_{T}(t)|.$
(46)

We choose $\beta:=\min\big{\\{}2k/\rho c,\,4(l^{2}+\mu^{2})/\tau_{q}\big{\\}}$
such that,

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\mathcal{L}(t)\leq-\beta
E(t)+|C_{T}(t)|.$

Utilizing the second inequality in (45), we obtain

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\mathcal{L}(t)\leq-\omega\mathcal{L}(t)+|C_{T}(t)|,\quad\mbox{for
all}\quad t\geq 0,$ (47)

where $\omega:=\frac{\beta}{3l^{2}+3\mu^{2}+2\tau_{q}k/\rho c}$. Multiplying
$(\ref{est.L})$ by $e^{\omega t}$ yields

$\displaystyle\frac{\textrm{d}}{\textrm{d}t}\Big{(}\mathcal{L}(t)e^{\omega
t}\Big{)}\leq|C_{T}(t)|e^{\omega t}.$ (48)

Integrating (48) over $[0,t]$,

$\displaystyle\mathcal{L}(t)\leq\mathcal{L}(0)e^{-\omega
t}+\int_{0}^{t}|C_{T}(r)|e^{-\omega(t-r)}\textrm{d}r,\quad\mbox{for all}\quad
t\geq 0,$ (49)

where $r\in[0,t]$. Choosing $M:=\frac{3l^{2}+3\mu^{2}+2\tau_{q}k/\rho
c}{l^{2}+\mu^{2}}>1$, $\gamma_{0}:=(l^{2}+\mu^{2})^{-1}$ and considering (45),
we finally obtain

$\displaystyle E(t)\leq ME(0)e^{-\omega
t}+\gamma_{0}\int_{0}^{t}|C_{T}(r)|e^{-\omega(t-r)}\textrm{d}r.$ (50)

Since $C_{T}\in L^{\infty}\big{(}(0,\infty),L^{2}(0,l)\big{)}$ and
$T_{0}\notin L_{*}^{2}(0,l)$ we have

$\displaystyle E(t)$ $\displaystyle\leq$ $\displaystyle ME(0)e^{-\omega
t}+(1-e^{-\omega t})\frac{\gamma_{0}}{\omega}\|C_{T}\|_{\infty}$ (51)
$\displaystyle\leq$ $\displaystyle ME(0)e^{-\omega
t}+2\frac{\gamma_{0}}{\omega}\|C_{T}\|_{\infty}.$

This proves $(i)$. On the other hand, if $T_{0}\in L_{*}^{2}(0,l)$, implies
that $|C_{T}(t)|=0$. From (50), it immediately leads to the proof of $(ii)$.

######  Remark 4.2.

The dynamics of the system (1)–(4) determined by the GK equation, can be
reduced to the dynamics of the system determined by the Fourier equation when
$(\tau_{q},\mu^{2})\rightarrow(0,0)$, i.e.,

$\displaystyle\rho cT_{t}+q_{x}$
$\displaystyle=0\quad\mbox{in}\quad(0,l)\times(0,\infty),$ (52) $\displaystyle
q+kT_{x}$ $\displaystyle=0\quad\mbox{in}\quad(0,l)\times(0,\infty).$ (53)

This means that the results of Theorems 3.1 and 4.1 are still valid in the
limit as $(\tau_{q},\mu^{2})\rightarrow(0,0)$.

##  5 A fully discrete finite difference model

In this section, we consider an implicit time integration method by finite
differences applied to (1)–(4) system in $(0,l)\times(0,\mathbf{T})$. We are
mainly concerned with showing that numerical energy preserves the properties
of the energy functional (13) from being positive and non-increasing. For this
reason, consider $J,N\in\mathbb{N}$.

We set $\Delta x=\displaystyle\frac{l}{J+1}$ and $\Delta
t=\displaystyle\frac{\mathbf{T}}{N+1}$ and introduce the nets

$\displaystyle 0=x_{0}<x_{1}=\Delta x<...<x_{j}=j\Delta
x<x_{J}<x_{J+1}=(J+1)\Delta x=l,$ (54) $\displaystyle 0=t_{0}<t_{1}=\Delta
t<...<t_{n}=n\Delta t<t_{N}<t_{N+1}=(N+1)\Delta t=\mathbf{T}.$ (55)

We consider the following finite-diference discretization of (1)–(4). For
$j=1,2,...,J$ and $n=1,2,...,N$ we have

$\displaystyle\rho c\frac{T_{j}^{n}-T_{j}^{n-1}}{\Delta
t}+\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta x}=0,$ (56)
$\displaystyle\tau_{q}\frac{q_{j}^{n}-q_{j}^{n-1}}{\Delta
t}+q_{j}^{n}-\mu^{2}\frac{q_{j+1}^{n}-2q_{j}^{n}+q_{j-1}^{n}}{\Delta
x^{2}}+k\frac{T_{j}^{n}-T_{j-1}^{n}}{\Delta x}=0.$ (57)

We denote by $(T_{j}^{n},q_{j}^{n})$ the numerical approximation for the
solutions $(T,q)$ in points $(x_{j},t_{n})$ of the mesh. For boundary
conditions, we use

$\displaystyle q_{0}^{n}=q_{J+1}^{n}=0,\quad\mbox{for all}\quad n=0,1,...,N+1$
(58)

and for the for the initial conditions, we adopted

$\displaystyle T_{j}^{0}=T_{0j},\quad q_{j}^{0}=q_{0j},\quad\mbox{for
all}\quad j=0,1,...,J+1.$ (59)

###  5.1 Discrete functional energy

In this subsection we analyze in detail some properties of the energy of the
discrete system (56)–(59). We define the functional energy of (56)–(59) by

$\displaystyle\displaystyle\mathbf{E^{n}}:=\rho c\frac{\Delta
x}{2}\sum_{j=0}^{J}|T_{j}^{n}|^{2}+\frac{\tau_{q}}{k}\frac{\Delta
x}{2}\sum_{j=0}^{J}|q_{j}^{n}|^{2},\quad\mbox{for all}\quad n=1,...,N+1.$ (60)

######  Theorem 5.1.

Let $(T_{j}^{n},q_{j}^{n})$ the numerical solution of the problem (56)–(59).
Thus, for all $\Delta x,\Delta t\in(0,1)$ the rate of change of the numerical
energy $\mathbf{E^{n}}$ in $(\ref{numerical.energy})$ at the instant $t_{n}$
is given by

$\displaystyle\frac{\mathbf{E^{n}}-\mathbf{E^{n-1}}}{\Delta
t}\leq-\frac{1}{k}\Delta
x\sum_{j=0}^{J}|q_{j}^{n}|^{2}-\frac{\mu^{2}}{k}\Delta
x\sum_{j=0}^{J}\bigg{|}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta x}\bigg{|}^{2}<0,$
(61)

for all $n=1,2,...,N+1$.

Proof:  Multiplying the Eq. $(\ref{Eq.num1})$ by $\Delta xT_{j}^{n}$ and
adding to $j=0,1,...,J$ we have

$\displaystyle\rho c\Delta
x\sum_{j=0}^{J}\bigg{(}\frac{T_{j}^{n}-T_{j}^{n-1}}{\Delta
t}T_{j}^{n}\bigg{)}+\Delta
x\sum_{j=0}^{J}\bigg{(}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta
x}T_{j}^{n}\bigg{)}=0.$ (62)

Using Young’s inequality we have

$\displaystyle\rho c\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|T_{j}^{n}|^{2}-\rho c\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|T_{j}^{n-1}|^{2}+\Delta
x\sum_{j=0}^{J}\bigg{(}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta
x}T_{j}^{n}\bigg{)}\leq 0.$ (63)

Similarly, multiplying the Eq. $(\ref{Eq.num2})$ by $k^{-1}\Delta xq_{j}^{n}$,
adding to $j=1,...,J$ and using the boundary conditions $q_{0}^{n}=0$ we have

$\displaystyle\frac{\tau_{q}}{k}\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|q_{j}^{n}|^{2}-\frac{\tau_{q}}{k}\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|q_{j}^{n-1}|^{2}+\frac{1}{k}\Delta
x\sum_{j=0}^{J}|q_{j}^{n}|^{2}$
$\displaystyle\underbrace{-\frac{\mu^{2}}{k}\Delta
x\sum_{j=1}^{J}\bigg{(}\frac{q_{j+1}^{n}-2q_{j}^{n}+q_{j-1}^{n}}{\Delta
x^{2}}q_{j}^{n}\bigg{)}+\Delta
x\sum_{j=1}^{J}\bigg{(}\frac{T_{j}^{n}-T_{j-1}^{n}}{\Delta
x}q_{j}^{n}\bigg{)}}_{J:=}\leq 0.$ (64)

Using the boundary conditions (58) we have

$\displaystyle J$ $\displaystyle:=$ $\displaystyle-\frac{\mu^{2}}{k}\Delta
x\sum_{j=1}^{J}\bigg{(}\frac{q_{j+1}^{n}-2q_{j}^{n}+q_{j-1}^{n}}{\Delta
x^{2}}q_{j}^{n}\bigg{)}+\Delta
x\sum_{j=1}^{J}\bigg{(}\frac{T_{j}^{n}-T_{j-1}^{n}}{\Delta
x}q_{j}^{n}\bigg{)}$ (65) $\displaystyle=$
$\displaystyle\frac{\mu^{2}}{k}\Delta
x\sum_{j=0}^{J}\bigg{|}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta
x}\bigg{|}^{2}-\Delta
x\sum_{j=0}^{J}\bigg{(}T_{j}^{n}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta
x}\bigg{)}.$

From (5.1) and (65) we obtain

$\displaystyle\frac{\tau_{q}}{k}\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|q_{j}^{n}|^{2}-\frac{\tau_{q}}{k}\frac{\Delta x}{2\Delta
t}\sum_{j=0}^{J}|q_{j}^{n-1}|^{2}+\frac{1}{k}\Delta
x\sum_{j=0}^{J}|q_{j}^{n}|^{2}+\frac{\mu^{2}}{k}\Delta
x\sum_{j=0}^{J}\bigg{|}\frac{q_{j}^{n+1}-q_{j}^{n}}{\Delta x}\bigg{|}^{2}$
$\displaystyle-\Delta
x\sum_{j=0}^{J}\bigg{(}T_{j}^{n}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta
x}\bigg{)}\leq 0.$ (66)

Adding the inequalities (63) and (5.1) and using the energy $\mathbf{E^{n}}$
in $(\ref{numerical.energy})$ we obtain

$\displaystyle\frac{\mathbf{E^{n}}-\mathbf{E^{n-1}}}{\Delta
t}\leq-\frac{1}{k}\Delta
x\sum_{j=0}^{J}|q_{j}^{n}|^{2}-\frac{\mu^{2}}{k}\Delta
x\sum_{j=0}^{J}\bigg{|}\frac{q_{j+1}^{n}-q_{j}^{n}}{\Delta x}\bigg{|}^{2},$
(67)

and therefore $\mathbf{E^{n}}\leq\mathbf{E^{0}}$ for all $n=1,2,...,N+1$.

##  6 Numerical simulation

Let us rewrite the scheme (56)–(59) in an equivalent vectorial form. We denote
by L and $\textbf{I}_{m}$ the tridiagonal matrix associated to the
approximation of the Laplacian and $m\times m$ identity matrix respectively,
i.e.,

$\textbf{L}:=\left(\begin{array}[]{cccccc}-2&1&0&0&\cdots&0\\\
1&-2&1&\ddots&\ddots&\colon\\\ 0&1&\ddots&\ddots&\ddots&0\\\
0&\ddots&\ddots&\ddots&1&0\\\ \colon&\ddots&\ddots&1&-2&1\\\
0&\cdots&0&0&1&-2\end{array}\right)_{J\times
J},\qquad\textbf{I}_{m}:=\left(\begin{array}[]{cccccc}1&0&0&0&\cdots&0\\\
0&1&0&\ddots&\ddots&\colon\\\ 0&0&\ddots&\ddots&\ddots&0\\\
0&\ddots&\ddots&\ddots&0&0\\\ \colon&\ddots&\ddots&0&1&0\\\
0&\cdots&0&0&0&1\end{array}\right)_{m\times m}.$

On the other hand, $\textbf{A}_{q}$ and $\textbf{A}_{T}$ are trapezoidal
matrices given by

$\textbf{A}_{q}:=\left(\begin{array}[]{cccccc}1&0&0&0&\cdots&0\\\
-1&1&0&\ddots&\ddots&\colon\\\ 0&-1&\ddots&\ddots&\ddots&0\\\
0&\ddots&\ddots&\ddots&0&0\\\ \colon&\ddots&\ddots&-1&1&0\\\
0&\cdots&0&0&-1&1\\\ 0&\cdots&0&0&0&-1\end{array}\right)_{J+1\times
J},\qquad\textbf{A}_{T}:=\left(\begin{array}[]{ccccccc}-1&1&0&0&\cdots&0&0\\\
0&-1&1&\ddots&\ddots&\colon&0\\\ 0&0&\ddots&\ddots&\ddots&0&0\\\
0&\ddots&\ddots&\ddots&1&0&0\\\ \colon&\ddots&\ddots&0&-1&1&0\\\
0&\cdots&0&0&0&-1&1\end{array}\right)_{J\times J+1}.$

Hence we define the mass matrices

$\displaystyle\textbf{B}:=\textbf{I}_{J}-\frac{\mu^{2}}{\tau_{q}+\Delta
t}\frac{\Delta t}{\Delta
x^{2}}\textbf{L},\qquad\textbf{C}:=\textbf{I}_{J+1}+\frac{k}{\rho
c(\tau_{q}+\Delta t)}\frac{\Delta t}{\Delta
x^{2}}\textbf{A}_{q}\textbf{B}^{-1}\textbf{A}_{T},\qquad\textbf{D}:=\textbf{A}_{q}\textbf{B}^{-1}.$
(68)

If we denote $\mathbb{T}^{n}=(T_{0}^{n},T_{1}^{n},...,T_{J}^{n})^{\top}$ and
$\mathbb{Q}^{n}=(q_{1}^{n},q_{2}^{n},...,q_{J}^{n})^{\top}$, then scheme
(56)–(59) takes the following vectorial form:

$\displaystyle\begin{cases}\displaystyle\mathbb{T}^{n}=\textbf{C}\mathbb{T}^{n-1}-\frac{\tau_{q}}{\rho
c(\tau_{q}+\Delta t)}\frac{\Delta t}{\Delta x}\textbf{D}\mathbb{Q}^{n-1},\quad
n=1,2,...,N+1,\\\ \displaystyle\mathbb{Q}^{n}=\frac{\tau_{q}}{\tau_{q}+\Delta
t}\textbf{B}^{-1}\mathbb{Q}^{n-1}-\frac{k}{\tau_{q}+\Delta t}\frac{\Delta
t}{\Delta x}\textbf{B}^{-1}\textbf{A}_{T}\mathbb{T}^{n-1},\quad
n=1,2,...,N+1,\\\
\displaystyle\mathbb{T}^{0}=(T_{0}^{0},T_{1}^{0},...,T_{J}^{0})^{\top},\quad\mathbb{Q}^{0}=(q_{1}^{0},q_{2}^{0},...,q_{J}^{0})^{\top}.\end{cases}$
(69)

Let us present a numerical example, especially with regard to the behavior of
discrete energy $\textbf{E}^{n}$ given by (60), with aiming to confirm the
analytical results established in previous sections. For our purposes, let us
consider $l=1\cdot 10^{-1}$ m and for the computational mesh we take $\Delta
x=2\cdot 10^{-4}$ m, $\Delta t=1.2\cdot 10^{-2}$ s. We use $\rho=2\cdot
10^{3}$ kg/m3, $c=5\cdot 10^{2}$ J/(kg K), $\tau_{q}=8\cdot 10^{-3}$ s,
$\mu^{2}=2.8\cdot 10^{-3}$ m2, $k=2\cdot 10^{3}$ W/(m K) and we take the
following initial conditions associated to the discrete system (69), following
[17]:

$\displaystyle T_{j}^{0}:=T_{b}+\frac{T_{f}}{2}\cos\Big{(}\frac{\pi
x_{j}}{l}\Big{)},\quad{j=0,1,...,J}.$ (70)

Below we provide the simulation of the solution and the energy functional
obtained from the numerical scheme (69), considering two cases in which we
solve the Fourier and Guyer-Krumhansl heat equations in regard to the Theorem
4.1 (item (i)). Since Theorem 4.1 (item (ii)) fals beyond the physical
interest, and also appears as a special case of item (i), we do not deal with
this situation. Although from the thermodynamic point of view the following
simulations are straightforward, we present them in order to support the
mathematical proofs, which properties are also preserved in the discrete case.
Naturally, the solutions of the Guyer-Krumhansl and Fourier heat equations
converges to the same equilibrium state with the same energy content, that
stands as a thermodynamic requirement.

Case I: We consider $T_{j}^{0}:=T_{b}+\frac{T_{f}}{2}\cos\big{(}\frac{\pi
x_{j}}{l}\big{)}$ with $T_{b}=1.5\cdot 10^{1}$ ∘C and $T_{f}=3\cdot 10^{1}$
∘C. This implies that $T_{j}^{0}\notin L_{*}^{2}(0,l)$. According to the
Theorem 4.1 (item (i)) the energy functional decays exponentially to a
positive constant. This is present in the simulation below.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 1: Guyer–Krumhansl heat equation. We use $\tau_{q}=8\cdot 10^{-3}$s,
$\mu^{2}=2.8\cdot 10^{-3}$ m2 and
$Z^{n}:=1+(M_{1}\|C_{T}\|_{\infty}\big{/}M_{0}E^{0})e^{\omega t_{n}}$. As
$T_{0}\notin L_{*}^{2}(0,l)$ the temperature history reaches an equilibrium
point (see (a) and (b)), exponential decay tends to a positive constant
$y\approx 1.24\cdot 10^{7}$.

The simulations for the Fourier heat equation are also obtained from the
numerical scheme (69) by choosing $(\tau_{q},\mu^{2})=(0,0)$.

Case II: Here, we consider the data from Case I with $\tau_{q}=\mu^{2}=0$.
According to the Theorem 4.1 (item (i)) the energy functional decays
exponentially to a positive constant.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 2: Solutions for the Fourier heat equation, presenting the agreement
between the physical and mathematical requirements in regard to the Theorem
4.1 (item (i)).

##  7 Discussion

It is known for decades that Fourier’s law has limitations, thus it is
inevitable to find its proper extension and discover the possibilities for
practical implementations. In this respect, the Guyer-Krumhansl equation is a
promising candidate. While the GK equation has a more complex structure than
Fourier’s law, its modeling capabilities can cover wide range of phenomena
from low-temperature to room temperature problems. Numerous experiments show
the validity of the GK equation in these cases. In order to take a step
forward to have the GK equation to be the next standard model, we presented
three crucial properties of that model.

First, we have seen that for temperature-dependent situations, the
thermodynamic background must be exploited as it connects the coefficients and
they are not independent of each other. It influences the entire model, and
introduces new terms into the consitutive equation. The discovery and
understanding of these new terms will help the experimental studies in the
future how and in what aspect to focus on the measurements and how to evaluate
them. Second, we have proved the well-posedness of the GK equation, that is,
its solution exist, unique and continously depends on the initial data. This
is crucial to have a strong mathematical and physical basis, and avoid the
issues which are characteristic for the DPL concept or any other models
lacking the proper background. Moreover, we also proved the uniform
stabilization for a practically relevant case in which the initial temperature
distribution $T_{0}\notin L_{*}^{2}(0,l)$, i.e., the initial energy content of
the system is not zero. Our statements are demonstrated through numerical
solutions, the proved properties are preserved in the discrete space as well.
These results establish a strong basis for the Guyer-Krumhansl heat equation
and for its later applications.

## Funding

A. J. A. Ramos was partially supported by CNPq Grant 310729/2019-0.

M. M. Freitas was partially supported by CNPq Grant 313081/2021-2.

D. S. Almeida Júnior was partially supported by CNPq Grant 314273/2020-4.

R. Kovács: The research reported in this paper is part of project no. BME-
NVA-02, implemented with the support provided by the Ministry of Innovation
and Technology of Hungary from the National Re-search, Development and
Innovation Fund, financed under the TKP2021 funding scheme. This paper was
supported by the János Bolyai Research Scholarship of the Hungarian Academy of
Sciences. The research reported in this paper and carried out at BME has been
supported by the grants National Research, Development and Innovation Office-
NKFIH FK 134277.

## Declarations

Conflict of interest The author declares no competing interests.

## References

  * [1] R. A. Guyer, J. A. Krumhansl, Thermal Conductivity, Second Sound, and Phonon Hydrodynamic Phenomena in Nonmetallic Crystals, Physical Review 148 (1966) 778–788.  doi:10.1103/PhysRev.148.778.
  * [2] R. A. Guyer, J. A. Krumhansl, Dispersion relation for second sound in solids, Physical Review 133 (5A) (1964) A1411.
  * [3] W. Dreyer, H. Struchtrup, Heat pulse experiments revisited, Continuum Mechanics and Thermodynamics 5 (1993) 3–50.
  * [4] P. Ván, T. Fülöp, Universality in heat conduction theory – weakly nonlocal thermodynamics, Annalen der Physik (Berlin) 524 (8) (2012) 470–478.
  * [5] M. Grmela, Generic guide to the multiscale dynamics and thermodynamics, Journal of Physics Communications 2 (3) (2018) 032001.
  * [6] M. Szücs, M. Pavelka, R. Kovács, T. Fülöp, P. Ván, M. Grmela, A case study of non-fourier heat conduction using Internal Variables and GENERIC, Journal of Non-Equilibrium Thermodynamics (2021).
  * [7] G. Lebon, From classical irreversible thermodynamics to extended thermodynamics, Acta Physica Hungarica 66 (1-4) (1989) 241–249.
  * [8] D. Jou, J. Casas-Vázquez, G. Lebon, Extended Irreversible Thermodynamics, Springer Verlag, Berlin, 2010 (fourth Edition), 2010.
  * [9] A. Fehér, N. Lukács, L. Somlai, T. Fodor, M. Szücs, T. Fülöp, P. Ván, R. Kovács, Size effects and beyond-fourier heat conduction in room-temperature experiments, Journal of Non-Equilibrium Thermodynamics 46 (2021) 403–411.
  * [10] R. Kovács, P. Ván, Generalized heat conduction in heat pulse experiments, International Journal of Heat and Mass Transfer 83 (2015) 613 – 620.
  * [11] S. Both, B. Czél, T. Fülöp, G. Gróf, Á. Gyenis, R. Kovács, P. Ván, J. Verhás, Deviation from the Fourier law in room-temperature heat pulse experiments, Journal of Non-Equilibrium Thermodynamics 41 (1) (2016) 41–48.
  * [12] P. Ván, A. Berezovski, T. Fülöp, G. Gróf, R. Kovács, Á. Lovas, J. Verhás, Guyer-Krumhansl-type heat conduction at room temperature, EPL 118 (5) (2017) 50005, arXiv:1704.00341v1.
  * [13] R. Kovács, P. Rogolino, Numerical treatment of nonlinear Fourier and Maxwell-Cattaneo-Vernotte heat transport equations, International Journal of Heat and Mass Transfer 150 (2020) 119281.
  * [14] M. Dreher, R. Quintanilla, R. Racke, Ill-posed problems in thermomechanics, Applied Mathematics Letters 22 (9) (2009) 1374–1379.
  * [15] S. A. Rukolaine, Unphysical effects of the dual-phase-lag model of heat conduction, International Journal of Heat and Mass Transfer 78 (2014) 58–63.
  * [16] S. A. Rukolaine, Unphysical effects of the dual-phase-lag model of heat conduction: higher-order approximations, International Journal of Thermal Sciences 113 (2017) 83–88.
  * [17] R. Kovács, Analytical treatment of nonhomogeneous initial states for non-Fourier heat equations. International Communications in Heat and Mass Transfer, 134, 106021, (2022).
  * [18] A. Pazy, Semigroups of Linear Operators and Applications to Partial Differential Equations. Springer, New York, (1983).
  * [19] C. A. Swanson, Wirtinger’s Inequality. SIAM Journal on Mathematical Analysis, 9(3), 484-491, (1978)