# Period two solution for a class of distributed delay differential equations

Yukihiko Nakata  Department of Mathematical Sciences, Aoyama Gakuin
University, 5-10-1 Fuchinobe, Chuo-ku, Sagamihara, Kanagawa 252-5258 Japan

###### Abstract.

We study the existence of a periodic solution for a differential equation with
distributed delay. It is shown that, for a class of distributed delay
deferential equations, a symmetric period 2 solution, where the period is
twice the maximum delay, is given as a periodic solution of a Hamiltonian
system of ordinary differential equations. Proof of the idea is based on
(Kaplan & Yorke, 1974, J. Math. Anal. Appl.) for a discrete delay differential
equation with an odd nonlinear function. To illustrate the results, we present
distributed delay differential equations that have periodic solutions
expressed in terms of the Jacobi elliptic functions.

###### Key words and phrases:

delay differential equation, periodic solution, Hamiltonian system,

##  1\. Introduction

In the study of delay differential equations, it has been shown that delay is
responsible to cause oscillation of the solution. Periodic solutions of delay
differential equations (DDEs) have been extensively studied, see e.g. [8] and
references therein.

In the paper [11], Kaplan and Yorke consider the existence of a symmetric
periodic solution of period 4 to the following DDE:

(1.1) $x^{\prime}(t)=-f\left(x(t-1)\right),$

where $f$ is a certain odd function from a subset of $\mathbb{R}$ to
$\mathbb{R}$. The authors pay attention to a periodic solution that satisfies
the symmetry:

$x(t)=-x(t+2).$

Thus the period of the solution is $4$. It is shown that the following planar
system of ordinary differential equations (ODEs):

(1.2) $x^{\prime}=-f(y),\ y^{\prime}=f(x)$

yields such a symmetric periodic solution to DDE (1.1).

The study [11] has been extended to a scalar DDE, including equation (1.1) as
a special case, in [9, 10]. The proof of Kaplan and Yorke in [11] is
elaborated by Nussbaum in [18]. Several authors extend the result to, for
example, a scalar differential equation with multiple discrete delays [4, 5]
and a system of differential equation with multiple discrete delays [6, 7]. In
these papers, variational approach is an essential tool to investigate the
existence of symmetric periodic solutions for Hamiltonian systems of higher
order. Recently, the author in [13] shows that any periodic solution of a
scalar DDE with monotone and symmetric feedback function, including (1.1), is
induced by the corresponding planar system of ODEs.

In the paper [15], the author studies the existence of a period 2 solution for
a logistic type differential equation with distributed delay. A similar ansatz
to the one in Kaplan and Yorke [11] yields a system of ODEs, from which we
obtain a periodic solution to the logistic DDE. The DDE turns out to be an
example that has a periodic solution expressed in terms of the Jacobi elliptic
functions ([14]).

It is promised to generalize the idea to a class of differential equations
with distributed delay, which we call distributed DDE. Our concern in this
paper is the following distributed DDE:

(1.3) $x^{\prime}(t)=-\int_{0}^{1}f(x(t-s))ds,$

where $f:D\subseteq\mathbb{R}\to\mathbb{R}$ is a Lipschitz continuous
function. If $f$ is a certain odd function, as in [11], it is shown that a
symmetric periodic solution is given as the solution of the following system
of ODEs:

(1.4a) $\displaystyle x^{\prime}(t)$ $\displaystyle=-y(t),$ (1.4b)
$\displaystyle y^{\prime}(t)$ $\displaystyle=2f(x(t)).$

The system of ODEs is a planar Hamiltonian system with the following
Hamiltonian function

(1.5) $y^{2}+4\int_{0}^{x}f(\xi)d\xi$

(see [19] and reference therein). For some particular nonlinear function
$f(x)$, explicit form of the periodic solution for (1.3) can be given in terms
of the Jacobi elliptic function.

The paper is organized as follows. In Section 2, we first characterize a
symmetric period $2$ solution of (1.3). It is shown that period 2 solution can
yield a family of periodic solutions of for a distributed DDE with a
multiplicative parameter. In Section 3, we introduce the notion of “Special
Symmetric Periodic Solution (SSPS)”, which is analogue to the symmetric
periodic solution of (1.1) in [11] for a discrete DDE. We show that SSPS of
(1.3) is related to a periodic solution of Hamiltonian system of ODEs (1.4),
see Proposition 5. In Section 3, we introduce 2 particular examples where the
period 2 solution is expressed in terms of the Jacobi elliptic functions.
Detailed computations and materials are summarized in Appendix.

##  2\. Periodic solution of period $2$

The following Proposition may motivate the study of the period 2 solution for
(1.3) as a special periodic solution.

######  Proposition 1.

Let $x:\mathbb{R}\to\mathbb{R}$ be a solution of equation (1.3). The following
statements are equivalent.

  1. (1)

There exists $t_{0}\in\mathbb{R}$ such that $x(t_{0}+t)=x(t_{0}-t),\ \forall
t\in\mathbb{R}$, i.e., $x$ is an even function.

  2. (2)

There exists $c\in\mathbb{R}$ such that

(2.1) $x(t)+x(t-1)=c,\ \forall t\in\mathbb{R}.$

  3. (3)

$x(t)=x(t+2),\ \forall t\in\mathbb{R}$ i.e., the period of $x$ is $2$.

###### Proof.

(1) $\Rightarrow$ (2) From the time translation invariance of the equation
(1.3), without loss of generality, one can assume that $t_{0}=0$. Suppose that
$x(t)=x(-t),\ \forall t\in\mathbb{R}$. From the definition, one obtains

$\displaystyle x^{\prime}(t)$ $\displaystyle=-x^{\prime}(-t)$
$\displaystyle=\int_{0}^{1}f(x(-t-s))ds.$
$\displaystyle=\int_{0}^{1}f(x(t+s))ds\quad(\because x(t)=x(-t))$
$\displaystyle=\int_{0}^{1}f(x(t+1-s))ds$ $\displaystyle=-x^{\prime}(t+1),$

which implies that $x(t)+x(t-1)$ is a constant.

(2) $\Rightarrow$ (1) Let $x(t)+x(t-1)=c$ (Const). Since $x(t-1)=c-x(t)$, from
$f(x(t-1))=f(c-x(t))$, the solution of (1.3) satisfies the following system of
ODEs

(2.2a) $\displaystyle x^{\prime}$ $\displaystyle=-y,$ (2.2b) $\displaystyle
y^{\prime}$ $\displaystyle=f(x)-f(c-x).$

Consider the solution of (2.2) with the initial condition
$(x(t_{0}),y(t_{0}))=(a,0),\ a>0$. One sees that
$\left(\tilde{x}(t+t_{0}),\tilde{y}(t+t_{0})\right)=\left(x(t_{0}-t),-y(t_{0}-t)\right)$
solves (2.2) through
$\left(\tilde{x}(t_{0}),\tilde{y}(t_{0})\right)=(x(t_{0}),y(t_{0}))=(a,0)$.
From the uniqueness of the solution, which is ensured by the assumption that
$f$ is a Lipschitz continuous function, we obtain the conclusion that
$x(t_{0}+t)=x(t_{0}-t)$.

(2) $\Leftrightarrow$ (3) It is obvious that (2) implies (3). We show that (3)
implies (2). From (3), we have $f(x)=f(x(t-2))$, thus

$\int_{t-2}^{t}f(x(s))ds=\text{Const},\ t\in\mathbb{R}.$

Therefore $x^{\prime}(t)+x^{\prime}(t-1)=\text{Const},\ t\in\mathbb{R}$. The
constant must be $0$, otherwise it contradicts to that $x$ is a periodic
solution, which is a bounded solution. ∎

######  Remark 2.

In (3) in Proposition 1, the minimal period of $x$ is not necessary to be $2$.

Next we show that the period 2 solution of (1.3) gives rise to a family of
periodic solutions to the equation of the following form

(2.3) $x^{\prime}(t)=-r\int_{0}^{1}f(x(t-s))ds,$

where $r>0$ is a multiplicative parameter.

######  Proposition 3.

Suppose that equation (1.3) has a periodic solution $x(t)$ of period $2$. Then
for every $n=1,2,\dots$ equation (2.3) with $r=\left(2n-1\right)^{2}$ has a
periodic solution $x((2n-1)t)$ of period $2/(2n-1)$.

###### Proof.

Denote by $x_{1}(t)$ the period 2 solution of equation (1.3). For
$n\in\mathbb{N}$, we let $x_{n}(t):=x_{1}(\alpha_{n}t)$, where
$\alpha_{n}:=2n-1$. From the direct computation,

(2.4)
$x_{n}^{\prime}(t)=\alpha_{n}x_{1}^{\prime}(\alpha_{n}t)=-\alpha_{n}\int_{0}^{1}f(x_{1}(\alpha_{n}t-s))ds.$

Since $x_{1}$ is the period $2$ solution, one obtains

(2.5)
$\int_{0}^{1}f(x_{1}(\alpha_{n}t-s))ds=\int_{0}^{\alpha_{n}}f(x_{1}(\alpha_{n}t-s))ds=\alpha_{n}\int_{0}^{1}f(x_{n}(t-s))ds.$

From 2.5 and 2.4, we see that $x_{n}(t)$ solves equation 2.3 with
$r=\alpha_{n}^{2}=\left(2n-1\right)^{2}$. ∎

##  3\. Distributed DDE (1.3) and Hamiltonian system (1.4)

In this section, we show that distributed DDE (1.3) and Hamiltonian system of
ODEs (1.4) are related in terms of the period 2 solution.

Let $D$ be an open interval of $\mathbb{R}$, containing $0$. We here suppose
the following conditions for the function $f(x)$:

  * (H)

$f(x)=-f(-x),\ xf(x)>0\ \left(x\in D\right)$

Then, analogously to [9, 10], we define a special symmetric periodic solution
(SSPS) for equation (1.3).

######  Definition 4.

Let $x:\mathbb{R}\to\mathbb{R}$ be a solution of equation (1.3). $x$ is called
a special symmetric periodic solution (SSPS) if $x$ is a periodic solution
with the minimal period $2$, i.e.,

(3.1) $x(t)=x(t-2),\ \forall t\in\mathbb{R}$

and satisfies

(3.2) $x(t)+x(t-1)=0,\ \forall t\in\mathbb{R}.$

Alike the well-known connection concerning the symmetric periodic solution
between DDE (1.1) and the Hamiltonian system (1.2), established in Kaplan and
Yorke [11], we here show the following connection between distributed DDE
(1.3) and Hamiltonian system (1.4).

######  Theorem 5.

Suppose that (H) holds. The following statements are true.

  1. (1)

Let us suppose that $x:\mathbb{R}\to\mathbb{R}$ is a SSPS of equation (1.3).
Then

(3.3) $\left(x(t),y(t)\right)=\left(x(t),\int_{0}^{1}f(x(t-s))ds\right),\
t\in\mathbb{R}$

is a periodic solution in $D\times\mathbb{R}$ with minimal period $2$ for
Hamiltonian system (1.4).

  2. (2)

Let $\left\\{\left(x(t),y(t)\right):t\in\mathbb{R}\right\\}\subset
D\times\mathbb{R}$ be a closed orbit with minimal period $2$ around the origin
for Hamiltonian system (1.4). Then $x(t),\ t\in\mathbb{R}$ is a SSPS of
equation (1.3).

Thus equation (1.3) has a SSPS if and only if system (1.4) has a closed orbit
in $D\times\mathbb{R}$ with minimal period $2$.

###### Proof.

(1) Let $x$ be a SSPS of equation (1.3). Define

$y(t):=\int_{0}^{1}f(x(t-s))ds.$

Then

$y^{\prime}(t)=f(x(t))-f(x(t-1)).$

From the condition (H) and the property (3.2), one has

$-f(x(t-1))=f(-x(t-1))=f(x(t)).$

Thus, one sees that $\left(x(t),y(t)\right)$ is a solution of Hamiltonian
system (1.4). Then it is obvious that $\left(x(t),y(t)\right)$ is a closed
symmetric orbit around the origin of minimal period $2$.

(2) Let $\left\\{\left(x(t),y(t)\right):t\in\mathbb{R}\right\\}$ be a closed
symmetric trajectory around the origin of minimal period $2$ for system (1.4).
Since $\left\\{\left(-x(t),-y(t)\right):t\in\mathbb{R}\right\\}$ is a same
trajectory, one can find that

$\left(x(t+1),y(t+1)\right)=\left(-x(t),-y(t)\right),\ t\in\mathbb{R}.$

Since $x(t)$ is of period $2$, the property $x(t-1)=-x(t)$ holds. Therefore,
it follows that, using the condition (H), $f(x(t))=-f(-x(t))=-f(x(t-1))$, thus

$y^{\prime}(t)=2f(x(t))=f(x(t))-f(x(t-1)),$

from which we obtain $y(t)=\int_{t-1}^{t}f(x(s))ds+\text{Const}$. We show that
the constant is zero. Consider a closed orbit around the origin for system
(1.4) through $x(\frac{1}{2})=0,\ y(\frac{1}{2})>0$. From Lemma 9 in Appendix
A, $\left(x,y\right)=\left(-x\left(-t\right),y(-t)\right)$ also solve system
(1.4). Therefore, from the uniqueness of the solution of (1.4), we have
$\left(x\left(\frac{1}{2}+t\right),y\left(\frac{1}{2}+t\right)\right)=\left(-x\left(\frac{1}{2}-t\right),y\left(\frac{1}{2}-t\right)\right)$.
Thus $x$ is an odd function about $t=\frac{1}{2}$, i.e.,
$x(t+\frac{1}{2})=-x(t-\frac{1}{2})$. It then follows that

$y(t)=\int_{t-1}^{t}f(x(s))ds.$

Consequently, the first component of the solution of equation (1.4) yields an
SSPS $x(t)$ for distributed DDE (1.3). ∎

##  4\. Examples

In this section we study distributed DDE (1.3) with two particular nonlinear
functions, namely $f(x)=r\sin x$ and $f(x)=r(e^{x}-1)$, where $r>0$. We
express the periodic solutions in terms of the Jacobi elliptic functions. See
[2, 14] for the properties of the Jacobi elliptic functions.

###  4.1. $f(x)=r\sin(x)$

In [16], we consider the existence of a period 2 solution to the following
equation

(4.1) $x^{\prime}(t)=-f\left(\int_{0}^{1}x(t-s)ds\right),$

with $f(x)=r\sin x,$ where $r>0$. In Appendix B, we show that equation (4.1)
can be transformed to the form of equation (1.3).

Let us consider the existence of period $2$ solution of (1.3) with $f(x)=r\sin
x$. The corresponding ODE (1.4) with $f(x)=r\sin x$ is a well-known pendulum
equation. We have the following result, see Appendix C for detail.

######  Theorem 6.

Let $f(x)=r\sin x,\ r>0$. Equation (1.3) has a SSPS if and only if
$r>\pi^{2}/2$. The SSPS is expressed as

$x(t)=2\arcsin\left(\sqrt{a}\text{sn}\left(\sqrt{2r}t,\sqrt{a}\right)\right),$

where $a$ is uniquely determined by the solution of the following equation

$\frac{2K\left(\sqrt{a}\right)}{\sqrt{2r}}=1.$

###  4.2. $f(x)=r\left(e^{x}-1\right)$

In [15] we study the following logistic equation with distributed delay

$x^{\prime}(t)=rx(t)\left(1-\int_{0}^{1}x(t-s)ds\right),$

where $r>0$. Changing the variable $x\mapsto\log x$, we obtain equation (1.3)
with $f(x)=r(e^{x}-1)$. Note that $f(x)$ is not an odd function, thus we may
not directly apply Theorem 5 to obtain SSPS of (1.3) with $f(x)=r(e^{x}-1)$.

Consider the existence of periodic solutions with period $2$, satisfying

$x(t)+x(t-1)=\text{Const}.$

Denote by $c$ the constant. The periodic solution satisfies the following
system of ODEs:

(4.2a) $\displaystyle x^{\prime}$ $\displaystyle=-y,$ (4.2b) $\displaystyle
y^{\prime}$ $\displaystyle=f(x)-f(c-x),$

(see also the proof of Proposition 1). Note that $c$ is a parameter to be
determined.

Let $w=x-\frac{1}{2}c$. Then system (4.2) becomes

$\displaystyle w^{\prime}$ $\displaystyle=-y,$ $\displaystyle y^{\prime}$
$\displaystyle=f\left(w+\frac{1}{2}c\right)-f\left(\frac{1}{2}c-w\right).$

For $f(x)=r(e^{x}-1)$, we see

$\displaystyle f\left(w+\frac{1}{2}c\right)-f\left(\frac{1}{2}c-w\right)$
$\displaystyle=r\left(e^{w+\frac{1}{2}c}-1\right)-r\left(e^{-w+\frac{1}{2}c}-1\right)$
$\displaystyle=2re^{\frac{1}{2}c}\sinh w.$

Thus we obtain the following system of ODEs:

(4.3a) $\displaystyle w^{\prime}$ $\displaystyle=-y,$ (4.3b) $\displaystyle
y^{\prime}$ $\displaystyle=2re^{\frac{1}{2}c}\sinh w.$

First we would like to find periodic solutions with period 2 for (4.3).
Consider periodic solutions of (4.3) with the following initial condition

(4.4) $\left(w(0),y(0)\right)=\left(a,0\right),\ a>0.$

Let

$\gamma=re^{\frac{1}{2}c},\quad\alpha=\sinh\left(\frac{a}{2}\right)>0.$

Define

(4.5) $\beta=\sqrt{2\gamma(1+\alpha^{2})},\quad
k=\frac{\alpha}{\sqrt{1+\alpha^{2}}}\in\left(0,1\right).$

Note that

(4.6) $\beta k=\sqrt{2\gamma}\alpha,\quad\alpha=\frac{k}{\sqrt{1-k^{2}}}$

hold. In terms of $\beta$ and $k$, the solution of (4.3) is expressed as
follows.

######  Lemma 7.

The solution of (4.3) with the initial condition (4.4) is expressed as

(4.7a) $\displaystyle w(t)=\log\left(\frac{\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)}{\text{dn}\left(\beta
t,k\right)-k\text{cn}\left(\beta t,k\right)}\right),$ (4.7b) $\displaystyle
y(t)=2\beta k\text{sn}\left(\beta t,k\right).$

The period of the solution (4.7) for (4.3) is 2 if and only if

(4.8) $\displaystyle\beta=2K(k),$

Furthermore, $x=w+\frac{1}{2}c$ satisfies (3.3) with $y(0)=0$ if and only if

$0=-\int_{0}^{1}f(x(-s))ds.$

Finally, we obtain the following result.

######  Theorem 8.

Let $f(x)=r(e^{x}-1),\ r>0$. Equation (1.3) has an SSPS if and only if
$r>\pi^{2}/2$. The SSPS is expressed as

(4.9)
$x(t)=\log\left(\frac{K(k)\left(\text{dn}\left(2K(k)t,k\right)+k\text{cn}\left(2K(k)t,k\right)\right)^{2}}{2E(k)-K(k)(1-k^{2})}\right),$

where $k\in\left(0,1\right)$ is a unique solution of

(4.10) $r=2K(k)\left(2E(k)-K(k)(1-k^{2})\right).$

### Acknowledgement

The author was supported by JSPS Kakenhi Grant Number 20K03734 and 19K21836.

## References

  * [1] D. Breda, O. Diekmann, D. Liessi, F. Scarabel, Numerical bifurcation analysis of a class of nonlinear renewal equations. Electron. J. Qual. Theo. Diff. Equ. 65 (2016) 1–24.
  * [2] P. F. Byrd and M. D. Friedmann, Handbook of Elliptic Integrals for Engineers and Physicists. Springer Berlin (1954)
  * [3] O. Diekmann, S. A. van Gils, S. M. Verduyn Lunel, H. O. Walther, Delay Equations: Functional-, Complex- and Nonlinear Analysis. Springer, New York (1995)
  * [4] G. Fei, G. Multiple periodic solutions of differential delay equations via Hamiltonian systems (I). Nonlinear Anal. Theory Methods Appl. 65 (2006) pp. 25–39.
  * [5] G. Fei, Multiple periodic solutions of differential delay equations via Hamiltonian systems (II). Nonlinear Anal. Theory Methods Appl. 65 (2006) pp. 40–58.
  * [6] B. Zheng, and Z. Guo, Multiplicity results on periodic solutions to higher-dimensional differential equations with multiple delays. Rocky Mt. J. Math., 44 (2014) pp. 1715-1744
  * [7] Z. Guo and J. Yu, Multiplicity results on period solutions to higher dimensional differential equations with multiple delays. J. Dyn. Differ. Equ., 23 (2011) pp. 1029–1052.
  * [8] J. K. Hale, S. M. Verduyn Lunel, Introduction to Functional Differential Equations. Springer, New York (1993)
  * [9] P. Dormayer, A. F. Ivanov, Symmetric periodic solutions of a delay differential equation. Conf. Publ. (Special) (1998) 220–230.
  * [10] P. Dormayer, A.F. Ivanov, Stability of symmetric periodic solutions with small amplitude of $x^{\prime}(t)=\alpha f(x(t),x(t-1))$. Discrete Contin. Dyn. Syst. 5 (1999) pp. 61–82
  * [11] J. L. Kaplan, J. A. Yorke, Ordinary differential equations which yield periodic solutions of differential delay equations. J. Math. Anal. Appl. 48 (1974) pp. 317–324
  * [12] B. Kennedy, Symmetric periodic solutions for a class of differential delay equations with distributed delay. Electron. J. Qual. Theo. Diff. Equ. 4 (2014) pp. 1–18.
  * [13] A. López Nieto, Periodic orbits of delay equations with monotone feedback and even-odd symmetry. arXiv:2002.01313 (2020)
  * [14] K. R. Meyer, Jacobi elliptic functions from a dynamical systems point of view. Am. Math. Mon. 108 (2001) pp. 729–737
  * [15] Y. Nakata, An explicit periodic solution of a delay differential equation. J. Dyn. Diff. Equ. 32 (2020) pp. 163–179
  * [16] Y. Nakata, Existence of a period two solution of a delay differential equation, Disc. Cont. Dyn. Syst. S. 14 (2021) pp. 1103–1110.
  * [17] R. D. Nussbaum, Periodic solutions of some nonlinear autonomous functional differential equations. Ann. di Mat. Pura ed Appl. 101 (1974) pp. 263–306.
  * [18] R. D. Nussbaum, Uniqueness and nonuniqueness for periodic solutions of $x^{\prime}(t)=-g(x(t-1))$. J. Differential. Equ., 34 (1979) pp. 25–54.
  * [19] R. Schaaf, Global Solution Branches of Two Point Boundary Value Problems. Lecture Notes in Mathematics Springer (1990)
  * [20] H. O. Walther, Topics in delay differential equations. Jahresber. Dtsch. Math. Ver., 116 (2014) pp. 87–114.

##  Appendix A Symmetry of the solution of Hamiltonian system (1.4) with odd
nonlinearity

In this section we suppose that (H) holds. Consider Hamiltonian system of ODEs
(1.4).

######  Lemma 9.

Suppose that (H) holds. Let $(x(t),y(t))$ be a solution of (1.4). Then
$\left(-x(t),-y(t)\right),\ \left(x(-t),-y(-t)\right),\
\left(-x\left(-t\right),y(-t)\right)$ are also solutions of (1.4).

By direct substitution, one can prove Lemma 9. We here omit the detail.

The solution of Hamiltonian system (1.4) with the initial condition

(A.1) $(x(0),y(0))=(a,0),\ a\in D$

is a closed orbit in $D\times\mathbb{R}$, satisfying

$y(t)^{2}+4F(x(t))=4F(a),$

where $F(x)=\int_{0}^{x}f(\xi)d\xi$. It is easy to see that the closed orbit
given by (A.1) is symmetric about the $x$-axis and the $y$-axis. The period of
the solution is given as $T,$ where

$T=2\int_{0}^{a}\frac{dx}{\sqrt{F(a)-F(x)}},$

i.e., $\left(x(t),y(t)\right)=\left(x(t+T),y(t+T)\right)$.

######  Proposition 10.

Suppose that (H) holds. Let $(x(t),y(t))$ be a periodic solution of (1.4) with
the initial condition (A.1). Then

  1. (1)

$x(t)$ is an even function of $t$ and $y(t)$ is an odd function of $t$.

  2. (2)

$x\left(\frac{1}{4}T+t\right)$ is an odd function of $t$ and
$y\left(\frac{1}{4}T+t\right)$ is an odd function of $t$.

In particular, one has

$\displaystyle\left(x\left(t\right),y(t)\right)$
$\displaystyle=\left(-x\left(\frac{1}{4}T-t\right),y\left(\frac{1}{4}T-t\right)\right)$
$\displaystyle=\left(-x\left(\frac{1}{4}T+t\right),-y\left(\frac{1}{4}T+t\right)\right)$
$\displaystyle=\left(x\left(T-t\right),-y\left(T-t\right)\right).$

##  Appendix B Equivalent formulation

Denote by $x$ a solution of (1.3). For a solution of (1.3), we define a
function $v:\left[-1,\infty\right)\to\mathbb{R}$ by

(B.1) $v^{\prime}=-f(x)$

with

$x(0)=\int_{0}^{1}v(-s)ds.$

Integrating (B.1), it follows that

$v(t)-v(t-1)=-\int_{t-1}^{t}f(x(s))dt.$

From equation (1.3), one sees that

$v(t)-v(t-1)=x^{\prime}(t).$

Therefore, we find the following equality:

(B.2) $x(t)=\int_{0}^{1}v(t-s)ds.$

Substituting (B.2) into equation (B.1), we obtain the following equation:

$v^{\prime}(t)=-f\left(\int_{0}^{1}v(t-s)ds\right).$

##  Appendix C Explicit periodic solutions

In this section, we derive expressions of periodic solutions of Hamiltonian
system (1.4) for $f(x)=r\sin x,\ r>0$ in Section C.1 and for
$f(x)=r\left(e^{x}-1\right),\ r>0$ in Section C.2. See also [2, 14] for the
fundamental properties of the Jacobi elliptic functions.

###  C.1. $f(x)=r\sin x$

Let $f(x)=r\sin x$. We consider the periodic solution of (1.4) with the
initial condition

(C.1) $\left(x(0),y(0)\right)=\left(a,0\right),\ a\in\left(0,\pi\right).$

From (1.5) we have

(C.2) $y(t)^{2}+8r\sin^{2}\frac{x(t)}{2}=8r\sin^{2}\frac{a}{2}.$

Let $u=\sin\frac{x}{2}$. Using (2.2a) we obtain

(C.3) $\displaystyle\left(u^{\prime}\right)^{2}$
$\displaystyle=2r\left(1-u^{2}\right)\left(k^{2}-u^{2}\right),$

where

$k=\sin\frac{a}{2}\in\left(0,1\right).$

######  Lemma 11.

Let $f(x)=r\sin x,\ r>0$. The periodic solution of (1.4) with the initial
condition (C.1) is expressed as

(C.4) $x(t)=2\arcsin\left(k\text{sn}\left(\sqrt{2r}t+K(k),k\right)\right).$

The period of the solution is given as $T$ where

(C.5) $T=\frac{4K(k)}{\sqrt{2r}}.$

###### Proof.

Equation (C.3) has the solution of the form

$u(t)=k\text{sn}\left(\sqrt{2r}t+c,k\right),$

where $c$ is determined by the initial condition $u(0)=\sin\frac{a}{2}=k$, so
that $\text{sn}\left(c,k\right)=1$. Then we see $c=K(k)$. From the definition,
we have $x=2\arcsin u$. Thus we obtain (C.4). ∎

###  C.2. $f(x)=r\left(e^{x}-1\right)$

From (1.5) we have

(C.6)
$y^{2}+8\gamma\sinh^{2}\left(\frac{w}{2}\right)=8\gamma\sinh^{2}\left(\frac{a}{2}\right),$

using $\cosh(x)=1+2\sinh^{2}x$. Let

(C.7) $u:=\sinh\left(\frac{w}{2}\right).$

Then, from (C.6), we obtain

(C.8)
$\left(u^{\prime}\right)^{2}=\frac{1}{4}\cosh^{2}\left(\frac{w}{2}\right)y^{2}=2\gamma\left(1+u^{2}\right)\left(\alpha^{2}-u^{2}\right).$

######  Lemma 12.

Let $f(x)=r(e^{x}-1)\ \left(x\in\mathbb{R}\right)$. The periodic solution of
(C.8) with the initial condition $u(0)=\alpha$ is expressed as

(C.9) $u(t)=\alpha\text{cn}\left(\beta t,k\right),$

where the period is given as

(C.10) $T=\frac{4K\left(k\right)}{\beta}.$

Here $\beta$ and $k$ are defined in (4.5).

###### Proof.

We consider the solution of the following form:

$u(t)=\sqrt{c_{1}}\text{cn}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right).$

By a direct computation, one sees

$\displaystyle\left(u^{\prime}\right)^{2}$
$\displaystyle=c_{1}c_{2}\text{sn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)\text{dn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)$
$\displaystyle=c_{1}c_{2}\left(1-\text{cn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)\right)\left(1-c_{3}\text{sn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)\right)$
$\displaystyle=c_{1}c_{2}\left(1-\text{cn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)\right)\left(1-c_{3}\left(1-\text{cn}^{2}\left(\sqrt{c_{2}}t,\sqrt{c_{3}}\right)\right)\right)$
$\displaystyle=c_{1}c_{2}\left(1-\frac{1}{c_{1}}u^{2}\right)\left(1-c_{3}+\frac{c_{3}}{c_{1}}u^{2}\right).$

From (C.8), we get

(C.11a) $\displaystyle c_{1}c_{2}\left(1-c_{3}\right)=2\gamma\alpha^{2},$
(C.11b) $\displaystyle
c_{2}\left(c_{3}-\left(1-c_{3}\right)\right)=2\gamma\left(\alpha^{2}-1\right)$
(C.11c) $\displaystyle\frac{c_{2}c_{3}}{c_{1}}=2\gamma.$

From (C.11), eliminating $c_{1}$ and $c_{2}$, we have

$\frac{\alpha^{2}}{\left(\alpha^{2}-1\right)^{2}}=\frac{c_{3}(1-c_{3})}{\left(2c_{3}-1\right)^{2}},$

from which we obtain

$0=c_{3}^{2}-c_{3}+\frac{\alpha^{2}}{\left(\alpha^{2}+1\right)^{2}}=\left(c_{3}-\frac{\alpha^{2}}{1+\alpha^{2}}\right)\left(c_{3}-\frac{1}{1+\alpha^{2}}\right).$

Hence

$c_{3}=\frac{\alpha^{2}}{1+\alpha^{2}},\
\frac{1}{1+\alpha^{2}}\in\left[0,1\right].$

Then, from (C.11b), we have

$c_{2}=\begin{cases}2\gamma(1+\alpha^{2})&\text{if
}c_{3}=\frac{\alpha^{2}}{1+\alpha^{2}},\\\ -2\gamma(1+\alpha^{2})&\text{if
}c_{3}=\frac{1}{1+\alpha^{2}}.\end{cases}$

From positivity of $c_{2}$, one obtains
$c_{3}=\frac{\alpha^{2}}{1+\alpha^{2}}$. Therefore, we now have

$\displaystyle\left(c_{1},c_{2},c_{3}\right)$
$\displaystyle=\left(\alpha^{2},2\gamma(1+\alpha^{2}),\frac{\alpha^{2}}{1+\alpha^{2}}\right).$
$\displaystyle=\left(\alpha^{2},\beta^{2},k^{2}\right).$

Hence we obtain (C.4). The period of the solution is given as in (C.10), from
the property of the Jacobi elliptic function cn. ∎

Let us prove Lemma 7 in the main text.

###### Proof.

(Proof of Lemma 7) From (C.9) in Lemma 12, we have

$\displaystyle w(t)$
$\displaystyle=2\sinh^{-1}\left(\alpha\text{cn}\left(\beta t,k\right)\right)$
$\displaystyle=2\log\left(\alpha\text{cn}\left(\beta
t,k\right)+\sqrt{1+\alpha^{2}\text{cn}^{2}\left(\beta t,k\right)}\right),$

where $\beta$ and $k$ are defined in (4.5). Then, using 4.6, we obtain

$\displaystyle w(t)$ $\displaystyle=2\log\left(\frac{k\text{cn}\left(\beta
t,k\right)+\sqrt{1-k^{2}+k^{2}\text{cn}^{2}\left(\beta
t,k\right)}}{\sqrt{1-k^{2}}}\right)$
$\displaystyle=2\log\left(\frac{k\text{cn}\left(\beta
t,k\right)+\text{dn}\left(\beta t,k\right)}{\sqrt{1-k^{2}}}\right)$
$\displaystyle=\log\left(\frac{\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)}{\text{dn}\left(\beta
t,k\right)-k\text{cn}\left(\beta t,k\right)}\right).$

Hence we obtain (4.7a). From (C.6)

$\displaystyle y(t)^{2}$
$\displaystyle=8\gamma\left(\alpha^{2}-\sinh^{2}\left(\frac{w}{2}\right)\right)$
$\displaystyle=8\gamma\alpha^{2}\left(1-\text{cn}^{2}\left(\beta
t,k\right)\right)$ $\displaystyle=8\gamma\alpha^{2}\text{sn}^{2}\left(\beta
t,k\right).$

Hence we obtain (4.7b). ∎

###### Proof.

(Proof of Theorem 8) From the solution obtained in Lemma 7, we can express
$x(t)$ as

(C.12) $x(t)=\log\left(\frac{\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)}{\text{dn}\left(\beta
t,k\right)-k\text{cn}\left(\beta t,k\right)}\right)+\frac{1}{2}c.$

Eq. (4.8) holds, if and only if the period of (4.7a) and (4.7b) is 2. We thus
assume that (4.8) holds. Then one has

(C.13) $re^{\frac{1}{2}c}=2K^{2}(k)(1-k^{2}).$

Furthermore we require that

$y(t)=\int_{0}^{1}f(x(t-s))ds$

holds for $t\in\mathbb{R}$. Since, from (4.8), $y(1)=2\beta
k\text{sn}\left(\beta,k\right)=0$ follows, we must have

(C.14) $0=\int_{0}^{1}\left(e^{x(t)}-1\right)dt.$

Substituting (C.12) into (C.14), we determine the parameters requiring

$\displaystyle 0$
$\displaystyle=e^{\frac{1}{2}c}\int_{0}^{1}\frac{\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)}{\text{dn}\left(\beta
t,k\right)-k\text{cn}\left(\beta t,k\right)}dt-1.$

We now compute

$\displaystyle\frac{\text{dn}\left(\beta t,k\right)+k\text{cn}\left(\beta
t,k\right)}{\text{dn}\left(\beta t,k\right)-k\text{cn}\left(\beta t,k\right)}$
$\displaystyle=\frac{\left(\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)\right)^{2}}{1-k^{2}}$
$\displaystyle=\frac{2\text{dn}^{2}\left(\beta
t,k\right)+2k\text{dn}\left(\beta t,k\right)\text{cn}\left(\beta
t,k\right)}{1-k^{2}}-1,$

Then, we obtain

$\int_{0}^{1}\frac{\text{dn}\left(\beta t,k\right)+k\text{cn}\left(\beta
t,k\right)}{\text{dn}\left(\beta t,k\right)-k\text{cn}\left(\beta
t,k\right)}dt=\frac{2}{1-k^{2}}\frac{E(k)}{K(k)}-1,$

thus

(C.15)
$e^{\frac{1}{2}c}\left(\frac{2}{1-k^{2}}\frac{E(k)}{K(k)}-1\right)-1=0.$

Here, we use

$\int_{0}^{1}\text{dn}^{2}\left(\beta
t,k\right)dt=\frac{2}{\beta}E(k)=\frac{E(k)}{K(k)},$

from the condition (4.8) (See also 314.02 in [2] for the former identity) and

$\int_{0}^{1}\text{dn}\left(\beta t,k\right)\text{cn}\left(\beta
t,k\right)dt=\frac{1}{\beta}\left[\text{sn}\left(\beta
t,k\right)\right]_{t=0}^{t=1}=0.$

From (C.13) and (C.15), we obtain the equation (4.10). From Lemma 6 in [15],
one can see that (4.10) has a unique solution for $r>\pi^{2}/2$. Then we
obtain that

$\frac{1}{2}c=\log\frac{2K^{2}(k)(1-k^{2})}{r}=\log\frac{K(k)(1-k^{2})}{2E(k)-K(k)(1-k^{2})}.$

Then, from (C.12), it follows that

$\displaystyle x(t)$ $\displaystyle=\log\left(\frac{\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta t,k\right)}{\text{dn}\left(\beta
t,k\right)-k\text{cn}\left(\beta
t,k\right)}\frac{K(k)(1-k^{2})}{2E(k)-K(k)(1-k^{2})}\right)$
$\displaystyle=\log\left(\frac{K(k)\left(\text{dn}\left(\beta
t,k\right)+k\text{cn}\left(\beta
t,k\right)\right)^{2}}{2E(k)-K(k)(1-k^{2})}\right).$

Using (4.8), we obtain the expression (4.9). ∎


# Descriptors of intrinsic hydrodynamic thermal transport:
screening a phonon database in a machine learning approach

Pol Torres  Department of Mechanical Engineering, The University of Tokyo,
7-3-1 Hongo, Bunkyo, Tokyo, 113-8656, Japan  EURECAT, Technology Center of
Catalonia, Applied Artificial Intelligence, 08290 Cerdanyola, Barcelona, Spain
Departament de Física, Universitat Autònoma de Barcelona (UAB), Campus de
Bellaterra, 08193 Bellaterra, Barcelona, Spain  Stephen Wu  Research
Organization of Information and Systems, The Institute of Statistical
Mathematics (ISM), 10-3 Midori-cho, Tachikawa, Tokyo 190-8562, Japan
Shenghong Ju  China-UK Low Carbon Collage, Shanghai Jiao Tong University,
Shanghai 201306, China  Department of Mechanical Engineering, The University
of Tokyo, 7-3-1 Hongo, Bunkyo, Tokyo, 113-8656, Japan  Chang Liu  Center for
Materials Research by Information Integration (CMI2), Research and Services
Division of Materials Data and Integrated System (MaDIS), National Institute
for Materials Science (NIMS), 1-2-1 Sengen, Tsukuba, Ibaraki 305-0047, Japan
Research Organization of Information and Systems, The Institute of Statistical
Mathematics (ISM), 10-3 Midori-cho, Tachikawa, Tokyo 190-8562, Japan  Terumasa
Tadano  Research Center for Magnetic and Spintronic Materials, National
Institute for Materials and Science, Tsukuba, Japan  Ryo Yoshida  Research
Organization of Information and Systems, The Institute of Statistical
Mathematics (ISM), 10-3 Midori-cho, Tachikawa, Tokyo 190-8562, Japan  Center
for Materials Research by Information Integration (CMI2), Research and
Services Division of Materials Data and Integrated System (MaDIS), National
Institute for Materials Science (NIMS), 1-2-1 Sengen, Tsukuba, Ibaraki
305-0047, Japan  Junichiro Shiomi  Department of Mechanical Engineering, The
University of Tokyo, 7-3-1 Hongo, Bunkyo, Tokyo, 113-8656, Japan
<EMAIL_ADDRESS>
###### Abstract

Machine learning techniques are used to explore the intrinsic origins of the
hydrodynamic thermal transport and to find new materials interesting for
science and engineering. The hydrodynamic thermal transport is governed
intrinsically by the hydrodynamic scale and the thermal conductivity. The
correlations between these intrinsic properties and harmonic and anharmonic
properties, and a large number of compositional (290) and structural (1224)
descriptors of 131 crystal compound materials are obtained, revealing some of
the key descriptors that determines the magnitude of the intrinsic
hydrodynamic effects, most of them related with the phonon relaxation times.
Then, a trained black-box model is applied to screen more than 5000 materials.
The results identify materials with potential technological applications.
Understanding the properties correlated to hydrodynamic thermal transport can
help to find new thermoelectric materials and on the design of new materials
to ease the heat dissipation in electronic devices.

##  I INTRODUCTION

From the Boltzmann transport equation (BTE) for phonons it is well known that
the thermal conductivity of a crystalline material is mainly determined by the
phonon scattering. The stronger the scattering is, the higher the thermal
resistance is. Phonons can scatter with each other anharmonically, against
boundaries, with impurities and defects. The BTE states that when the phonon
distribution is perturbed from its equilibrium, it relaxes to the equilibrium
through the scattering events Guyer and Krumhansl (1966a); Esfarjani and
Stokes (2008); Chaput (2013); Torres (2018).

When studying the thermal transport in bulk materials and under slow heating
conditions, the Fourier’s law provides a simple expression for the diffusive
thermal transport:

$\mathbf{Q}=-\mathbf{\kappa}\boldsymbol{\nabla}T\;,$ (1)

where $\mathbf{Q}$ is the heat flux, $\mathbf{\kappa}$ the thermal
conductivity and $T$ the temperature. In recent experiments at the micro/nano
scale and/or using fast and/or large thermal gradient, it has been observed
that not all thermal transport can be explained by the Fourier’s law Johnson
et al. (2013); Siemens et al. (2010); Minnich et al. (2011); Regner et al.
(2013); Wilson and Cahill (2014); Hu et al. (2015); Hoogeboom-Pot et al.
(2015). In such situations, other transport phenomena such as ballistic
transport and Lévy flights Vermeersch et al. (2014); Vermeersch and Shakouri
(2014) or hydrodynamic transport appear Cepellotti et al. (2015); Lee et al.
(2015); Torres et al. (2018). In this article we will focus on the latter
transport regime.

In the second half of the last century, Guyer and Krumhansl introduced an
expression of the thermal transport beyond the Fourier’s law to describe the
hydrodynamic thermal transportGuyer and Krumhansl (1966a, b):

$\tau\dot{\textbf{Q}}+\textbf{Q}=-\kappa\boldsymbol{\nabla}T+\ell^{2}\nabla^{2}\textbf{Q}\;,$
(2)

where $\tau$ is the phonon mean free time and $\ell$ is the hydrodynamic
scale. Both of these are the values averaged over the entire phonon spectrum
Torres (2018); Torres et al. (2018). Eq. (2) is a general hydrodynamic heat
flux equation that includes memory (time derivative term) and non-local
effects (Laplacian term), and it has been used successfully to describe recent
experiments (Ziabari et al., 2018; Torres et al., 2018) that could not be
reproduced only by Fourier’s law.

It is known that the thermal hydrodynamic effects are related with the
relative importance of Normal (N) scattering (where phonon-phonon collisions
conserve momentum) in front of restive (R) scattering (where phonon-phonon
collisions do not conserve momentum, like Uumklapp scattering or mass
deffect), but little is known about the intrinsic origin of these phenomena.
In addition, as described by Eq. (2), the temporal scale and the geometry will
also contribute to this effect. To elucidate the key structural/chemical
factors that determine the hydrodynamic thermal transport, leaving aside the
temporal and geometrical scales, in this paper we employ machine learning
techniques. A Python open source platform of materials informatics that we
have developed called XenonPy R (2020) is used to generate compositional and
sctructural descriptors. By providing only the Materials Project ID of a
certain material Jain et al. (2013); Ong et al. (2015), the tool provides 290
compositional descriptors and 1224 structural descriptors like the
polarizability, lattice constant, and Van Der Waals radius, etc.

##  II METHODOLOGY

In this work we perform full first-principles calculations for 131 materials.
The interatomic force constants from 94 of those materials are provided by A.
Togo Seko et al. (2015), 18 are from the AlmaBTE database Carrete et al.
(2014) and 20 are calculated by ourselves. Our DFT calculations are performed
using the Vienna Ab initio Simulation Package (VASP) in the framework of the
projector augmented wave (PAW) potentials in the Perdew-Burke-Ernzerhof (PBE)
approximation. The energy cutoff was fixed at a value 10% higher than the
value specified in the pseudopotential file. The atom displacements to
calculate the harmonic and anharmonic interatomic force constants are provided
by the Phonopy Togo and Tanaka (2015) and Phono3py Togo et al. (2015) codes
respectively.

The parameters that describe the hydrodynamic thermal transport from Eq. (2),
mainly $\kappa$, $\tau$ and $\ell$ are computed from the first-principles
calculations beyond the relaxation-time approximation (RTA) by using Phono3py
and the KCM code Torres et al. (2020). All the calculations are performed at
300 K. From the calculations, we have obtained extra first-principles
descriptors like the averaged Normal and resistive (Umklapp plus mass
variance/impurity) relaxation times ($\tau_{N}$, $\tau_{R}$), specific heat
$C$, and the group velocity $v$. In recent articles, it has been described
that other parameters such as the maximum frequency of the phonon dispersion
relations ($\omega_{max}$), the Grüneisen parameter, the volume of the unit
cell ($V$), and the scattering phase space ($SPS$) show a good correlation
with the thermal conductivity Juneja et al. (2019); Ju et al. (2020). These
are harmonic and quasi-harmonic properties but not all of them are equally
easy to calculate. While the maximum frequency is obtained with a single
harmonic calculation, evaluation of the Grüneisen parameter requires a few
calculations with different volumes. The $SPS$ is defined as:

$P_{3}^{\pm}(\mathbf{q}j)=\frac{1}{N_{\mathbf{q}}}\sum_{q1,q2,j1,j2}\delta(\omega_{q,j}\pm\omega_{q1,j1}-\omega_{q2,j2})\delta_{q\pm
q1,q2+G}\;,$ (3)

where $\omega$ and $\mathbf{q}$ are the phonon frequency and the crystalline
momentum of each phonon involved in the three-phonon interaction. Therefore,
although $SPS$ is a harmonic property, it requires exploration of all the
three-phonon scattering processes, which takes a some time.

The recent works have found good correlation between $w_{max}$ and $\kappa$
Juneja et al. (2019); Ju et al. (2020), and from our experience we know that
in several cases $\kappa$ shows a good correlation with $\ell$. The
correlation of $w_{max}$ and $\kappa$ can be broken in materials with high
frequency optical bands. Those bands are usually flat and therefore the
velocity is small and does no contribute significantly to the thermal
conductivity. For instance, lithium hydride (LiH) and lithium fluoride (LiF)
have $\omega_{max}\sim 32$ THz ($\sim 20$ THz) and $\kappa=22$ W/mK ($14$
W/mK) respectively Lindsay (2016) while silicon (Si) has $\omega_{max}\sim 16$
THz and $\kappa=145$ W/mK at 300 K Ward and Broido (2010).

In the effort of finding a new descriptor including the velocity to avoid the
effect of high frequency optical bands, we propose two candidates. The first
one is equivalent to the Landauer thermal conductance expressed in the
$\mathbf{q}$ space considering the transmission $T=1$ for all the modes:

$G_{max}=\frac{1}{N_{\mathbf{q}}V}\sum_{\mathbf{q},j}\frac{\hbar\omega}{2\pi
k_{B}T^{2}}\frac{e^{\hbar\omega/k_{B}T}}{(e^{\hbar\omega/k_{B}T}-1)^{2}}\textbf{v}\;,$
(4)

where $N$ and $V$ are number of q-points in the mesh sampling and the volume
of the unit cell, respectively. The second descriptor is a pseudo thermal
conductivity defined in the following. We know that the thermal conductivity
under the RTA can be expressed as:

$\kappa^{\alpha\beta}=\frac{1}{N_{\mathbf{q}}V}\sum_{\mathbf{q},j}C_{\mathbf{q},j}v_{\mathbf{q},j}^{\alpha}v_{\mathbf{q},j}^{\beta}\tau_{\mathbf{q},j}\;,$
(5)

where $\alpha,\beta$ are cartesian cordinates, $C_{\mathbf{q},j}$ the specific
heat and $\tau_{\mathbf{q},j}$ the relaxation time. The latter is calculated
from the phonon-phonon scattering rate as $\tau_{\mathbf{q},j}=\hbar/2\Gamma$,
where

$\displaystyle\Gamma_{\mathbf{q}j}=\frac{\pi\hbar}{2N_{\mathbf{q}j}}\sum_{\mathbf{q}^{\prime}j^{\prime}\mathbf{q}^{\prime\prime}j^{\prime\prime}}\frac{|\Phi(\mathbf{q}j\mathbf{q}^{\prime}j^{\prime}\mathbf{q}^{\prime\prime}j^{\prime\prime})|^{2}}{8\omega_{\mathbf{q}j}\omega_{\mathbf{q^{\prime}}j^{\prime}}\omega_{\mathbf{q^{\prime\prime}}j^{\prime\prime}}}\cdot$
(6)
$\displaystyle[(1+n_{\mathbf{q}^{\prime}j^{\prime}}+n_{\mathbf{q}^{\prime\prime}j^{\prime\prime}})\delta(\omega_{\mathbf{q}j}-\omega_{\mathbf{q}^{\prime}j^{\prime}}-\omega_{\mathbf{q}^{\prime\prime}j^{\prime\prime}})$
$\displaystyle+2(n_{\mathbf{q}^{\prime}j^{\prime}}-n_{\mathbf{q}^{\prime\prime}j^{\prime\prime}})\delta(\omega_{\mathbf{q}j}+\omega_{\mathbf{q}^{\prime}j^{\prime}}-\omega_{\mathbf{q}^{\prime\prime}j^{\prime\prime}})]\;,$

and $|\Phi(\mathbf{q}\mathbf{q}^{\prime}\mathbf{q}^{\prime\prime})|^{2}$ is
the three-phonon scattering matrix. If we approximate
$|\Phi(\mathbf{q}\mathbf{q}^{\prime}\mathbf{q}^{\prime\prime})|^{2}/8\omega_{\mathbf{q}j}\omega_{\mathbf{q^{\prime}}j^{\prime}}\omega_{\mathbf{q^{\prime\prime}}j^{\prime\prime}}\rightarrow
1/m^{2}$, then

$\displaystyle\Gamma_{\mathbf{q}}=\frac{\pi\hbar}{2m^{2}}(W^{-}_{\mathbf{q},j}+2W^{+}_{\mathbf{q},j})\;.$
(7)

The magnitude $W^{\pm}_{\mathbf{q},j}$ is the temperature dependent $SPS$
($TSPS$), which as well as the $SPS$, can be obtanied directly from Phono3py
Togo et al. (2015). With these definitions, we introduce the pseudo thermal
conductivity as:

$\kappa^{\alpha\beta}_{pseudo}=\frac{1}{N_{\mathbf{q}}V}\sum_{\mathbf{q},j}C_{\mathbf{q},j}v_{\mathbf{q},j}^{\alpha}v_{\mathbf{q},j}^{\beta}\left[\frac{\pi}{m^{2}}(W^{-}_{\mathbf{q},j}+2W^{+}_{\mathbf{q},j})\right]^{-1}\;.$
(8)

##  III RESULTS AND DISCUSSION

The values of $\ell$, $\kappa$ and all the descriptors described in the
previous section has been calculated for the 131 materials from first-
principles and with XenonPy.

The first-principles calculations show that the materials with stronger
intrinsic hydrodynamic behaviour (larger $\ell$) are those made of carbon,
such as diamond (mp-66), lonsdaleite (mp-47), and other carbon allotropes
(mp-611426 and mp-616440), as well as BAs in wurtzite(WZ) and zinc-blend(ZB)
structures. The former carbon based materials have a thermal conductivity of
around 2500 W/mK while the hydrodynamic scale ranges from 800-1000 nm. In
comparison with BAs, the $\kappa$ is smaller while $\ell$ is larger. More
explicity, WZ(ZB) BAs have a thermal conductivity of $\sim$1250
W/mK($\sim$1500 W/mK) and $\ell\sim$1100 nm($\sim$1400 nm). Another case is
BN, which has a large thermal conductivity (1200 W/mK) but smaller $\ell$ (319
nm) than BAs in wurtzite (WZ) and zinc-blend (ZB) phases. Those materials with
large $\ell$ are interesting from the viewpoint of studying the experimental
signatures of the hydrodynamic thermal transport, as they also exhibit a high
relative importance of N versus resistive scattering processes, as reported in
several works (Ward et al., 2009; Broido et al., 2013). Also, these materials
can be usefull to study the improvement of heat release in hot spots in
electronic devices due to fast operation rates, as hydrodynamic effects near
to the edges of electronic circuits create a vorticity that helps to increase
the heat dissipation due to an effective change of the thermal conductivity
Torres et al. (2018).

Table 1: Materials with $\ell>\gamma\kappa$. $\kappa$ is expressed in W/mK and $\ell$ in nm. More materials fulfilling this condition can be found in the APPENDIX. Material | $\mathbf{\kappa}$ | $\mathbf{\ell}$ |  | Material | $\mathbf{\kappa}$ | $\mathbf{\ell}$
---|---|---|---|---|---|---
AlP(mp-1550) | 101.5 | 151.7 |  | GaN(mp-830) | 201.4 | 202.3
BeS(mp-422) | 164.4 | 173.2 |  | BeTe(mp-252) | 200 | 351.7
Bi2Te3 | 1.5 | 10 |  | NaBr(mp-22916) | 3.4 | 11.6
(mp-1227339) |  |  |  |  |  |
KBr(mp-23251) | 3.3 | 17.1 |  | BaSe(mp-1253) | 11.5 | 53.2
BaTe(mp-1000) | 13.4 | 63.8 |  | AgI(mp-22925) | 1.2 | 42.6
CdS(mp-2469) | 25.6 | 111.3 |  | CdSe(mp-2691) | 8.2 | 59.4
CdTe(mp-406) | 6 | 42.6 |  | CuBr(mp-22913) | 1 | 11.2
GaAs(mp-2534) | 36.5 | 94.3 |  | GaP(mp-2490) | 84.7 | 180.0
GaSb(mp-1156) | 40.8 | 145.4 |  | InAs | 24.6 | 97.0
|  |  |  | (mp-20305) |  |
InP(mp-20351) | 77.0 | 426.2 |  | InSb | 15.6 | 77.8
|  |  |  | (mp-20012) |  |
MgTe(mp-13033) | 15.8 | 88.8 |  | ZnS | 23.4 | 75.6
|  |  |  | (mp-10695) |  |
ZnSe(mp-1190) | 17.5 | 52.4 |  | ZnTe(mp-2176) | 15.9 | 96.2
AgI(mp-22894) | 1 | 18.6 |  | AlAs(mp-8881) | 82.3 | 198.6
AlP(mp-8880) | 83.1 | 126.5 |  | AlSb | 84.7 | 302
|  |  |  | (mp-1018100) |  |
CdS(mp-672) | 19.4 | 82.1 |  | CdSe(mp-1070) | 6.8 | 42.8
CdTe(mp-12779) | 4.6 | 42.4 |  | CuI(mp-569346) | 6.2 | 46.1
GaP(mp-8882) | 79.2 | 190.3 |  | GaSb | 26.8 | 118
|  |  |  | (mp-1018059) |  |
InAs(mp-1007652) | 18.7 | 93.6 |  | InP | 67.6 | 314
|  |  |  | (mp-966800) |  |
InSb(mp-1007661) | 10.8 | 72 |  | MgTe(mp-1039) | 11.8 | 75.3
ZnTe(mp-8884) | 13.1 | 82.6 |  | Sn(mp-117) | 28.3 | 108.8
TlMgLa | 9.3 | 54.2 |  | InN(mp-22205) | 110.1 | 226.3
(mp-962063) |  |  |  |  |  |
AlSb | 131.5 | 417.3 |  | InN | 108.7 | 235.4
(mp-2624) |  |  |  | (mp-20411) |  |
AlAs(mp-2172) | 114.7 | 237 |  |  |  |
Be2C(mp-1569) | 149.6 | 213 |  | B2AsP | 433 | 576
|  |  |  | (mp-1008528) |  |

The list of 131 materials used in this work corresponds to some of the most
studied and useful insulators and seminconductors for thermal management and
thermoelectrics Muñoz Rojo et al. (2016, 2017); Koga et al. (1998); Maurya et
al. (2019); Hoshino and Shimoji (1970); Xiao and Zhao (2018). From the first-
principles calculations it can be observed that most of them fulfil the
condition $\ell>\gamma\kappa$, where $\gamma=1$ nm$\cdot$mK/W (Table 1).
$\gamma$ is just a factor to evaluate the difference in the order of magnitude
and units between $\kappa$ and $\ell$. Taking into account this observation,
from the data of the studied 131 materials, we use machine learning models to
look for materials with some interesting hydrodynamic properties, such as
materials with $\ell$ values higher than 1000 nm (similar to C-based materials
or BAs) or materials fulfilling the condition $\ell>\gamma\kappa$.

Recently various machine learning techniques have been used to predict the
thermal conductivity, like Bayesian optimization Seko et al. (2015), Gaussian
process regression Juneja et al. (2019) or transfer learning Ju et al. (2020).
In this work we have examined different machine learning models trained using
the scikit-learn library in Python. Three types of ensemble learning with
regression tree models and neural networks with different combination of
layers were considered. The corresponding packages are ExtraTreesRegressor
Geurts et al. (2006), RandomForestRegressor Breiman (2001), and
GradientBoostingRegressor Friedman (2002) and XenonPy. Hyperparameters for
each type of model were selected using 5-fold cross-validation based on the
GridSearchCV package, using in all cases 25$\%$ of data for test. We selected
negative mean squared error as the scoring metric. Finally, the model with the
lowest root mean squared error (RMSE) on the test data was selected for the
high-throughput screening of $\ell$ and $\kappa$.

Figure 1: Scatter plot of relevant correlations found in the present work for
the 131 materials set. a) $SPS$ as function of $\omega_{max}$. b) $\kappa$ as
function of $\omega_{max}$. c) $\ell$ as function of $\kappa$. d) $\ell$ as
function of Min. Mendeleev number. For $\ell$ and $\kappa$ the median value of
the three components have been considered.

Considering all the descriptors described in the previous section, those
calculated from first-principles and those obtained directly from XenonPy, we
have calculated the Pearson’s correlation matrix and obtained the
corresponding absolute correlation coefficients. By examining the correlation
matrix for all the 131 materials we can observe a few interesting things.
First of all, we see that the $SPS$ have a correlation of $0.996$ with
$w_{max}$ (Fig. 1a). From the harmonic point of view, we have also confirmed
that $w_{max}$ strongly correlates with $\kappa$ (Fig. 1b), with a correlation
coefficient of $0.86$, while $v$ and $G_{max}$ with values of $0.78$ and
$0.71$ respectively are better than $V$ ($0.53$). Those correlations for the
hydrodynamic scale $\ell$ are worse. In the case of $w_{max}$ is $0.59$, while
for $v$ and $G_{max}$ are $0.48$ and $0.33$ respectively. Unfortunately the
correlation coefficient between $\ell$ and the pseudo thermal conductivity is
lower than 0.5 (0.48). On the other hand, the correlation between $\kappa$ and
$\ell$ is quite high ($0.88$) (Fig. 1c), as both parameters have a strong
correlation with the phonon relaxation times according to their definitions
Torres (2018); Torres et al. (2018). We also observe that $\ell$ strongly
correlates with the total relaxation time $\tau$ ($0.92$) and the Normal
relaxation time $\tau_{N}$ ($0.90$), as by definition
$\ell=\ell(\tau,\tau_{N})$ Torres (2018). For $\kappa$ these correlations are
weaker, $0.72$ and $0.73$ respectively. In addition to the Pearson’s
correlation coefficients, which accounts for linear correlations, the Maximal
Information Coefficient (MIC) has been also calculated, which also consider
non-linear correlations.

Table 2: Top descriptors with MIC correlation with $\ell$ higher than $0.7$. Descriptor |  Correlation with $\mathbf{\ell}$
---|---
Min. Mendeleev number | 0.79
Ave. dipole polarizability | 0.76
Max. dipole polarizability | 0.76
Max. polarizability | 0.72
Min. Allen energy | 0.71
Var. Mendeleev number | 0.7
Min. atomic number | 0.7
Min. atomic weight | 0.7
Max. Van Der Waals radius (Alvarez) | 0.7

The descriptors with a MIC with $\ell$ higher than $0.7$ are summarized in
Table 2. Regarding those descriptors we can observe that none of the harmonic
or quasi-harmonic descriptors are in the list, and that the highest
correlation with $\ell$ is found for the minimum Mendeleev number ($0.79$)
(Fig. 1d). While the meaning of this compositional descriptors is self-
explained by the names, the Min, Max corresponds to Min/Max values of the
descriptor for each unique atom of the compound, and the Ave to the average of
the descriptor of all the atoms of the compound. This information provides
valuable physical insight that is supported by previous works. As detailed
before, the hydrodynamic thermal transport is strongly related to the phonon
relaxation times and to several compositional descriptors. Several works have
determined analytical expressions for the phonon relaxation times Han and
Klemens (1993); Han (1996); Herring (1954); Tamura (1983); Klemens (1955,
1960); Ward and Broido (2010), showing that those magnitudes are directly
related with the material density $\rho$, group velocity, and frequency for
the case of Normal and Umklapp scattering, in agreement with some of the
higher correlated descriptors for $\kappa$. In the case of impurity/mass
defect scattering the analytic expresions show the relation to the atomic
masses of the atoms of a compound Inyushkin et al. (2003, 2004). This also
agrees with the high correlation of $\kappa$ and $\ell$ with the atomic number
and atomic weight. Our calculations have also identified the polarizability as
a key factor to determine the intrinsic hydrodynamic effects. A recent work
have shown that by strain engineering it is possible to tune the
polarizability of a material and modify the thermal conductivity Torres et al.
(2019), revealing the high correlation between both magnitudes. In our case
the correlation of the Max. polarizability with $\kappa$ is 0.63, while for
$\ell$ is 0.76. Therefore it is expected that modifications in the
polarizability will have a big influence on $\ell$.

The descriptors for $\kappa$ showing MIC higher than $0.7$ can be found in
Table 3. Note that here we have omitted the $SPS$, showing also a correlation
of $0.7$, as for our calculated materials it is almost equivalent to
$\omega_{max}$ but it is much more expensive computationally (Fig. 1a). In
materials with large band gaps this might not be true.

Table 3: Top descriptors with MIC correlation with $\kappa$ higher than $0.7$. Descriptor |  Correlation with $\mathbf{\kappa}$
---|---
$w_{max}$ | 0.86
Min. atomic weight | 0.79
Min. atomic number | 0.79
$v$ | 0.78
Ave. heat of formation | 0.77
Ave. GS energy | 0.71
$G_{max}$ | 0.71

By using the descriptors from Table 2 and Table 3, we first create machine
learning models for $\ell$ and $\kappa$ with our known 131 materials
calculated from first-principles.

A neural network with 6 layers (with 9-69-59-30-13-1 neurons in each layer)
has been selected for $\ell$, and a neural network with 5 layers (with
7-59-52-17-1 neurons in each layer) for $\kappa$ as the models with lower
RMSE. The best models showed an RMSE value of 2.64 nm for $\ell$ and 5.47 W/mK
for $\kappa$ respectively. In Fig. 2 the first-principles values of $\kappa$
(top) and $\ell$ (bottom) are represented versus the prediction obtained from
informatics techniques. It can be observed that the predictions are slightly
better for high $\kappa$ values than for the lower ones. Low $\kappa$
materials are highly anharmonic and the present formalism might not be
accurate enough. The straight line with slope of 1 represents an exact
prediction. The calculated values of $\ell$ and $\kappa$ are averaged over the
three axial components.

Now, we proceed to screen all the materials available in the phonon database
of the Kyoto University Jain et al. (2013); Ong et al. (2013, 2015); Togo
(2020). In the phonon database we can find the harmonic information for 10035
materials. From the initial list, we have removed those showing negative
phonon bands, reducing the number to 5059. This harmonic information allows us
to compute harmonic and quasi-harmonic properties such as $\omega_{max}$,
$G_{max}$, $C$, $v$, $SPS$ and $\kappa_{pseudo}^{\alpha,\beta}$, in addition
to the material density $\rho$ and volume $V$, for all the 5059 materials.

Figure 2: Values of $\kappa$ and $\ell$ calculated by first-principles (FP) versus the prediction obtained by machine learning algorithms for the 25$\%$ test set. The straight line represents and ideal exact prediction. Table 4: Validation set of materials with $\ell>\gamma\kappa$. $\kappa$ is expressed in W/mK and $\ell$ in nm. Material | $\mathbf{\kappa_{x}}$ | $\mathbf{\kappa_{y}}$ | $\mathbf{\kappa_{z}}$ | $\mathbf{\ell_{x}}$ | $\mathbf{\ell_{y}}$ | $\mathbf{\ell_{z}}$
---|---|---|---|---|---|---
GeAs2(mp-17524) | 8 | 2.6 | 2.4 | 56 | 102.9 | 82.9
SiP(mp-2798) | 41.9 | 19.4 | 12.9 | 93.5 | 41 | 133.9
BeP2(mp-27148) | 73.5 | 71.5 | 57.5 | 92.9 | 91.4 | 105
GeO2(mp-470) | 63.6 | 63.6 | 101 | 93.8 | 93.8 | 95.5

After screening over the 5059 materials, the models have found 147 compounds
fulfilling the $\ell>\gamma\kappa$ criteria. From the whole list we have
selected 5 materials to check the validity of the predictions from first-
principles: GeAs2, SiP, BeP2 and GeO2. The values of thermal conductivity and
hydrodynamic scale calculated from first-principles for each axial component
can be found in Table 4. Due to the limited resources, all the 147 materials
can not be calculated from first-principles. The full list of materials
fulfilling the condition $\ell>\gamma\kappa$ from the phonon database
according to our calculations can be found in the APPENDIX.

From Table 4 it can be observed that all the materials fullfil the condition
$\ell>\gamma\kappa$ also from first-principles calculations. It is expected
that, according to this condition, the materials could be useful for thermal
management and/or thermoelectrics. In several works it has been analysed that
the dominance of Normal phonon-phonon and electron-phonon scattering in front
of Umklapp processes produces a fluid-like behaviour (hydrodynamics) that
generates phonon drag and increases the Seebeck coefficient, and consequently
an enhancement of the thermoelectric efficiency (ZT) Zhou et al. (2015);
Kuleev et al. (2003); Levchenko and Schmalian (2020); Gooth et al. (2018);
Kuleyev et al. (2016); Huang and Lucas (2021); Glazov (2019). In the case of
SiP, E. Jiang et. al Jiang et al. (2019) showed that in its monolayer form can
achieve a ZT$\sim$0.8 and a Seebeck coefficient of 3 mV/K, opening the door to
further improvements for nanoscale P-based materials. In the case of BeP2, its
bulk ZT is around 1, and in addition several works have highlighted the
greatly tunable electronic properties that exhibit in its monolayer and penta-
monolayer form Azarmi et al. (2019); Kansara et al. (2020). GeAs2 has been
reported to be a material with ultralow thermal conductivity and high
thermoelectric efficiency (Zhao et al., 2017), showing a Seebeck coefficent
higher than 0.3 mV/K with high carrier concentration. Finally, GeO2 has been
examined in the last decades as material used to improve the thermoelectric
efficiency of other materials by alloying with Cu, SiO2 and CdO to create thin
films Lucy (2006); Snedaker (2015); Narushima et al. (2000) or coating
conducting polymers with GeO2 nanoparticles Shiraishi et al. (2017).
Unfortunately, the models have not found any new material with high $\ell$
value above 1000 nm.

Figure 3: Values of $\kappa$ and $\ell$ for the materials corresponding to
Table 1 (blue circles), Table 4 (red squares) and C and BN allotropes (green
diamonds). The insets show the predicted VS the first-principles (FP)
calculations for $\kappa$ (top) and $\ell$ (bottom).

The results from Table 1 and Table 4 are represented in Fig. 3, which seems to
indicate that thermoelectric materials tend to exhibit a $\ell>\gamma\kappa$
behaviour. The insets show, for the materials of Table 4, the predicted values
by the machine learning models as function of their values calculated from
first-principles. The straight line represents the prefect prediction. It can
be observed that for low values of $\kappa$ the prediction is worse than for
high values. This might be because the present formalism used to calculate the
thermal conductivity is not suitable to capture the high anharmonicity of low
$\kappa$ materials, affecting the values used to train the models and
therefore the predictions. In addition, materials with $\ell>$ 1000 nm and
below the $\ell=\gamma\kappa$ line would be the best candidates to improve
heat dissipation in thermal management of electronics.

##  IV CONCLUSIONS

In conclusion, we have used XenonPy’s structural and compositional material
descriptors together with first-principles calculations to explore the
hydrodynamic thermal transport. By using machine learning calculations, we
have found that the structural descriptors are not relevant to study the
hydrodynamic thermal transport, while there are 9 compositional descriptors
with a correlation higher than $0.7$ with $\ell$, the key parameter expressing
the scale of phonon hydrodynamics. The calculations have also revealed that
diamond and BAs and their allotropes have both high $\ell$. Those materials
are the best candidates to identify hydrodynamic signatures from an
experimental point of view. The second group of materials that have been
considered are those having $\ell>\gamma\kappa$, which seems to be a good
indicator to explore new materials to be used for thermal management and
thermoelectricity. Here several materials have been found in the original data
set, while the machine learning techniques identified 147 materials. From the
identified samples, BeP2, SiP, GeAs2 and GeO2 have been calculated from first-
principles and used to validate the models, which according to recent works
seem to have high thermoelectric performance, in agreement with our original
guess. Although these materials have been already reported for this propose,
the model found them without knowing of their existence a priori. Finally, the
full list of values of $\kappa$ and $\ell$ of the materials from the original
dataset and the 147 found materials with $\ell>\gamma\kappa$ detailed in the
APPENDIX can be useful to look for new compounds for thermal management and
thermoelectrics.

## ACKNOWLEDGMENTS

We acknowledge financial support of The Canon Foundation. This work was
partially supported by JSPS KAKENHI (19H00744, 19H01132, 19H50820, 18K18017,
19K14902) and JST CREST (JPMJCR20Q3, JPMJCR19I2, JPMJCR19I1).

## APPENDIX

In the next table all the values of thermal conductivity $\kappa$ and non-
local length $\ell$ for the materials used for the training are displayed. The
values that have been already showed in the main text are not included.

Table 5: Values of $\kappa$ and $\ell$ of the full data set used for training
the models. The averaged values of the three components are considered.

Material |  $\mathbf{\kappa}$ (W/mK) |  $\mathbf{\ell}$ (nm)
---|---|---
Ge(mp-32) | 58.6 | 116.0
C(mp-66) | 2552.3 | 880.2
LiCl(mp-22905) | 5.1 | 4.9
LiBr(mp-23259) | 3.5 | 5.0
LiI(mp-22899) | 3.3 | 4.6
NaF(mp-682) | 50.9 | 26.7
NaCl(mp-22862) | 9.7 | 18.3
NaI(mp-23268) | 1.7 | 9.1
KF(mp-463) | 7.9 | 15.8
KCl(mp-23193) | 9.5 | 26.1
KI(mp-22898) | 1.6 | 12.3
RbF(mp-11718) | 2.8 | 9.6
RbF(mp-11718) | 2.8 | 9.6
RbCl(mp-23295) | 2.2 | 12.9
RbBr(mp-22867) | 5.3 | 28.8
RbI(mp-22903) | 2.5 | 20.3
AgCl(mp-22922) | 0.7 | 0.9
AgBr(mp-23231) | 0.4 | 1.2
CsF(mp-1784) | 2.4 | 8.2
MgO(mp-1265) | 57.4 | 35.2
CaO(mp-2605) | 24.7 | 22.8
CaS(mp-1672) | 38.0 | 57.6
CaSe(mp-1415) | 23.1 | 46.7
CaTe(mp-1519) | 12.4 | 41.6
SrO(mp-2472) | 10.3 | 13.8
BaO(mp-1342) | 5.2 | 5.5

BaS(mp-1500) | 4.8 | 21.4
---|---|---
CdO(mp-1132) | 11.1 | 8.9
PbS(mp-21276) | 3.7 | 2.8
PbSe(mp-2201) | 1.3 | 2.8
PbTe(mp-19717) | 1.2 | 4.3
AlN(mp-1700) | 267.2 | 177.9
BAs(mp-10044) | 1537.8 | 1389.0
BeO(mp-1778) | 527.5 | 194.2
BeSe(mp-1541) | 239.6 | 255.2
BeTe(mp-252) | 200.0 | 351.7
BN(mp-1639) | 1192.7 | 319.1
BP(mp-1479) | 449.8 | 301.3
CuCl(mp-22914) | 1.3 | 7.7
CuI(mp-22895) | 6.1 | 39.8
GaN(mp-830) | 201.4 | 202.4
SiC(mp-8062) | 471.5 | 317.4
ZnO(mp-1986) | 49.9 | 70.6
AlAs(mp-8881) | 83.0 | 199.0
CdO(mp-1132) | 11.1 | 8.9
PbS(mp-21276) | 3.7 | 2.8
PbSe(mp-2201) | 1.3 | 2.8
PbTe(mp-19717) | 1.2 | 4.3
AlN(mp-1700) | 267.2 | 177.9
BAs(mp-10044) | 1537.8 | 1389.0
BeO(mp-1778) | 527.5 | 194.2
BeSe(mp-1541) | 239.6 | 255.2
BeTe(mp-252) | 200.0 | 351.7
BN(mp-1639) | 1192.7 | 319.1
BP(mp-1479) | 449.8 | 301.3
CuCl(mp-22914) | 1.3 | 7.7
CuI(mp-22895) | 6.1 | 39.8
GaN(mp-830) | 201.4 | 202.4
SiC(mp-8062) | 471.5 | 317.4
ZnO(mp-1986) | 49.9 | 70.6
AlAs(mp-8881) | 83.0 | 199.0
AlN(mp-661) | 243.4 | 173.2
BAs(mp-984718) | 1169.7 | 1043.6
BeO(mp-2542) | 392.9 | 165.5
BN(mp-2653) | 965.8 | 277.1
BP(mp-1008559) | 355.5 | 240.3
CdS(mp-672) | 19.8 | 81.6
CdTe(mp-12779) | 4.8 | 42.5
CuCl(mp-1184046) | 1.1 | 7.8
GaAs(mp-8883) | 31.2 | 91.0
GaN(mp-804) | 204.4 | 222.0
GaP(mp-8882) | 80.7 | 190.8
SiC(mp-7140) | 389.0 | 293.3
ZnO(mp-2133) | 43.7 | 63.5
ZnS(mp-560588) | 23.0 | 76.9
ZnSe(mp-380) | 16.1 | 57.9
AgSbBa(mp-984720) | 0.4 | 2.7
CoAsZr(mp-961689) | 36.5 | 42.8
CoSbZr(mp-22377) | 39.8 | 44.9
NiScSb(mp-3432) | 4.7 | 6.8
BeO(mp-1778) | 527.5 | 194.2
BeSe(mp-1541) | 239.6 | 255.2
BeTe(mp-252) | 200.0 | 351.7
BN(mp-1639) | 1192.7 | 319.1
BP(mp-1479) | 449.8 | 301.3

CuCl(mp-22914) | 1.3 | 7.7
---|---|---
CuI(mp-22895) | 6.1 | 39.8
GaN(mp-830) | 201.4 | 202.4
SiC(mp-8062) | 471.5 | 317.4
ZnO(mp-1986) | 49.9 | 70.6
AlAs(mp-8881) | 83.0 | 199.0
CdO(mp-1132) | 11.1 | 8.9
PbS(mp-21276) | 3.7 | 2.8
NiSnTi(mp-924130) | 7.7 | 12.0
NiSnTi(mp-623646) | 7.7 | 12.0
NiSnTi(mp-22782) | 7.7 | 12.0
NiSnZr(mp-924129) | 28.1 | 39.9
NiSnZr(mp-30806) | 28.1 | 39.9
PdSrTe(mp-961663) | 0.5 | 3.2
PtGaTa(mp-961677) | 47.1 | 84.2
PtGeTi(mp-1008680) | 7.8 | 11.7
PtGeTi(mp-961671) | 7.8 | 11.7
PtNbIn(mp-961708) | 18.6 | 49.4
RhNbSi(mp-1100400) | 24.2 | 26.5
RhSbHf(mp-10367) | 29.8 | 52.3
RuAsV(mp-1100404) | 11.6 | 11.5
SiSrCd(mp-962076) | 4.4 | 7.7
SnSrBa(mp-962062) | 2.4 | 14.2
GeC(mp-1002164) | 776.3 | 376.9
Lonsdaleite(mp-47) | 2159.6 | 821.4
BC2N(mp-30148) | 1112.1 | 679.2
SiSn(mp-1009813) | 151.1 | 317.0
BSb(mp-997618) | 379.5 | 494.7
B2AsP(mp-1008528) | 433.5 | 576.6
C3N4(mp-571653) | 334.4 | 387.9
BC2N(mp-629458) | 993.8 | 581.5
BeCN2(mp-15703) | 416.3 | 207.2
C3N4(mp-2852) | 300.6 | 124.1
SnC(mp-1009820) | 261.7 | 308.1
B(mp-160) | 62.5 | 44.3
LiBC(mp-9244) | 130.0 | 39.6

The full list of materials identified by the machine learning models
fulfilling the condition $\ell>\gamma\kappa$ are:

Se(mp-147), Sm2O3(mp-218), YbSe(mp-(286), B2O3(mp-306), MgP4(mp-384),
Sc2S3(mp-401), CdP2(mp-402), AlF3(mp-468),GeO2(mp-470), NiP2(mp-486),
Ga2S3(mp-539), TeO2(mp-557), GaF3(mp-588), B6As(mp-624), SnSe2(mp-665),
GeSe(mp-700), P2Pt(mp-730), Li3P(mp-736), Li3As(mp-757), Ga2O3(mp-886),
As2Se3(mp-909), CuP2(mp-927), GeTe(mp-938), NbZn3(mp-953), K2O(mp-971),
SrF2(mp-981), BaF2(mp-1029), SrS(mp-1087), Al2O3(mp-1143), ThS2(mp-1146),
YbO(mp-1216), MgS(mp-1315), Ga2Se3(mp-1340), ZnP2(mp-1392), Rb2O(mp-1394),
Sn2S3(mp-1509), As2O3(mp-1581), Al4C3(mp-1591), Tm2O3(mp-1767), YbTe(mp-1779),
SiF4(mp-1818), YbS(mp-1820), SiAs(mp-1863), ZnF2(mp-1873), SrTe(mp-1958),
Mg3As2(mp-1990), Li3Sb(mp-2074), Rb2P3(mp-2079), Si5Ir3(mp-2084),
KGe(mp-2146), Er2S3(mp-2234), GeS(mp-2242), Tm2S3(mp-2309), Dy2O3(mp-2345),
Ga2Te5(mp-2371), NaSi(mp-2402), P4Se5(mp-2447), Al2S3(mp-2654), CaF2(mp-2741),
SrSe(mp-2758), SiP(mp-2798), Lu2S3(mp-2826), InF3(mp-6949), Be3N2(mp-6977),
Rb2Se3(mp-7447), GeF2(mp-7595), K2Se3(mp-7670), ThSe3(mp-7951), PF5(mp-8511),
AsF5(mp-8723), Al2Te5(mp-9254), GeAs(mp-9548), ZnS(mp-9946), MgSe(mp-10760),
RbSn(mp-11054), Sr3As4(mp-15339), Ge3Os2(mp-16610), GeAs2(mp-17524),
InS(mp-19795), InSe(mp-20485), Si3Ru2(mp-22192), In2O3(mp-22323),
CsCl(mp-22865), YbBr2(mp-22882), XeF2(mp-22885), CaBr2(mp-22888),
ZnCl2(mp-22889), CaCl2(mp-22904), CsBr(mp-22906), PtBr3(mp-23165),
IO2(mp-23170), SrI2(mp-23181), BaCl2(mp-23199), InI(mp-23202), MgI2(mp-23205),
SrCl2(mp-23209), MgCl2(mp-23210), AsI3(mp-23218), YbCl2(mp-23220),
BaI2(mp-23260), BeCl2(mp-23267), BeP2(mp-27148), ZnI2(mp-27161),
BaBr2(mp-27456), Te3Cl2(mp-27628), ThI4(mp-27697), PdI2(mp-27747),
GeI2(mp-27922), CdI2(mp-28248), P2Pd(mp-28266), PtI2(mp-28319),
Pt3I8(mp-28320), GaBr2(mp-28384), Br2O(mp-28460), In7Cl9(mp-28730),
Ca5P8(mp-28879), Cl2O(mp-29537), NaGe(mp-29657), Tl2Te3(mp-29711),
LaAs2(mp-29815), CaI2(mp-30031), MgBr2(mp-30034), BeI2(mp-30140),
GaCl3(mp-30952), Si9Te8(mp-31135), NbI5(mp-31487), Tm2Se3(mp-32850),
AuBr(mp-505366), Cs2S(mp-540703), KSn(mp-542374), ThCl4(mp-567431),
SrBr2(mp-567744), HfI4(mp-569059), SiB3(mp-569128), TeI(mp-569766),
ThBr4(mp-570229), YbI2(mp-570418), SeBr(mp-570589), Ga3Os(mp-570844),
ZrI4(mp-571235), Ga2I3(mp-636675), In7Te${}_{1}0$(mp-669311),
Lu2Se3(mp-673650), Sc2Se3(mp-684690), CsSi(mp-866482), SiGe(mp-978534).

## References

  * Guyer and Krumhansl (1966a) R. A. Guyer and J. A. Krumhansl, Phys. Rev. 148, 766 (1966a).
  * Esfarjani and Stokes (2008) K. Esfarjani and H. Stokes, Phys. Rev. B 77, 144112 (2008), ISSN 1098-0121, URL http://link.aps.org/doi/10.1103/PhysRevB.77.144112.
  * Chaput (2013) L. Chaput, Phys. Rev. Lett. 110, 265506 (2013), eprint 1303.4062.
  * Torres (2018) P. Torres, _Thermal Transport in Semiconductors: First Principles and Phonon Hydrodynamics_ (Springer Nature, Switzerland, 2018).
  * Johnson et al. (2013) J. A. Johnson, A. A. Maznev, J. Cuffe, J. K. Eliason, A. J. Minnich, T. Kehoe, C. M. S. Torres, G. Chen, and K. A. Nelson, Phys. Rev. Lett. 110, 025901 (2013), ISSN 00319007, eprint 1204.4735, URL http://link.aps.org/doi/10.1103/PhysRevLett.110.025901.
  * Siemens et al. (2010) M. E. Siemens, Q. Li, R. Yang, K. K. A. Nelson, E. H. Anderson, M. M. Murnane, and H. C. Kapteyn, Nat. Mater. 9, 26 (2010), ISSN 1476-1122, URL http://www.nature.com/nmat/journal/v9/n1/abs/nmat2568.html.
  * Minnich et al. (2011) A. J. Minnich, J. A. Johnson, A. J. Schmidt, K. Esfarjani, M. S. Dresselhaus, K. A. Nelson, and G. Chen, Phys. Rev. Lett. 107 (2011), ISSN 00319007, URL {%}3CGoto.
  * Regner et al. (2013) K. T. Regner, D. P. Sellan, Z. Su, C. H. Amon, A. J. H. McGaughey, and J. Malen, Nature Comm. 4, 1640 (2013), ISSN 2041-1723, URL http://www.ncbi.nlm.nih.gov/pubmed/23535661.
  * Wilson and Cahill (2014) R. B. Wilson and D. G. Cahill, Nature Comm. 5, 5075 (2014), ISSN 2041-1723, URL http://www.ncbi.nlm.nih.gov/pubmed/25270394.
  * Hu et al. (2015) Y. Hu, L. Zeng, A. J. Minnich, M. S. Dresselhaus, and G. Chen, Nat. Nanotechnol. 10, 701 (2015), ISSN 1748-3387, URL http://www.scopus.com/inward/record.url?eid=2-s2.0-84930350912{&}partnerID=tZOtx3y1{%}5Cnhttp://www.nature.com/doifinder/10.1038/nnano.2015.109.
  * Hoogeboom-Pot et al. (2015) K. M. Hoogeboom-Pot, J. N. Hernandez-Charpak, X. Gu, T. D. Frazer, E. H. Anderson, W. Chao, R. W. Falcone, R. Yang, M. M. Murnane, H. C. Kapteyn, et al., Proc. Natl. Acad. Sci. 112, 201503449 (2015), ISSN 0027-8424, URL http://www.pnas.org/lookup/doi/10.1073/pnas.1503449112.
  * Vermeersch et al. (2014) B. Vermeersch, A. M. S. Mohammed, G. Pernot, Y. R. Koh, and A. Shakouri, Phys. Rev. B 90, 1 (2014), ISSN 1550235X.
  * Vermeersch and Shakouri (2014) B. Vermeersch and A. Shakouri, arXiv:1412.6555 pp. 1–21 (2014), eprint 1412.6555, URL http://arxiv.org/abs/1412.6555.
  * Cepellotti et al. (2015) A. Cepellotti, G. Fugallo, L. Paulatto, M. Lazzeri, F. Mauri, and N. Marzari, Nat. Commun. 6, 6400 (2015).
  * Lee et al. (2015) S. Lee, D. Broido, K. Esfarjani, and G. Chen, Nat. Commun. 6, 6290 (2015).
  * Torres et al. (2018) P. Torres, A. Ziabari, A. Torelló, J. Bafaluy, J. Camacho, X. Cartoixà, A. Shakouri, and F. X. Alvarez, Phys. Rev. Materials 2, 076001 (2018).
  * Guyer and Krumhansl (1966b) R. A. Guyer and J. A. Krumhansl, Phys. Rev. 148, 778 (1966b).
  * Ziabari et al. (2018) A. Ziabari, P. Torres, B. Vermeersch, X. Cartoixà, A. Torelló, Y. Xuan, J.-H. Bahk, Y. Koh, Y. P., F. Alvarez, et al., Nat. Commun. 9, 255 (2018).
  * R (2020) Y. R, _XenonPy_ (2020), URL https://github.com/yohida-lab/XenonPy.
  * Jain et al. (2013) A. Jain, Y. J. Yu, and A. J. H. McGaughey, Phys. Rev. B 87, 1 (2013), ISSN 10980121.
  * Ong et al. (2015) S. P. Ong, S. Cholia, A. Jain, M. Brafman, D. Gunter, G. Ceder, and K. A. Persson, Computational Materials Science 97, 209 (2015), ISSN 0927-0256, URL http://www.sciencedirect.com/science/article/pii/S0927025614007113.
  * Seko et al. (2015) A. Seko, A. Togo, H. Hayashi, K. Tsuda, L. Chaput, and I. Tanaka, Phys. Rev. Lett. 115, 205901 (2015), URL https://link.aps.org/doi/10.1103/PhysRevLett.115.205901.
  * Carrete et al. (2014) J. Carrete, W. Li, N. Mingo, S. Wang, and S. Curtarolo, Phys. Rev. X 4, 011019 (2014), URL https://link.aps.org/doi/10.1103/PhysRevX.4.011019.
  * Togo and Tanaka (2015) A. Togo and I. Tanaka, Scr. Mater. 108, 1 (2015).
  * Togo et al. (2015) A. Togo, L. Chaput, and I. Tanaka, Phys. Rev. B 91, 094306 (2015).
  * Torres et al. (2020) P. Torres, F. Alvarez, and X. Cartoixa, _Kinetic Collective Model: BTE-based hydrodynamic model for thermal transport_ (2020), URL https://physta.github.io/.
  * Juneja et al. (2019) R. Juneja, G. Yumnam, S. Satsangi, and A. K. Singh, Chem. Mater. 31, 5145 (2019).
  * Ju et al. (2020) S. Ju, R. Yoshida, C. Liu, K. Hongo, T. Tadano, and J. Shiomi (2020), URL https://arxiv.org/abs/1909.11234.
  * Lindsay (2016) L. Lindsay, Phys. Rev. B 94, 174304 (2016), URL https://link.aps.org/doi/10.1103/PhysRevB.94.174304.
  * Ward and Broido (2010) A. Ward and D. A. Broido, Phys. Rev. B 81, 085205 (2010), ISSN 1098-0121, eprint arXiv:1107.5288v1, URL http://link.aps.org/doi/10.1103/PhysRevB.81.085205.
  * Ward et al. (2009) A. Ward, D. A. Broido, D. A. Stewart, and G. Deinzer, Phys. Rev. B 80, 125203 (2009).
  * Broido et al. (2013) D. A. Broido, L. Lindsay, and T. L. Reinecke, Phys. Rev. B 88, 214303 (2013).
  * Muñoz Rojo et al. (2016) M. Muñoz Rojo, Y. Zhang, C. V. Manzano, R. Alvaro, J. Gooth, M. Salmeron, and M. Martin-Gonzalez, Sci. Rep. 6, 19014 (2016), ISSN 2045-2322, URL http://www.nature.com/articles/srep19014.
  * Muñoz Rojo et al. (2017) M. Muñoz Rojo, B. Abad, C. V. Manzano, P. Torres, X. Cartoixà, F. X. Alvarez, and M. S. Martín-González, Nanoscale 9, 6741 (2017), ISSN 2040-3364, URL http://pubs.rsc.org/en/Content/ArticleLanding/2017/NR/C7NR02173A.
  * Koga et al. (1998) T. Koga, X. Sun, S. Cronin, and M. Dresselhausb, Applied Physics Letters 73 (1998).
  * Maurya et al. (2019) V. Maurya, U. Paliwal, G. Sharma, and K. B. Joshi, RSC Adv. 9, 13515 (2019), URL http://dx.doi.org/10.1039/C9RA01573F.
  * Hoshino and Shimoji (1970) H. Hoshino and M. Shimoji, Journal of Physics and Chemistry of Solids 31, 1553 (1970), ISSN 0022-3697, URL https://www.sciencedirect.com/science/article/pii/0022369770900405.
  * Xiao and Zhao (2018) Y. Xiao and L.-D. Zhao, npj Quantum Materials 3, 2397 (2018), ISSN 1.
  * Geurts et al. (2006) P. Geurts, D. Ernst, and L. Wehenkel, Machine Learning. 63, 3 (2006).
  * Breiman (2001) L. Breiman, Machine Learning. 54, 5 (2001).
  * Friedman (2002) J. H. Friedman, Computational Statistics and Data Analysis 38, 367 (2002).
  * Han and Klemens (1993) Y. J. Han and P. G. Klemens, Phys. Rev. B 48, 6033 (1993), ISSN 01631829.
  * Han (1996) Y.-J. Han, Phys. Rev. B 54, 8977 (1996), ISSN 0163-1829.
  * Herring (1954) C. Herring, Phys. Rev. 95, 954 (1954), ISSN 0031899X, URL http://link.aps.org/doi/10.1103/PhysRev.95.954.
  * Tamura (1983) S. I. Tamura, Phys. Rev. B 27, 858 (1983), ISSN 01631829.
  * Klemens (1955) P. G. Klemens, Proc. Phys. Soc. A 68, 1113 (1955), ISSN 0370-1298.
  * Klemens (1960) P. Klemens, Phys. Rev. 119, 507 (1960), ISSN 0031-899X, URL http://link.aps.org/doi/10.1103/PhysRev.119.507.
  * Inyushkin et al. (2003) A. V. Inyushkin, A. N. Taldenkov, A. Y. Yakubovsky, A. V. Markov, L. Moreno-Garsia, and B. N. Sharonov, Semicond. Sci. Tech. 18, 685 (2003), ISSN 02681242, URL http://stacks.iop.org/0268-1242/18/i=7/a=315.
  * Inyushkin et al. (2004) A. V. Inyushkin, A. N. Taldenkov, A. M. Gibin, A. V. Gusev, and H.-J. Pohl, Phys. Status Solidi (C) 1, 2995 (2004), ISSN 16101634, URL http://doi.wiley.com/10.1002/pssc.200405341.
  * Torres et al. (2019) P. Torres, J. Íñiguez, and R. Rurali, Phys. Rev. Lett. 123, 185901 (2019), URL https://link.aps.org/doi/10.1103/PhysRevLett.123.185901.
  * Ong et al. (2013) S. P. Ong, W. D. Richards, A. Jain, G. Hautier, M. Kocher, S. Cholia, D. Gunter, V. L. Chevrier, K. A. Persson, and G. Ceder, Computational Materials Science 68, 314 (2013), ISSN 0927-0256, URL http://www.sciencedirect.com/science/article/pii/S0927025612006295.
  * Togo (2020) A. Togo, _Phonon database at Kyoto university_ (2020), URL http://phonondb.mtl.kyoto-u.ac.jp/index.html.
  * Zhou et al. (2015) J. Zhou, B. Liao, B. Qiu, S. Huberman, K. Esfarjani, M. S. Dresselhaus, and G. Chen, Proceedings of the National Academy of Sciences 112, 14777 (2015), ISSN 0027-8424, eprint https://www.pnas.org/content/112/48/14777.full.pdf, URL https://www.pnas.org/content/112/48/14777.
  * Kuleev et al. (2003) I. Kuleev, I. Kuleev, A. Taldenkov, A. Inyushkin, V. Ozhogin, K. Itoh, and E. Haller, Journal of Experimental and Theoretical Physics 96, 1078 (2003).
  * Levchenko and Schmalian (2020) A. Levchenko and J. Schmalian, Annals of Physics 419, 168218 (2020), ISSN 0003-4916, URL https://www.sciencedirect.com/science/article/pii/S0003491620301524.
  * Gooth et al. (2018) J. Gooth, F. Menges, N. Kumar, V. SüB, C. Shekhar, Y. Sun, U. Drechsler, R. Zierold, C. Felser, and B. Gotsmann, Nat. Commun. 9, 4093 (2018).
  * Kuleyev et al. (2016) I. Kuleyev, I. Kuleyev, S. Bakharev, and V. Ustinov, J. Exp. Theor. Phys. 123, 489 (2016), ISSN 0003-4916.
  * Huang and Lucas (2021) X. Huang and A. Lucas, Phys. Rev. B 103, 155128 (2021), URL https://link.aps.org/doi/10.1103/PhysRevB.103.155128.
  * Glazov (2019) M. M. Glazov, Phys. Rev. B 100, 045426 (2019), URL https://link.aps.org/doi/10.1103/PhysRevB.100.045426.
  * Jiang et al. (2019) E. Jiang, X. Zhu, T. Ouyang, C. Tang, J. Li, C. He, C. Zhang, and J. Zhong, Journal of Applied Physics 126, 185106 (2019), eprint https://doi.org/10.1063/1.5121517, URL https://doi.org/10.1063/1.5121517.
  * Azarmi et al. (2019) Z. Azarmi, M. Naseri, and S. Parsamehr, Chemical Physics Letters 728, 136 (2019), ISSN 0009-2614, URL http://www.sciencedirect.com/science/article/pii/S0009261419303902.
  * Kansara et al. (2020) S. Kansara, Y. Sonvane, P. N. Gajjar, and S. K. Gupta, RSC Adv. 10, 26804 (2020), URL http://dx.doi.org/10.1039/D0RA03599H.
  * Zhao et al. (2017) T. Zhao, Y. Sun, Z. Shuai, and D. Wang, Chemistry of Materials 29, 6261 (2017), eprint https://doi.org/10.1021/acs.chemmater.7b01343, URL https://doi.org/10.1021/acs.chemmater.7b01343.
  * Lucy (2006) I. B. Lucy, Journal of Materials Science 41, 2517 (2006).
  * Snedaker (2015) M. L. Snedaker, _Alternative Approaches to Group IV Thermoelectric Materials_ (PhD dissertation, University of California, Santa Barbara, 2015).
  * Narushima et al. (2000) S. Narushima, H. Hosono, J. Jisun, T. Yoko, and K. Shimakawa, Journal of Non-Crystalline Solids 274, 313 (2000), ISSN 0022-3093, physics of Non-Crystalline Solids 9, URL http://www.sciencedirect.com/science/article/pii/S0022309300001940.
  * Shiraishi et al. (2017) Y. Shiraishi, S. Hata, Y. Okawauchi, K. Oshima, H. Anno, and N. Toshima, Chemistry Letters 46, 933 (2017), eprint https://doi.org/10.1246/cl.170218, URL https://doi.org/10.1246/cl.170218.