# The normal Casimir–Lifshitz force for laterally moving graphene

Mauro Antezza, Nataliia Emelianova, Nail Khusnutdinov

► To cite this version:

Mauro Antezza, Nataliia Emelianova, Nail Khusnutdinov. The normal Casimir–Lifshitz force for laterally moving graphene. Nanotechnology, 2024, 35 (23), pp.235001. 10.1088/1361-6528/ad2f1c. hal-04608544

HAL Id: hal-04608544

https://hal.science/hal-04608544v1

Submitted on 24 Mar 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Casimir-Lifshitz force for moving graphene

Mauro Antezza*

Laboratoire Charles Coulomb (L2C), UMR 5221 CNRS-University of Montpellier, F-34095 Montpellier, France and
Institut Universitaire de France, 1 rue Descartes, Paris Cedex 05 F-75231, France

Natalia Emelianova† and Nail Khusnutdinov‡

CMCC, Universidade Federal do ABC, Avenida dos Estados 5001, CEP 09210-580, SP, Brazil

(Dated: March 7, 2023)

We consider the system of two parallel sheets of graphene which are moving with relative parallel velocity $v$ and calculate the Casimir energy by using the scattering approach. The energy has real and imaginary parts. The real part contains a velocity contribution to the force perpendicular to the sheets. The imaginary part is related to a friction force parallel to the graphene sheets. We prove that the friction appears if the relative velocity becomes greater than the velocity of Fermi and if the modulo of reflection coefficient is greater than the one which is significant for virtual photon production. We analyze in detail the real contribution to the Casimir energy for two systems – graphene/graphene and ideal metal/graphene. In the non-relativistic case $v \ll v_F$, the relative correction to the Casimir energy $(\varepsilon_v - \varepsilon_0)/\varepsilon_0$ is proportional to the $(v/c)^2$ (the maximum value is 0.0033 for the gapless case and $v = v_F$) for the first system, and it is zero up to the Fermi velocity $v = v_F$ for system ideal metal/graphene.

# I. INTRODUCTION

Nowadays, the Casimir effect [1] is a well-known and exper-
imentally observed phenomenon – the macroscopic evidence
of the microscopic quantum origin of fields and matter (see, for
example, recent review [2] and book [3]). The Casimir effect for
real materials depends on the physical properties of solids and
shapes of their boundaries, temperatures of bodies and fields,
chemical potential, etc. [3]. The different temperatures of the
solids lead to heat transfer, and the non-equilibrium Casimir
effect [4] takes place. The new materials such as graphene [5],
Weyl semi-metals [6], and topological matter [7] demonstrate
the great diversity of properties of the Casimir effect.

Graphene is a mono-atomic two-dimensional layer of Car-
bon atoms, and the main role of the Casimir effect for graphene
plays its 2D conductivity tensor. In the frameworks of Kubo ap-
proach it was calculated in Ref. [8] and by using (2 + 1)D quan-
tum electrodynamics in Ref. [9].

The Casimir force depends on the relative velocity of slabs, too. The uniform motion perpendicular to the slabs was considered in Ref. [10] for scalar field and electromagnetic fields (see also the recent review [11] on the dynamical Casimir effect). The relative movement with velocity parallel to the slabs gives two different contributions to the Casimir force – perpendicular and parallel to the slabs. The first one, for dielectric layers, was considered in Ref. [12], where repulsive Casimir force was observed. The second contribution is called van der Waals/Casimir (quantum or non-contact) friction. It was investigated for a long time (see, for example, [13] and reference therein and the recent review [14]).

The Casimir friction for graphene was considered in Ref. [15]

in the framework of effective action. It was obtained that the
friction appears as an imaginal part of effective action and the
threshold for the friction force appears – it is zero if the relative
velocity is smaller than the velocity Fermi. This statement is in
agreement with Ref. [16] where it was shown that a frictional
force arises when the relative velocity exceeds that of light in the
medium. The Dirac electron is described by the Dirac equation
with velocity Fermi instead of the speed of light.

In this paper, we use the following approach to obtain the
Casimir energy of graphene moving with relative parallel veloc-
ity. By using the Lorentz transformation of graphene's conduc-
tivity tensor [9] we obtain the tensor of the moving graphene
and then use the formulas obtained in Ref. [17] to obtain the
real and imaginary parts of the Casimir energy and the force. To
do this we rotate the contour of integration from real frequen-
cies to imaginary ones. The integration over the imaginary axis
gives the real contribution to the Casimir force perpendicular to
the sheets of graphene. The specific complex frequencies pro-
duce an imaginary contribution which is the Casimir friction –
the force parallel to the sheets. In this paper, we consider the
real contribution, only.

The remainder of the paper is organized as follows. In Sec. II we consider the Lorentz transformation for the action. In Sec. III the tensor of conductivity of moving graphene in the laboratory frame is derived. We use the definition of conductivity as a boundary condition in the form of Ohm low and the Lorentz transformation to the laboratory frame.

The Sec. IV is devoted to the derivation of the Casimir energy and force. We start from the expression for the energy at the real frequencies and then move the integration to the imaginary axes. Two parts of the energy appeared. The part of the energy with integration along the imaginary frequencies $\xi = -i\omega$ (Sec. IVA) gives the real contribution, only to the energy. This statement is the consequence of the general relations for the conductivity tensor. These relations are proved in the Appendix A. The imaginary contribution to the energy and

<footnote>* mauro.antezza@umontpellier.fr</footnote>

<footnote>† natalia7emelianova@gmail.com</footnote>

<footnote>‡ nail.khusnutdinov@gmail.com</footnote>
force may come from the last part, only (Sec. IV B) which has the form of contour integrals along close contours in the upper and down half-planes of complex frequencies $\omega$ which corresponds to the opposite direction of scattering. The non-zero contributions from these integrals exist if some specific condition is satisfied. This condition in the scalar form has been noted by Lifshitz in Ref. [1] and also in Ref. [18]). In this section, we proved the threshold for velocity – the imaginary contribution to the energy is zero up to the velocity of Fermi. We proved also that this contribution takes pales if and only if the modulo of the reflections coefficient is greater than 1, that is the photon production has to appear.

In Sec. V we consider in detail the velocity correction to the perpendicular force, due to the real part of the energy, and for two systems – graphene/graphene and graphene/ideal metal cases.

The direct calculations of the boosted polarization tensor are shown in Appendix B.

We use the natural units $\hbar = c = 1$ and the following notations for indexes: $\mu, \nu = 0, 1, 2, 3 = t, x, y, z$; $a, b = 0, 1, 2 = t, x, y$, and $i, j = 1, 2 = x, y$. The wave vector, $k^a = (\omega, k)$, and $k'^a = (\omega', k')$.

## II. GENERAL REMARKS

The presence of the velocity of Fermi instead of the velocity of light in the 3D Dirac equation breaks the Lorentz invariance of the action. Indeed, in the laboratory frame, $K$ the action for a single graphene sheet reads

$$ \mathcal{A} = -\frac{1}{4} \int d^4 x F^2 + \int d^3 x \bar{\psi} [\tilde{\gamma}^b (i \partial_b - e A_b) - m] \psi |_{z=a}, \quad (1) $$

where $\tilde{\gamma}^a = \tilde{g}_b^a \gamma^b$ and $\tilde{g}_b^a = \text{diag}(1, v_F, v_F)$. The graphene is situated perpendicular to axis $z$ at point $z = a$. Then we make a boost in the plane of graphene with velocity $v$ to the moving frame $K'$: $x'^a = \Lambda_{.b}^{a'} x^b$, where the Lorentz matrix reads

$$ \Lambda = [\Lambda_{.b}^{a'}] = \begin{pmatrix} \gamma & -\gamma v \\ -\gamma v & I + (\gamma - 1) \frac{v \otimes v}{v^2} \end{pmatrix}, \quad (2) $$

and $\gamma = 1/\sqrt{1-v^2}$. The spinor transforms with matrix $S = \exp(-i/4 \sigma^{ab} \lambda_{ab})$: $\psi'(x') = S\psi(x)$, $\tilde{\psi}'(x') = \tilde{\psi}(x)S^{-1}$, where

$$ \sigma^{ab} = \frac{i}{2} [\gamma^a, \gamma^b], \quad \lambda_{ab} = \begin{pmatrix} 0 & -v^1 & -v^2 \\ v^1 & 0 & 0 \\ v^2 & 0 & 0 \end{pmatrix}. \quad (3) $$

The action in the frame $K'$ becomes

$$ \mathcal{A}' = -\frac{1}{4} \int d^4 x' F'^2 + \int d^3 x' \bar{\psi}' [\tilde{\Gamma}^b (i \partial_{b'} - e A_{b'}^c) - m] \psi' |_{z=a}, $$

where

$$ \begin{aligned} \tilde{\Gamma}^{a'} &= \Lambda_{.b}^{a'} S \tilde{\gamma}^b S^{-1} = \Lambda_{.b}^{a'} \Lambda_{.e'}^c \tilde{g}_c^b \gamma^{e'} = \tilde{g}_{e'}^{a'} \gamma^{e'} \\ &= v_F \gamma^a + u^a (1 - v_F) \psi, \end{aligned} \quad (4) $$

and $\psi = \gamma^a u_a$, $u_a = (u_0, u_1, u_2) = \gamma(1, -v)$. Here the matrix $\tilde{g}' = \Lambda^{-1}\tilde{g}\Lambda$ is the Lorentz transform of $\tilde{g}$. We observe that $\tilde{\Gamma}^{a'} \neq \tilde{\gamma}^{a'}$ if $v \neq 0$, and therefore, $\mathcal{A}' \neq \mathcal{A}$. For $v_F = 1$ we obtain the Lorentz invariance, $\tilde{\Gamma}^{a'} = \tilde{\gamma}^{a'} = \gamma^{a'}$, as should be the case.

## III. THE CONDUCTIVITY OF MOVING GRAPHENE

Let us consider a system of two parallel graphenes which are situated perpendicular to the axis $z$ and are separated by distance $a$ apart. The first graphene is at rest in the laboratory frame $K$ and the second graphene is moving parallel to the first graphene with the 3-velocity $u^a = (u^0, u)$. The first graphene is described by the tensor conductivity $\sigma_1$ and the second one – by the tensor conductivity $\sigma_2$. Both conductivities are for the laboratory frame $K$.

The expression for tensor conductivity of graphene has been calculated in the co-moving frame $K'$ in Ref. [9] for zero temperature and non-zero temperatures in terms of the polarization tensor in 3D quantum electrodynamics. To obtain the tensor conductivity of graphene $\sigma_2$ in the laboratory frame $K$ we make the Lorentz boost from the co-moving frame $K'$ to the laboratory frame $K$. The conductivity of the second graphene in its co-moving frame one denotes with prime, $\sigma_2'$.

Let us repeat the derivation of this tensor and extend it to the problem under consideration. The boundary conditions in the co-moving frame $K'$, where graphene is at rest, have the following form

$$ \begin{aligned} [A'_{\mu}(x')]_{z'=a} &= 0, \\ [F'^{zv}(x')]_{z'=a} &= \int d^3 \tilde{x}'' \tilde{\Pi}^{\prime v\rho} (x' - x'') A'_{\rho}(x'') |_{z'=z''=a}, \end{aligned} \quad (5) $$

where $[f]_{z'=a} = f_{z'=a-0} - f_{z'=a+0}$ and $\tilde{\Pi}^{\prime 3\mu} = \tilde{\Pi}^{\prime \mu 3} = 0$. In what follows, we will assume that $z = z' = z'' = a$. Then, we make the Fourier transform with respect to the three coordinates $x'^a$ and obtain

$$ \begin{aligned} [A'_{\mu}(k')] &= 0, \\ [F'^{zv}(k')] &= \tilde{\Pi}^{\prime v\rho}(k') A'_{\rho}(k'). \end{aligned} \quad (6) $$

Let us consider first the case $v = i$. By definition $F'_{i0} = \partial_{i'} A'_0 - \partial_{0'} A'_i = -E'_i(x')$ and therefore $A'_i = E'_i/i\omega' + A'_0 k_i/\omega'$. Then, by using the gauge invariance condition $\tilde{\Pi}^{\prime v\rho}(k') k'_\rho = 0$, the second boundary condition (6) is represented in the following form

$$ [F'^{zi}(k')] = \frac{1}{i\omega'} \tilde{\Pi}^{\prime ij} E'_j. \quad (7) $$

It has the form of boundary condition for tangent components of a magnetic field

$$ [\tilde{H}'] \times \tilde{n} = \tilde{J}', \quad (8) $$

with surface current

$$ J^{\prime i} = \sigma^{\prime ij} E'_{j}, \quad (9) $$
and normal vector $n^a = \delta_z^a$. Therefore,

$$ \sigma^{ij} = \frac{1}{i\omega'} \tilde{\Pi}^{ij}. \qquad (10) $$

By setting now the index $v=0$ in Eq. (6) we obtain the boundary condition for the normal component of the electric field

$$ [\tilde{E}'] \cdot \tilde{n} = -\frac{1}{i\omega'} \tilde{\Pi}^{0j} E'_j = -\rho'. \qquad (11) $$

To obtain the conductivity tensor, $\sigma$, in the laboratory frame $K$, we make boost, $K' \to K: x'^a = \Lambda_{.b}^{a'} x^b, k'_a = \Lambda_{.a'}^b k_b$, in the boundary equations (6) to the laboratory frame with matrix

$$ \Lambda = [\Lambda_{.b}^{a'}] = \begin{pmatrix} u^0 & -\mathbf{u} \\ -\mathbf{u} & \mathbf{I} + \frac{\mathbf{u} \otimes \mathbf{u}}{u^0+1} \end{pmatrix}, \qquad (12) $$

where $\mathbf{u} = (u^1, u^2) = v\gamma$ and $\mathbf{I}$ is 2 × 2 identical matrix. In matrix notations, $x' = \Lambda x$ and $k' = k\Lambda^{-1}$. Then, in the laboratory frame $K$ we obtain the boundary conditions

$$ \begin{aligned} [A_{\mu}(k)] &= 0, \\ [F^{zv}(k)] &= \tilde{\Pi}^{v\rho}(k) A_{\rho}(k), \end{aligned} \qquad (13) $$

where

$$ \tilde{\Pi}^{v\rho}(k) = \Lambda_{.v'}^{v} \Lambda_{.\rho'}^{\rho'} \tilde{\Pi}^{v'\rho'}(k'). \qquad (14) $$

The same expressions for $\tilde{\Pi}^{v\rho}(k)$ may be obtained by direct calculation of boosted Dirac equations (see Appendix B).

By doing the same steps (6)-(10) but with wave vector $k$ instead of $k'$, we obtain the conductivity tensor of graphene in the laboratory frame $K$

$$ \sigma^{ij}(k) = \frac{1}{i\omega} \tilde{\Pi}^{ij}(k) = \frac{1}{i\omega} \Lambda_{.v}^{i} \Lambda_{.\rho'}^{j} \tilde{\Pi}^{v'\rho'}(k'). \qquad (15) $$

The PT has two contributions namely, symmetric (parity-even) and antisymmetric (parity-odd) $\tilde{\Pi}^{\mu\nu} = \tilde{\Pi}_s^{\mu\nu} + \tilde{\Pi}_a^{\mu\nu}$, where $\tilde{\Pi}_a^{\mu\nu} = \xi e^{i\nu\alpha} k'_a$. Due to gauge invariance $\tilde{\Pi}^{\mu\nu} k'_\mu = \tilde{\Pi}_s^{\mu\nu} k'_\mu \equiv 0$. We use these relations to express time components

$$ \tilde{\Pi}_s^{00} = \frac{1}{\omega'^2} \tilde{\Pi}_s^{ij} k'_i k'_j, \quad \tilde{\Pi}_s^{0i} = -\frac{1}{\omega'} \tilde{\Pi}_s^{ij} k'_j. \qquad (16) $$

Taking these relations into account we obtain

$$ \tilde{\Pi}^{ij}(k) = \tilde{\Pi}^{i'n'm'} G_{.n'}^i G_{.m'}^j, \qquad (17) $$

where

$$ G_{.n'}^i = \delta_{n'}^i - \frac{u^i u_n}{u_0 + 1} - \frac{u^i k'_n}{\omega'}. \qquad (18) $$

By using manifest form of $k' = k\Lambda^{-1}$,

$$ \omega' = (ku), k' = k - \mathbf{u} \frac{\omega + (ku)}{1 + u_0}, \qquad (19) $$

we obtain

$$ G_{.n'}^i = \delta_{n'}^i + \frac{\omega u^i u_n}{(ku)(u_0 + 1)} - \frac{u^i k_n}{(ku)}, \qquad (20) $$

where $(ku) = k^a u_a$. In the matrix notations, one has

$$ \mathbf{G} = \mathbf{I} - \frac{\omega}{(ku)} \frac{\mathbf{u} \otimes \mathbf{u}}{u_0 + 1} + \frac{\mathbf{u} \otimes \mathbf{k}}{(ku)}, \qquad (21) $$

and

$$ \tilde{\Pi} = \mathbf{G} \tilde{\Pi}' \mathbf{G}^T. \qquad (22) $$

Therefore, the conductivity tensor $\eta_i = 2\pi\sigma_i$ in the laboratory frame is

$$ \eta_2 = \frac{\omega'}{\omega} \mathbf{G} \eta_2' \mathbf{G}^T, \qquad (23) $$

where $\eta_2'$ is the conductivity tensor of the second plane in its frame, which depends on the $\omega'$, $k'$. The straightforward calculations give

$$ \det \mathbf{G} = \frac{\omega}{\omega'} = \frac{\omega}{(ku)}, \qquad (24) $$

and therefore $\det \eta_2 = \det \eta_2'$.

A graphene at rest has the following conductivity [9]

$$ \eta_1 = \eta_{\text{gr}} \frac{\tilde{k}}{\omega} \left( \mathbf{I} + v_F^2 \frac{\mathbf{k} \otimes \mathbf{k}}{\tilde{k}^2} \right) \tilde{\Phi} \left( \frac{\tilde{k}}{2m} \right), \qquad (25) $$

where $\eta_{\text{gr}} = 2\pi\sigma_{\text{gr}} = \pi\alpha/2 = 0.0114$, and

$$ \tilde{\Phi}(y) = \frac{2i}{\pi y} \left\{ 1 - \frac{y^2 + 1}{y} \operatorname{arctanh} y \right\}, \\ \tilde{k} = \sqrt{\omega^2 - v_F^2 \mathbf{k}^2}. $$

The conductivity of the moving graphene has form (23) where

$$ \eta_2' = \eta_{\text{gr}} \frac{\tilde{k}'}{\omega'} \left( \mathbf{I} + v_F^2 \frac{\mathbf{k}' \otimes \mathbf{k}'}{\tilde{k}'^2} \right) \tilde{\Phi} \left( \frac{\tilde{k}'}{2m} \right). \qquad (26) $$

Here $\tilde{k}' = \sqrt{\omega'^2 - v_F^2 \mathbf{k}'^2}$, or by using the Lorentz boost (19) one obtains

$$ \tilde{k}' = \sqrt{\tilde{k}^2 + \frac{1 - v_F^2}{1 - v^2} (\omega^2 v^2 + (kv)^2 - 2\omegav)}. \qquad (27) $$

The tensor (23) is symmetric and depends on the vectors $\mathbf{k}$ and $\mathbf{v}$. Therefore, the general structure of this tensor is the following

$$ \eta_2 = A_1 \mathbf{I} + A_2 \mathbf{k} \otimes \mathbf{k} + A_3 \mathbf{v} \otimes \mathbf{v} + A_4 (\mathbf{k} \otimes \mathbf{v} + \mathbf{v} \otimes \mathbf{k}). \qquad (28) $$

Calculating different contractions of the above relation with $\mathbf{I}$, $\mathbf{k}$, and $\mathbf{v}$ we obtain

$$ \eta_2 = i_1 \mathbf{I} + i_2 \mathbf{k} \otimes \mathbf{k} + i_3 (\mathbf{k} \otimes \mathbf{v} + \mathbf{v} \otimes \mathbf{k}), \qquad (29) $$
where

$$
\begin{align}
i_1 &= \frac{\eta_{\text{gr}} \tilde{\Phi}'}{\omega k^2 \tilde{k}'} \left( k^2 \tilde{k}'^2 + \frac{1 - v_F^2}{1 - v^2} ((kv)^2 - k^2 v^2) k_3^2 \right), \nonumber \\
i_2 &= \frac{\eta_{\text{gr}} \tilde{\Phi}'}{\omega k^2 \tilde{k}'} \left( k^2 v_F^2 + \frac{1 - v_F^2}{1 - v^2} v^2 k_3^2 \right), \nonumber \\
i_3 &= \frac{\eta_{\text{gr}} \tilde{\Phi}'}{k^2 \tilde{k}'} \frac{1 - v_F^2}{1 - v^2} (k^2 - \omega), \tag{30}
\end{align}
$$

and

$$
\tilde{\Phi}' = \tilde{\Phi} \left( \frac{\tilde{k}'}{2m} \right), \quad k_3^2 = \omega^2 - k^2.
$$

For $v = 0$ we obtain conductivity of graphene at rest (25).

IV. THE CASIMIR ENERGY AND FORCE

In the frameworks of the scattering matrix approach [17] the
Casimir energy per unit area for real frequencies reads

$$
\mathcal{E} = -\frac{1}{2i} \iint \frac{d^2 k}{(2\pi)^3} (I_- - I_+), \quad (31)
$$

where

$$
I_{\pm} = \int_{0}^{\infty} d\omega \ln \det \left[ 1 - e^{\pm 2i\omega k_3} \mathcal{R}(\pm k_3) \right], \quad (32)
$$

$\mathcal{R}(\pm k_3) = r'_1(\pm k_3)r_2(\pm k_3)$, and $k_3 = \sqrt{\omega^2 - k^2}$. The Eq. (31) contains all contributions to the Casimir energy – due to evanescent and waveguide surface contributions ($\omega < k$) and propagating waves ($\omega > k$). The subscript 1 (2) means that all reflection matrices related to the left (right) subsystem ($z = 0, a$). Two terms $I_\pm$ correspond to opposite directions of scattering.

The scattering matrix for each part of the system has the fol-
lowing form

$$
\mathcal{S} = \begin{pmatrix} r & t' \\ t & r' \end{pmatrix},
$$

with the corresponding index and argument. This matrix for the
conductive plane at position $z = a$ was found in Ref. [17]

$$
r = -e^{2ik_3 a} r_0, \quad r' = -e^{-2ik_3 a} r_0, \quad t = t' = \mathbf{I} + r_0, \quad (33)
$$

where

$$
r_0 = \frac{\omega^2 \eta - k \otimes (k\eta) + \mathbf{I} \omega k_3 \det \eta}{\omega^2 \operatorname{tr} \eta - k \mathbf{k} \eta + \omega k_3 (1 + \det \eta)},
$$

is the reflection matrix for plane at position $z = 0$.

The pressure is the negative derivative of the energy per unit
area with respect to distance $\mathcal{P} = -\mathcal{E}'_a$. With simple algebra we
obtain

$$
\mathcal{P} = \iint \frac{d^2 k}{(2\pi)^3} (J_- + J_+), \quad (34)
$$

where

$$
J_{\pm} = \int_{0}^{\infty} d\omega k_{3} \frac{e^{\pm 2i\omega k_{3}} (\operatorname{tr} \mathcal{R}(\pm k_{3}) - 2e^{\pm 2i\omega k_{3}} \det \mathcal{R}(\pm k_{3}))}{\det [\mathbf{I} - e^{\pm 2i\omega k_{3}} \mathcal{R}(\pm k_{3})]} . \quad (35)
$$

The energy (31) and the pressure (34) have real and imaginary
contributions due to velocity. The Casimir friction has to appear
as an imaginary part of the energy. The Casimir energy for real
frequencies has been considered in Ref. [19]. To separate real
and imaginary parts we shift the integration from the real axis
to the imaginary one as shown in Fig. 1. The sign of phase of k₃
in the complex plain of ω coincides with the sign of phase of ω
and defines the regions where the exponents become damping.
Therefore,

<figure><img src="image_2.png" /><figcaption>FIG. 1. Contours of integrations.</figcaption></figure>

$$
I_{\pm} : \int_{0}^{\infty} \Rightarrow \oint_{\gamma_{\pm}} - \int_{C_{\pm} A_{\pm}} - \int_{B_{\pm} C_{\pm}} \Rightarrow \sum_{i=1}^{3} \mathcal{E}_{i}, \quad (36)
$$

where limits, $\epsilon, \epsilon \to 0$, and radius of arcs, $R \to \infty$, are assumed.
The contributions of arcs $B_\pm C_\pm$ are zeros in the last limit, $\mathscr{E}_3 \to 0$.
Therefore, the energy is a sum of two contributions $\mathscr{E} = \mathscr{E}_1 + \mathscr{E}_2$. Let us consider consecutively these terms separately.

A. $\mathscr{E}_2$: The real contribution

The second integral alongside the imaginary axis reads

$$
\mathscr{E}_2 = \frac{1}{2i} \iint \frac{d^2 k}{(2\pi)^3} \left\{ \int_{C_- A_-} d\omega \ln \det \left[ 1 - e^{-2i\omega k_3} \mathcal{R}(-k_3) \right] - \int_{C_+ A_+} d\omega \ln \det \left[ 1 - e^{2i\omega k_3} \mathcal{R}(k_3) \right] \right\}. \quad (37)
$$
One changes integrand variable $\omega = i\xi + \epsilon$. In this case $k_3 = ik_E$ sign($\xi$), where $k_E = \sqrt{\xi^2 + k^2}$. Therefore,

$$ \mathcal{E}_2 = \frac{1}{2} \iint \frac{d^2 k}{(2\pi)^3} \left\{ \int_{-\infty}^{0} d\xi \ln \det \left[ 1 - e^{-2ak_E} \mathcal{R}(ik_E) \right] - \int_{+\infty}^{0} d\xi \ln \det \left[ 1 - e^{-2ak_E} \mathcal{R}(ik_E) \right] \right\}, \quad (38) $$

and we arrive with formula [17]

$$ \mathcal{E}_2 = \iint \frac{d^2 k}{2(2\pi)^3} \int_{-\infty}^{+\infty} d\xi \ln \det \left[ 1 - e^{-2ak_E} \mathcal{R}(ik_E) \right]. \quad (39) $$

Appendix A shows that the conductivities of graphene obey the general relation [20, §82] which must be fulfilled for any realistic conductivity

$$ \begin{aligned} \eta_1(-\omega^*) &= \eta_1^*(\omega), \\ \eta_2(-\omega^*, v) &= \eta_2^*(\omega, -v). \end{aligned} \quad (40) $$

Taking into account these relations, only, it is easy to show that the $\mathcal{E}_2$ is real (see Appendix A) and

$$ \mathcal{E}_2 = \Re \iint \frac{d^2 k}{(2\pi)^3} \int_0^\infty d\xi \ln \det \left[ \mathbf{I} - e^{-2ak_E} \mathcal{R}(ik_E) \right]. \quad (41) $$

This expression has represented in Ref. [17] in the following manifest form

$$ \begin{aligned} \mathcal{E}_2 = & \int \frac{d^2 k}{(2\pi)^3} \int_0^\infty d\xi \ln \left( 1 + e^{-4ak_E} \frac{\xi^2 k_E^2}{b_1 b_2} \det \eta_1 \det \eta_2 \right) \\ & - e^{-2ak_E} \left[ \frac{\xi^2 k_E^2}{b_1 b_2} \left[ (1 - \det \eta_1)(1 - \det \eta_2) + \det(\eta_1 - \eta_2) \right] \right. \\ & \left. - \frac{\xi k_E}{b_1} - \frac{\xi k_E}{b_2} + 1 \right] \end{aligned} \quad (42) $$

where $b_i = \xi^2 \text{tr} \eta_i + (kk\eta_i) + \xi k_E (1 - \det \eta_i)$.

## B. $\mathcal{E}_1$: The Casimir friction

Let us proceed to the last contribution which comes from the integrals over closed contours $\gamma_\pm$

$$ \mathcal{E}_1 = -\frac{1}{2i} \iint \frac{d^2 k}{(2\pi)^3} \left\{ \oint_{\gamma_-} d\omega \ln \det \left[ \mathbf{I} - e^{-2iak_3} \mathcal{R}(-k_3) \right] - \oint_{\gamma_+} d\omega \ln \det \left[ \mathbf{I} - e^{2iak_3} \mathcal{R}(k_3) \right] \right\}. \quad (43) $$

This part is the source of the Casimir friction because it is imaginary. The contribution to pressure reads

$$ \mathcal{P}_1 = \iint \frac{d^2 k}{(2\pi)^3} \times \left\{ \oint_{\gamma_-} d\omega k_3 \frac{e^{-2iak_3} (\text{tr} \mathcal{R}(-k_3) - 2e^{-2iak_3} \det \mathcal{R}(-k_3))}{\det [\mathbf{I} - e^{-2iak_3} \mathcal{R}(-k_3)]} \right\}. \quad (44) $$

$$ + \oint_{\gamma_+} d\omega k_3 \frac{e^{2iak_3} (\text{tr} \mathcal{R}(k_3) - 2e^{2iak_3} \det \mathcal{R}(k_3))}{\det [\mathbf{I} - e^{2iak_3} \mathcal{R}(k_3)]} \right\}. \quad (44) $$

Therefore, the imaginary part of force may appear due to zeros of the denominator

$$ \det \left[ \mathbf{I} - e^{\pm 2iak_3} r'_1(\pm k_3) r_2(\pm k_3) \right] = 0, \quad (45) $$

as residues in these points.

For the scalar Fresnel coefficients these poles have been noted by Lifshitz in Ref. [1] when the integration over real axes is changed to the imaginary one (see also [18]). In the Ref. [19] (see also review [21]) noted that the roots of this relation in the scalar case may appear in the complex plane, only. Lifshitz [1] proved an absence of poles for normal dielectrics.

Let us suppose that the matrix $\mathcal{R}(\pm k_3) = r'_1(\pm k_3)r_2(\pm k_3)$ has eigenvalues $r_{\text{tm}}^2(\pm k_3)$ and $r_{\text{te}}^2(\pm k_3)$, which means that the following relation is valid

$$ \mathcal{R}(\pm k_3) = \mathbf{T}^{-1}\mathbf{D}_{\pm}\mathbf{T}, \quad (46) $$

where

$$ \mathbf{D}_{\pm} = \begin{pmatrix} r_{\text{tm}}^2(\pm k_3) & 0 \\ 0 & r_{\text{te}}^2(\pm k_3) \end{pmatrix}. \quad (47) $$

For $v = 0$, the $r_{\text{tm}}^2$ and $r_{\text{te}}^2$ are the squares of reflection coefficients of TM and TE modes respectively for a system of two planes with the corresponding sign of $k_3$.

Taking into account Eq. (46) and relation

$$ \oint_{\gamma_{-}} d\omega k_{3} + \oint_{\gamma_{+}} d\omega k_{3} = 0, $$

we obtain

$$ \mathcal{P}_1 = \iint \frac{d^2 k}{(2\pi)^3} \left\{ \oint_{\gamma_-} d\omega k_3 \left[ \frac{1}{1 - e^{-2iak_3} r_{\text{tm}}^2(-k_3)} + \frac{1}{1 - e^{-2iak_3} r_{\text{te}}^2(-k_3)} \right] + \oint_{\gamma_+} d\omega k_3 \left[ \frac{1}{1 - e^{-2iak_3} r_{\text{tm}}^2(k_3)} + \frac{1}{1 - e^{-2iak_3} r_{\text{te}}^2(k_3)} \right] \right\}. \quad (48) $$

Therefore, the contributions to these integrals come from the roots

$$ 1 - e^{\pm 2iak_3} r_{\text{tm}}^2(\pm k_3) = 0, \quad 1 - e^{\pm 2iak_3} r_{\text{te}}^2(\pm k_3) = 0, \quad (49) $$

as residues in these poles if the poles are inside the corresponding contour. For real frequencies, these relations describe the waveguide modes [22], which are the evanescent modes with $\omega < k$. A similar relation was considered in Ref. [23].

Let us denote the solutions of equations

$$ 1 - e^{-2iak_3} r_{\text{tm}}^2(-k_3) = 0, \quad 1 - e^{-2iak_3} r_{\text{te}}^2(-k_3) = 0, \quad (50) $$

by $\omega_{\text{tm}}^{n-}$, $\omega_{\text{te}}^{n-}$ with negative imaginary parts $\Im\omega_{\text{tm}}^{n-} < 0$, $\Im\omega_{\text{te}}^{n-} < 0$, and the solutions of equations

$$ 1 - e^{2iak_3} r_{\text{tm}}^2(k_3) = 0, \quad 1 - e^{2iak_3} r_{\text{tm}}^2(k_3) = 0, \quad (51) $$
by $\omega_{\text{tm}}^{n+}, \omega_{\text{te}}^{n+}$ with positive imaginary parts $\Im\omega_{\text{tm}}^{n+} > 0, \Im\omega_{\text{te}}^{n+} > 0$. Then, the pressure reads

$$ \mathcal{P}_1 = i \sum_n \iint \frac{d^2 k}{(2\pi)^2} \times \left[ \begin{aligned} & \operatorname{Res}_{\omega=\omega_{\text{tm}}^{n+}} \frac{k_3}{1 - e^{2iak_3} r_{\text{tm}}^2(k_3)} - \operatorname{Res}_{\omega=\omega_{\text{tm}}^{n-}} \frac{k_3}{1 - e^{-2iak_3} r_{\text{tm}}^2(-k_3)} \\ & + \operatorname{Res}_{\omega=\omega_{\text{te}}^{n+}} \frac{k_3}{1 - e^{2iak_3} r_{\text{te}}^2(k_3)} - \operatorname{Res}_{\omega=\omega_{\text{te}}^{n-}} \frac{k_3}{1 - e^{-2iak_3} r_{\text{te}}^2(-k_3)} \end{aligned} \right]. \quad (52) $$

Let us analyze the general properties of existing poles in the contours following by Lifshitz [1]. We represent the relations (49) in the following form

$$ e^{\pm 2iak_3} = r_{\text{tm}}^{-2}(\pm k_3), \quad e^{\pm 2iak_3} = r_{\text{te}}^{-2}(\pm k_3), \quad (53) $$

and consider the modulo of both sides.

It is easy to see that $\text{sign}(\Im k_3) = \text{sign}(\Im \omega)$. Therefore, the modulo of the left-hand side of Eqs. (53) reads

$$ |e^{\pm 2iak_3}| = (e^{\pm 2iak_3} e^{\mp 2iak_3^*})^{1/2} = e^{\mp 2a\Im k_3} = e^{-2a|\Im k_3|} < 1. \quad (54) $$

Therefore, the solutions exist if and only if the relations

$$ |r_{\text{tm}}^{-2}(\pm k_3)| < 1, \quad |r_{\text{te}}^{-2}(\pm k_3)| < 1, \quad (55) $$

are valid, or

$$ |r_{\text{tm}}^{2}(\pm k_3)| > 1, \quad |r_{\text{te}}^{2}(\pm k_3)| > 1. \quad (56) $$

It means that it has to be photon production.

Let us consider now the threshold of velocity for the Casimir friction. It was found in the Ref. [15] in the framework of effective action. Here we prove this threshold by using the simple observation that the origin for the Casimir friction is the mixing between positive and negative frequencies [16].

For graphene at rest, $v = 0$, the spectrum is defined by relation $\tilde{k} = \sqrt{\omega^2 - v_F^2 k^2} = 0$. We obtain two branches

$$ \omega_0^{\pm}/k = \pm v_F. \quad (57) $$

For moving graphene the spectrum is defined by relation $\tilde{k}' = \sqrt{\omega'^2 - v_F^2 k'^2} = 0$. We have a quadratic equation

$$ \omega^2 - v_F^2 k^2 + \frac{1 - v_F^2}{1 - v^2} (\omega^2 v^2 + (kv)^2 - 2\omega(kv)) = 0, \quad (58) $$

with solutions

$$ \begin{aligned} \omega_{\nu}^{\pm}/k = v & \frac{1 - v_F^2}{1 - v^2 v_F^2} \cos\varphi \\ & \pm v_F \frac{\sqrt{1 - v^2}}{1 - v^2 v_F^2} \sqrt{1 - v^2 v_F^2 - v^2(1 - v_F^2)\cos^2\varphi}. \end{aligned} \quad (59) $$

The expression under the square root is always positive. Here, the $\varphi$ is the angle between vectors $k$ and $v$. It is easy to verify that $\lim_{v \to 0} \omega_{\nu}^{\pm} = \omega_0^{\pm}$.

To find a condition for mixing the frequencies let us first for simplicity expand (59) over velocity $v < 1$. We set $v = v_F\beta$ and expand over $v_F \ll 1$:

$$ \omega_{\nu}^{\pm}/k = v_F(\pm 1 + \beta \cos\varphi) + O(v_F^3). \quad (60) $$

For $\beta < 1$ ($v < v_F$) we observe that $\omega_{\nu}^{+} > 0$ and $\omega_{\nu}^{-} < 0$ for all $\varphi$, that is there is no mixing between modes. In the opposite case $\beta > 1$ ($v > v_F$), there are domains of angle $\varphi$ where $\omega_{\nu}^{+} < 0$ or $\omega_{\nu}^{-} > 0$ – the mixing appears and therefore the Casimir friction appears.

This statement is valid without expansion over $v$. Indeed, the condition $\omega_{\nu}^{+} < 0$ reads

$$ v_F \frac{\sqrt{1 - v^2}}{1 - v^2 v_F^2} \sqrt{1 - v^2 v_F^2 - v^2(1 - v_F^2)\cos^2\varphi} < -v \frac{1 - v_F^2}{1 - v^2 v_F^2} \cos\varphi. $$

The first observation is that $\varphi \in [\frac{\pi}{2}, \frac{3\pi}{2}]$, because expression on right hand side should be positive. With this limitation, one can put in a square on both sides and one obtains a condition

$$ \cos^2 \varphi > 1 - \frac{v^2 - v_F^2}{v^2 (1 - v_F^2)}. \quad (61) $$

We observe that this relation can not be valid if $v < v_F$. For $v > v_F$ there are domains of $\varphi$ where solutions exist and therefore the mixing of modes exists, and the Casimir friction comes up.

The condition $\omega_{\nu}^{-} > 0$ reads

$$ v_F \frac{\sqrt{1 - v^2}}{1 - v^2 v_F^2} \sqrt{1 - v^2 v_F^2 - v^2(1 - v_F^2)\cos^2\varphi} < v \frac{1 - v_F^2}{1 - v^2 v_F^2} \cos\varphi. $$

Therefore, $\varphi \in [-\frac{\pi}{2}, \frac{\pi}{2}]$. With this limitation one can put in a square both sides and one obtains the same condition

$$ \cos^2 \varphi > 1 - \frac{v^2 - v_F^2}{v^2 (1 - v_F^2)}, \quad (62) $$

which has solutions for $v > v_F$, only.

<figure><img src="image_3.png" /><figcaption>FIG. 2. The spectrum (59). Left panel: $v/v_F = 0.5$ (there is no mixing). Right panel: $v/v_F = 2$ (there is mixing).</figcaption></figure>

Let us consider threshold for photon production. It happens if [16, 23]

$$ \omega_{\nu}^{+} + \omega_{0}^{+} < 0. \quad (63) $$
From Fig. 2 we observe that we have to consider this relation for $\varphi = \pm \pi$. It gives threshold for velocity

$$v > \frac{2v_F}{1 + v_F^2} \approx 2v_F, \qquad (64)$$

which is close to that observed in Refs. [16, 23].

The detailed analysis of the Casimir friction we will make elsewhere in a separate paper.

## V. THE CASIMIR ENERGY: VELOCITY CORRECTION

Let us consider the contribution $\mathcal{E}_2$ to the Casimir energy given by Eq. (41) for two systems – graphene/graphene (gg) and graphene/ideal metal (gm). We consider first the zero-velocity energy $\mathcal{E}_0 \equiv \mathcal{E}_2|_{v=0}$ and then the velocity correction which we will denote $\Delta_v \mathcal{E}_2 = \mathcal{E}_2 - \mathcal{E}_0$. Also, we will give expressions for velocity correction to the Casimir pressure $\Delta_v \mathcal{P}_2 = \mathcal{P}_2 - \mathcal{P}_0$.

### A. Two graphene

The zero-velocity Casimir energy reads

$$\mathcal{E}_0^{gg} = \int_0^\infty \frac{k_E^2 dk_E}{(2\pi)^2} \int_0^1 dx \\ \times \left\{ \ln\left(1 - e^{-2ak_E} r_{\text{tm}}^2\right) + \ln\left(1 - e^{-2ak_E} r_{\text{te}}^2\right) \right\}, \quad (65)$$

where

$$r_{\text{tm}}^{-1} = 1 + \frac{t}{\eta_{\text{gr}} \Phi \left( \frac{k_F t}{2ma} \right)}, \quad r_{\text{te}}^{-1} = -1 - \frac{1}{\eta_{\text{gr}} t \Phi \left( \frac{k_F t}{2ma} \right)}, \quad (66)$$

are the reflection coefficients for TM and TE modes [9]. Here, $\Phi(y) = \bar{\Phi}(iy)$, and we used spherical coordinates, $k_E = \sqrt{\xi^2 + k^2}$, $x = \cos\theta$, and $t = \sqrt{x^2 + v_F^2(1 - x^2)}$. Let us estimate the $\mathcal{E}_0^{gg}$ for large and small separations.

For $a \to \infty$, the function $\Phi(k_F t/2ma) \approx 4k_F t/3\pi ma$ and we obtain

$$\frac{\mathcal{E}_0^{gg}}{\mathcal{E}_{\text{id}}} \approx \frac{32\eta_{\text{gr}}^2 (9 + 2v_F^2 + 4v_F^4)}{\pi^6 (am)^2} \approx \frac{288\eta_{\text{gr}}^2}{\pi^6 (am)^2}, \\ \frac{\mathcal{P}_0^{gg}}{\mathcal{P}_{\text{id}}} \approx \frac{160\eta_g r^2 (9 + 2v_F^2 + 4v_F^4)}{3\pi^6 (am)^2} \approx \frac{480\eta_{\text{gr}}^2}{\pi^6 (am)^2}, \quad (67)$$

where $\mathcal{E}_{\text{id}} = -\pi^2/720a^3$ is the Casimir energy for two ideal plates. Therefore, for large separation $\mathcal{E}_0^{gg} \approx -2\eta_{\text{gr}}^2/5\pi^4 a^5 m^2$.

For $a \to 0$, the function $\Phi \to 1$ (all expressions see in Appendix C). Numerically, $\mathcal{E}_0^{gg}/\mathcal{E}_{\text{id}} \approx 0.0048$, for $v_F = 1/300$, and therefore $\mathcal{E}_0^{gg} \sim 1/a^3$.

Note, that for constant conductivity case [24] the energy is proportional to the first degree of $\eta_{\text{gr}}$ due to TM mode and falls as $1/a^3$ for any separations. The constant conductivity case is realized for $v_F = 0$ and $\Phi = 1$ (see Eq. (25)). But setting $v_F = 0$

does not lead to relation $\Phi = 1$. The dispersion over angle $\theta$ is still survived. The case $\Phi = 1$ is realized for small separation. By setting $v_F = 0$ in this relation we obtain $\mathcal{E}_0^{gg}/\mathcal{E}_{\text{id}} \to 45\eta_{\text{gr}}/\pi^4$ which coincides with result of Ref. [24] for constant conductivity case. The large separation case is completely different from the constant conductivity case – the dispersion over angle $\theta$ and non-zero gap play an important role.

The velocity correction term $\mathcal{E}_2^{gg}$ has a complicated form and no reason to reproduce it here. For large separation, $(am) \to \infty$ and arbitrary velocity

$$\frac{\Delta_v \mathcal{E}_2^{gg}}{\mathcal{E}_{\text{id}}} \approx \frac{112\eta_{\text{gr}}^2 (1 - v_F^2)^2}{\pi^6 (am)^2} \frac{v^2}{1 - v^2} \approx \frac{112\eta_{\text{gr}}^2}{\pi^6 (am)^2} \frac{v^2}{1 - v^2}, \\ \frac{\Delta_v \mathcal{P}_2^{gg}}{\mathcal{P}_{\text{id}}} \approx \frac{560\eta_{\text{gr}}^2 (1 - v_F^2)^2}{3\pi^6 (am)^2} \frac{v^2}{1 - v^2} \approx \frac{560\eta_{\text{gr}}^2}{3\pi^6 (am)^2} \frac{v^2}{1 - v^2}. \quad (68)$$

and for small separation, $(am) \to 0$ and small velocity $v \ll v_F$

$$\frac{\Delta_v \mathcal{E}_2^{gg}}{\mathcal{E}_{\text{id}}} \approx C_e v^2, \quad \frac{\Delta_v \mathcal{P}_2^{gg}}{\mathcal{P}_{\text{id}}} \approx C_p v^2. \quad (69)$$

For $v_F = 1/300$ one has $C_e = C_p = 1.396$. Therefore, $\mathcal{E}_2^{gg} \sim a^{-5}$ for large separations and $\mathcal{E}_2^{gg} \sim a^{-3}$ for small ones. In the non-relativistic case $v \ll v_F$ one has

$$\frac{\Delta_v \mathcal{E}_2^{gg}}{\mathcal{E}_0^{gg}} = E(am) v^2, \quad \frac{\Delta_v \mathcal{P}_2^{gg}}{\mathcal{P}_0^{gg}} = P(am) v^2. \quad (70)$$

Due to Eqs. (68), (69), the function $E(am)$ has the following limits

$$E(ma) = \begin{cases} \frac{7}{18}, & ma \to \infty \\ 290, & ma \to 0 \end{cases}. \quad (71)$$

The function $P$ has the same values for the limit cases.

Numerical evaluations of $E(ma)$ are shown in Fig. 3. The contribution $\mathcal{E}_2$ to the Casimir energy $\mathcal{E}$ is suppressed by factor $(v/C)^2$.

The numerical evaluation of exact expression (41) is shown in Fig. (4). We observe dependence $\sim v^2$ for domain $v \ll v_F$. In the ultrarelativistic case, the energy is divergent. The energy has a maximum at specific velocity $v_{cr}^{gg}$. The numerical evaluations reveal the following dependence

$$v_{cr}^{gg} - v_F \approx \frac{1}{2}(ma), \quad (72)$$

for $ma < 1$. The velocity Fermi defines the position of the $v_{cr}^{gg}$. In the gapless case $v_{cr}^{gg} = v_F$.

### B. Graphene and ideal metal

In this case, we have to take limit $\eta_{\text{gr}} \to \infty$ for graphene at rest in Eq. (65). In this limit $r_{\text{tm}} = 1$ and $r_{\text{te}} = -1$. Therefore, the zero-velocity term reads [9]

$$\mathcal{E}_0^{gm} = \int_0^\infty \frac{\kappa^2 d\kappa}{(2\pi)^2} \int_0^1 dx$$
<figure><img src="image_4.png" /><figcaption>FIG. 3. The functions E and F (70) describe the velocity correction. Here we used $v_F = 1/300$. The limits $a \to 0$ and $a \to \infty$ are in accordance with Eq. (71).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 4. The numerical evaluation the relative correction $\Delta_v \mathcal{E}_2^{gm}/\mathcal{E}_0^{gm} = (\mathcal{E}_2^{gm} - \mathcal{E}_0^{gm})/\mathcal{E}_0^{gm}$ for Casimir energy (41) for case of two graphene. Here we used $v_F = 1/300$ and $\eta_{gr} = 0.0114$. Left panel: $ma = 0.05$ ($m = 0.1$ (eV), $a = 100$ (nm)). The quadratic behavior (70) occurs for small velocities $v \ll v_F$. Right panel: $ma = 0, 0.05, 0.1, 0.5, 0.9, 2$ (from the bottom upwards). The maxima occurred for critical velocity (72).</figcaption></figure>

$$ \times \{\ln(1 - e^{-2ak}r_{\text{tm}}) + \ln(1 + e^{-2ak}r_{\text{te}})\}. \quad (73) $$

For $a \to \infty$, we obtain

$$ \frac{\mathcal{E}_0^{gm}}{\mathcal{E}_{\text{id}}} \Big|_{a \to \infty} \approx \frac{60\eta_{\text{gr}}(2 + v_F^2)}{\pi^5(am)} \approx \frac{120\eta_{\text{gr}}}{\pi^5(am)}. \quad (74) $$

Therefore, for large separation $\mathcal{E}_0^{gm} \sim \eta_{gr}/a^4 m$ [9].

For $a \to 0$ (see Appendix C) one has numerically $\mathcal{E}_0^{gm}/\mathcal{E}_{\text{id}} = 0.026$ for $v_F = 1/300$ and thus, $\mathcal{E}_0^{gm} \sim 1/a^3$.

The numerical evaluations of exact expression (41) are shown in Fig. 5. We observe that for small velocities the energy is constant, and does not depend on the velocity up to some critical

<figure><img src="image_6.png" /><figcaption>FIG. 5. The numerical evaluation the relative correction $\Delta_v \mathcal{E}_2^{gm}/\mathcal{E}_0^{gm} = (\mathcal{E}_2^{gm} - \mathcal{E}_0^{gm})/\mathcal{E}_0^{gm}$ for Casimir energy (41) for system ideal metal/graphene. Here we used $v_F = 1/300$ and $\eta_{gr} = 0.0114$. Left panel: $ma = 0.05$ ($m = 0.1$ (eV), $a = 100$ (nm)). The correction is zero for $v < v_F$ in accordance with (75). Right panel: $ma = 0, 0.05, 0.1, 0.5, 0.9, 2$ (from the bottom upwards).</figcaption></figure>

velocity $v_{cr}^{gm}$. Numerical analysis reveals the same linear dependence (72) for critical velocity

$$ v_{cr}^{gm} - v_F \approx \frac{1}{2}(ma), \quad (75) $$

for $ma < 1$.

## VI. DISCUSSION AND CONCLUSION

We considered here the Casimir energy for two parallel graphenes with relative velocity $v$ in the framework of scattering theory. To obtain the conductivity tensor of moving graphene in the laboratory frame the Lorentz transformation of the boundary conditions is used, which defines the conductivity tensor in the terms of the polarization tensor. The Casimir energy consists of two parts – real (Sec. IV A) and imaginary (Sec. IV B). The real part is a velocity correction to the usual Casimir energy and the imaginary one is the Casimir friction. For both of these energies, the Fermi velocity plays a crucial role. The Casimir friction is zero up to the Fermi velocity, for $v < v_F$. This fact is the simple consequence of the observation that the origin for the Casimir friction is the mixing between positive and negative frequencies [16] (see the end of Sec. IV B) and was considered in Ref. [15] in the framework of effective action approach. The Casimir friction appears if the Eqs. (49) have corresponding solutions with a non-zero imaginary part of the frequency. It takes place if and only if the modulo of reflection coefficients are greater than one (56). The detailed evaluations of the Casimir friction will be considered in a separate paper.

We considered here in detail the real part of the energy for two systems – graphene/graphene and ideal metal/graphene with a relative velocity $v$. In the first system, the velocity correction to the Casimir energy grows up as $v^2$ (70) and then falls after some maximum value (4). The position of this maximum is defined by the velocity Fermi and has linear dependence with (ma) (72). For the system ideal metal/graphene, the energy correction is zero up to some value which has the same form (75) and is defined by the velocity Fermi (see Fig. 5).
For laboratory experiments, the velocity Fermi $v_F \approx 1/300 \approx 10^6$ m/s and situation $v \ll v_F$ takes place. The function $E$ in the Eq. (70) has maximum value 290 for gapless case, $m = 0$. In this case and for $v = v_F$

$$ \frac{\Delta_v \mathcal{E}_2^{gg}}{\mathcal{E}_0^{gg}} = E(0) v_F^2 = 0.0033, \qquad (76) $$

for system graphene/graphene and it is zero for the ideal metal/graphene configuration.

## ACKNOWLEDGMENTS

NK was supported in part by the grants 2022/08771-5, 2021/10128-0 of São Paulo Research Foundation (FAPESP). One of us (NK) is grateful to M. Bordag, M. Silveirinha and D. Vasilevich for fruitful discussions. MA acknowledges support from the French Agence Nationale pour la Recherche – ANR ("CAT" project).

## Appendix A: Analytical properties of conductivity

To calculate the Casimir energy we have to consider the conductivity as a function of complex frequency $\omega = \omega_r + i\omega_i$. The general property of permittivity [20, §82] is

$$ \varepsilon(-\omega^*) = \varepsilon^*(\omega). \qquad (A1) $$

For conductors $\varepsilon(\omega) = I + 2i\eta(\omega)/\omega$, and

$$ \eta(-\omega^*) = \eta^*(\omega). \qquad (A2) $$

A graphene at rest has the following conductivity

$$ \eta_1 = \eta_{gr} \frac{\tilde{k}}{\omega} \left( \mathbf{I} + v_F^2 \frac{\mathbf{k} \otimes \mathbf{k}}{\tilde{k}^2} \right) \tilde{\Phi} \left( \frac{\tilde{k}}{2m} \right), \qquad (A3) $$

where $\eta_{gr} = 2\pi\sigma_{gr} = \pi\alpha/2 = 0.0114$, and

$$ \tilde{\Phi}(y) = \frac{2i}{\pi y} \left\{ 1 - \frac{y^2 + 1}{y} \operatorname{arctanh} y \right\}, \\ \tilde{k} = \sqrt{\omega^2 - v_F^2 \mathbf{k}^2}. $$

Obviously, that $y(-\omega^*) = y(\omega^*) = y^*(\omega)$, and therefore

$$ \operatorname{arctanh} y(-\omega^*) = \operatorname{arctanh} y^*(\omega) = \operatorname{arctanh}^* y(\omega), \qquad (A4) $$

and

$$ \tilde{\Phi}(y(-\omega^*)) = -\tilde{\Phi}^*(y(\omega)). $$

Therefore, the relation (A2) is valid: $\eta_1(-\omega^*) = \eta_1^*(\omega)$. Then, it is easy to verify that

$$ \eta_1(\omega^*) = -\eta_1^*(\omega). \qquad (A5) $$

At the imaginary axis $\omega = i\xi$ we obtain

$$ \eta_1(i\xi) = \eta_1^*(i\xi), \quad \eta_1(-i\xi) = -\eta_1^*(i\xi) = \eta_1(i\xi), \qquad (A6) $$

and conductivity is real function

$$ \eta_1(i\xi) = \eta_{gr} \frac{k_F}{i\xi} \left( \mathbf{I} - v_F^2 \frac{\mathbf{k} \otimes \mathbf{k}}{k_F^2} \right) \Phi \left( \frac{k_F}{2m} \right), \qquad (A7) $$

where

$$ \Phi(y) = \frac{2}{\pi y} \left\{ 1 + \frac{y^2 - 1}{y} \operatorname{arctan} y \right\}, \\ k_F = \sqrt{\xi^2 + v_F^2 \mathbf{k}^2}. $$

The conditions (A6) for function $\tilde{\eta}_1(\xi) = \eta_1(i\xi)$ read

$$ \tilde{\eta}_1(\xi) = \tilde{\eta}_1(\xi)^*, \quad \tilde{\eta}_1(-\xi) = -\tilde{\eta}_1(\xi)^* = -\tilde{\eta}_1(\xi). \qquad (A8) $$

The conductivity of the moving graphene has form (23)

$$ \eta_2 = \frac{\omega'}{\omega} \mathbf{G} \eta_2' \mathbf{G}^T, \qquad (A9) $$

where

$$ \mathbf{G} = \mathbf{I} - \frac{\omega}{(ku)} \frac{\mathbf{u} \otimes \mathbf{u}}{u_0 + 1} + \frac{\mathbf{u} \otimes \mathbf{k}}{(ku)}, \qquad (A10) $$

and

$$ \eta_2' = \eta_{gr} \frac{\tilde{k}'}{\omega'} \left( \mathbf{I} + v_F^2 \frac{\mathbf{k}' \otimes \mathbf{k}'}{\tilde{k}'^2} \right) \tilde{\Phi} \left( \frac{\tilde{k}'}{2m} \right). \qquad (A11) $$

Here,

$$ \tilde{k}' = \sqrt{\omega'^2 - v_F^2 \tilde{k}'^2} \\ = \sqrt{\tilde{k}^2 + \gamma^2(1 - v_F^2)(\omega^2 v^2 + (kv)^2 - 2\omega(kv))}, \qquad (A12) $$

and $\omega' = (ku) = \gamma(\omega - (kv))$.

We observe now that

$$ \tilde{k}'(-\omega^*, v) = \tilde{k}'^*(\omega, -v), \quad \omega'(-\omega^*, v) = -\omega'^*(\omega, -v), $$

and $(y' = \tilde{k}'/2m)$

$$ \tilde{\Phi}(y'(-\omega^*, v)) = -\tilde{\Phi}^*(y'(\omega, -v)). $$

Therefore,

$$ \eta_2'(-\omega^*, v) = \eta_2'^*(\omega, -v). \qquad (A13) $$

It is easy to see from (A10) that,

$$ \mathbf{G}(-\omega^*, v) = \mathbf{G}'^*(\omega, -v) $$

and therefore,

$$ \eta_2(-\omega^*, v) = \eta_2^*(\omega, -v). \qquad (A14) $$
The conductivity of the moving graphene obeys the relation (A2) but with a corresponding changing direction of the velocity. It is similar to the case with an external magnetic field.

Then, because

$$ 
\begin{aligned}
\tilde{k}'(\omega^*, \mathbf{v}) &= \tilde{k}'^*(-\omega, -\mathbf{v}) = \tilde{k}'^*(\omega, \mathbf{v}), \\
\omega'(\omega^*, \mathbf{v}) &= -\omega'^*(-\omega, -\mathbf{v}) = \omega'^*(\omega, \mathbf{v}),
\end{aligned}
$$

we have

$$ \tilde{\Phi}(y'(\omega^*, v)) = -\tilde{\Phi}^*(y'(\omega, v)), $$

and therefore,

$$ 
\begin{aligned}
\eta_2'(\omega^*, \mathbf{v}) &= \eta_2'^*(-\omega, -\mathbf{v}) = -\eta_2'^*(\omega, \mathbf{v}), \\
\eta_2(\omega^*, \mathbf{v}) &= \eta_2^*(-\omega, -\mathbf{v}) = -\eta_2^*(\omega, \mathbf{v}).
\end{aligned}
$$

At the imaginary axis $\omega = i\xi$, the conductivity $\tilde{\eta}_2(\xi, \mathbf{v}) = \eta(i\xi, \mathbf{v})$ is not real, but it obeys to relations

$$ \tilde{\eta}_2(\xi, \mathbf{v}) = \tilde{\eta}_2^*(\xi, -\mathbf{v}), \quad \tilde{\eta}_2(-\xi, \mathbf{v}) = -\tilde{\eta}_2^*(\xi, \mathbf{v}). \quad (A15) $$

The scattering matrix (33) at the imaginary axes reads

$$ r_i(ik_E) = - \frac{\xi^2 \tilde{\eta}_i + k \otimes (k \tilde{\eta}_i) + \mathbf{I} \xi k_E \det \tilde{\eta}_i}{\xi^2 \mathrm{tr} \tilde{\eta}_i + (k \tilde{\eta}_i) + \xi k_E (1 + \det \tilde{\eta}_i)}. \quad (A16) $$

Due to relations (A8) and (A15) we obtain

$$ 
\begin{aligned}
r_1(ik_E, -\xi) &= r_1(ik_E, \xi) = r_1^*(ik_E, \xi) \\
r_2(ik_E, -\xi, \mathbf{v}) &= r_1^*(ik_E, \xi, \mathbf{v}).
\end{aligned}
\quad (A17) $$

Let us consider the energy (41)

$$ \mathcal{E}_2 = \iint \frac{d^2 k}{2(2\pi)^3} \int_{-\infty}^{+\infty} d\xi L(\xi, \mathbf{v}), \quad (A18) $$

where

$$ L(\xi, \mathbf{v}) = \ln \det \left[ 1 - e^{-2ak_E} \mathbf{r}_1'(ik_E) \mathbf{r}_2(ik_E) \right]. $$

Taking into account the relations (A17) we obtain, $L(-\xi, \mathbf{v}) = L^*(\xi, \mathbf{v})$, and therefore,

$$ 
\begin{aligned}
\mathcal{E}_2 &= \iint \frac{d^2 k}{2(2\pi)^3} \int_{-\infty}^{\infty} d\xi L(\xi, \mathbf{v}) \\
&= \iint \frac{d^2 k}{(2\pi)^3} \int_0^{\infty} d\xi \Re L(\xi, \mathbf{v}).
\end{aligned}
\quad (A19) $$

Therefore, the $\mathcal{E}_2$ is a real quantity and does not give a contribution to the Casimir friction.

## Appendix B: Direct calculation of the boosted PT

The PT in the frame $K'$ has the following form

$$ \tilde{\Pi}^{ab}(k') = ie^2 \iiint \frac{d^3 p'}{(2\pi)^3} \mathrm{tr}[\tilde{\gamma}^a \mathcal{S}'(p' + k') \tilde{\gamma}^b \mathcal{S}'(p')], \quad (B1) $$

where

$$ \mathcal{S}'(p') = \frac{\tilde{\gamma}^a p'_a + m}{(\tilde{\gamma}^b p'_b)^2 - m^2}, \quad (B2) $$

and $\tilde{\gamma}^a = \tilde{g}_b^a \gamma^b$, $\tilde{g}_b^a = \mathrm{diag}(1, v_F, v_F)$. In the frame $K$ one has the same form of PT

$$ \tilde{\Pi}^{ab}(k) = ie^2 \iiint \frac{d^3 p}{(2\pi)^3} \mathrm{tr}[\tilde{\Gamma}^a \mathcal{S}(p + k) \tilde{\Gamma}^b \mathcal{S}(p)], \quad (B3) $$

but with matrices $\tilde{\Gamma}^a$ and the spinor Green function $\mathcal{S}$.

To obtain these matrices and Green function we make boost $x^a \to x'^a = \Lambda_b^{a'} x^b (x' = \Lambda x)$. Then the bispinors $\psi(x), \bar{\psi}(x)$ in new coordinates read $\psi'(x') = S\psi(x), \bar{\psi}'(x') = \bar{\psi}(x)S^{-1}$ (see, for example, [25]). Therefore,

$$ \tilde{\Gamma}^{a'} = S \tilde{\gamma}^b S^{-1} \Lambda^{a'}_b, \quad \mathcal{S} = \frac{\tilde{\Gamma}^a p_a + m}{(\tilde{\Gamma}^b p_b)^2 - m^2}. \quad (B4) $$

Obviously, that $\Gamma^{a'} = \gamma^{a'}$. By straightforward calculations we obtain

$$ 
\begin{aligned}
\tilde{\Gamma}^a &= v_F \gamma^a + u^a (1 - v_F) \psi, \\
\mathcal{S} &= \frac{v_F \mathbf{p} + (1 - v_F)(up) \psi + m}{p^2 v_F^2 + (up)^2 (1 - v_F^2) - m^2},
\end{aligned}
\quad (B5) $$

where $\psi = \gamma^a u_a$ and $(up) = u^a p_a$. Taking into account the above expressions we arrive at the following expression of polarization tensor (B3)

$$ \tilde{\Pi}^{\mu\nu}(k) = \frac{ie^2}{2\pi^3} \iiint \frac{t^{\mu\nu} d^3 p}{Q(q)Q(p)}, \quad (B6) $$

where

$$ 
\begin{aligned}
t^{ab} &= v_F^4 p^{(a} q^{b)} + v_F^2 (1 - v_F^2) (up) p^{(a} u^{b)} + v_F^2 (1 - v_F^2) (up) q^{(a} u^{b)} \\
&\quad + (1 - v_F^2) [m^2 - (qp) v_F^2 + (up)(uq)(1 - v_F^2)] u^a u^b \\
&\quad + v_F^2 [m^2 - (qp) v_F^2 - (up)(uq)(1 - v_F^2)] g^{ab}.
\end{aligned}
\quad (B7) $$

Here, $q = p + k$, $Q(p) = p^2 v_F^2 + (up)^2(1 - v_F^2) - m^2$, and $a^{(ab)} = a^{ab} + a^{ba}$. Then, taking into account the Pauli-Villars renormalization we arrive at (14).

## Appendix C: Zero velocity terms for short distance

The direct evaluations integrals in Eqs. (65), (73) where $a \to 0$ has taken, give the following expressions for zero-velocity terms
$$
\begin{align}
\left. \frac{\mathcal{E}_{\text{tm},0}^{\text{gg}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} &= \frac{45}{\pi^4} \frac{\eta_{\text{gr}}^2}{\eta_{\text{gr}}^2 - v_F^2} \left\{ \frac{\eta_{\text{gr}}}{\eta_{\text{gr}} + 1} - \frac{v_F^2 \operatorname{arctanh}\left(\frac{\sqrt{1-v_F^2}\sqrt{\eta_{\text{gr}}^2 - v_F^2}}{v_F^2 + \eta_{\text{gr}}}\right)}{\sqrt{1-v_F^2}\sqrt{\eta_{\text{gr}}^2 - v_F^2}} \right\} = 0.0047, \nonumber \\
\left. \frac{\mathcal{E}_{\text{te},0}^{\text{gg}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} &= \frac{45}{\pi^4 \sqrt{1-v_F^2}\eta_{\text{gr}}(1-v_F^2\eta_{\text{gr}}^2)} \left\{ \frac{2-3v_F^2\eta_{\text{gr}}^2}{\sqrt{1-v_F^2}\eta_{\text{gr}}^2} \operatorname{arctanh}\left(\frac{\sqrt{1-v_F^2}\sqrt{1-v_F^2\eta_{\text{gr}}^2}}{v_F^2\eta_{\text{gr}} + 1}\right) - 2(1-v_F^2\eta_{\text{gr}}^2)\operatorname{arctanh}\left(\sqrt{1-v_F^2}\right) \right. \nonumber \\
&\qquad \left. + \frac{\sqrt{1-v_F^2}\eta_{\text{gr}}}{\eta_{\text{gr}} + 1}(2+\eta_{\text{gr}}-v_F^2\eta_{\text{gr}}^2(\eta_{\text{gr}}+1)) \right\} = 0.000019, \tag{C1}
\end{align}
$$

for system of two graphene, and

$$
\begin{align}
\left. \frac{\mathcal{E}_{\text{tm},0}^{\text{gm}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} &= \frac{45 \eta_{\text{gr}}}{\pi^4} \left\{ \frac{\operatorname{arctanh}\left(\sqrt{1-v_F^2}\right)}{\sqrt{1-v_F^2}} - \frac{\eta_{\text{gr}} \operatorname{arctanh}\left(\frac{\sqrt{1-v_F^2}\sqrt{\eta_{\text{gr}}^2 - v_F^2}}{v_F^2 + \eta_{\text{gr}}}\right)}{\sqrt{1-v_F^2}\sqrt{\eta_{\text{gr}}^2 - v_F^2}} \right\} = 0.023, \nonumber \\
\left. \frac{\mathcal{E}_{\text{te},0}^{\text{gm}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} &= \frac{45}{\pi^4 \eta_{\text{gr}}} \left\{ \eta_{\text{gr}} - \frac{\operatorname{arctanh}\left(\sqrt{1-v_F^2}\right)}{\sqrt{1-v_F^2}} + \frac{\operatorname{arctanh}\left(\frac{\sqrt{1-v_F^2}\sqrt{1-v_F^2\eta_{\text{gr}}^2}}{v_F^2\eta_{\text{gr}}+1}\right)}{\sqrt{1-v_F^2}\sqrt{1-v_F^2\eta_{\text{gr}}^2}} \right\} = 0.0026, \tag{C2}
\end{align}
$$

for graphene/ideal metal configuration. Numerical evaluation
has been made for $v_F = 1/300$. For two graphene, the contribu-
tion of TE mode is 3 orders smaller than TM mode, while for the
second case, one has single order, only. Because $\eta_{\text{gr}}/v_F = 3.4$ we

can not expand this expression over $\eta_{\text{gr}}$. For $v_F = 0$ we obtain

$$
\left. \frac{\mathcal{E}_{\text{tm},0}^{\text{gg}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} \approx \frac{45\eta_{\text{gr}}}{\pi^4}, \quad \left. \frac{\mathcal{E}_{\text{te},0}^{\text{gg}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} \approx \frac{15\eta_{\text{gr}}^2}{\pi^4}, \qquad (\text{C3})
$$

for two graphene, and

$$
\left. \frac{\mathcal{E}_{\text{tm},0}^{\text{gm}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} \approx -\frac{45\eta_{\text{gr}}\ln\eta_{\text{gr}}}{\pi^4}, \quad \left. \frac{\mathcal{E}_{\text{te},0}^{\text{gm}}}{\mathcal{E}_{\text{id}}} \right|_{a \to 0} \approx \frac{45\eta_{\text{gr}}}{2\pi^4}, \qquad (\text{C4})
$$

for second case.

[1] H. B. G. Casimir, Kon. Ned. Akad. Wetensch. Proc. 51, 793-795 (1948); E. M. Lifshitz, Sov. Phys. JETP 2, 73 (1956); I. Dzyaloshinskii, E. Lifshitz, and L. Pitaevskii, Adv. Phys. 10, 165-209 (1961).

[2] L. M. Woods, D. A. R. Dalvit, A. Tkatchenko, P. Rodriguez-Lopez, A. W. Rodriguez, and R. Podgornik, Rev. Mod. Phys. 88, 45003 (2016).

[3] M. Bordag, G. L. Klimchitskaya, U. Mohideen, and V. M. Mostepanenko, *Advances in the Casimir Effect* (2009) p. 768.

[4] M. Antezza, L. P. Pitaevskii, and S. Stringari, Phys. Rev. A 70, 053619 (2004); Phys. Rev. Lett. 95, 113202 (2005); M. Antezza, J. Phys. A. Math. Gen. 39, 6117-6126 (2006); J. M. Obrecht, R. J. Wild, M. Antezza, L. P. Pitaevskii, S. Stringari, and E. A. Cornell, Phys. Rev. Lett.

98, 063201 (2007); M. Antezza, L. P. Pitaevskii, S. Stringari, and V. B. Svetovoy, Phys. Rev. A 77, 022901 (2008).

[5] G. L. Klimchitskaya, U. Mohideen, and V. M. Mostepanenko, Int. J. Mod. Phys. A 37, 22410032 (2022); C. Abbas, B. Guizal, and M. Antezza, Phys. Rev. Lett. 118, 126101 (2017).

[6] J. H. Wilson, A. A. Allocca, and V. Galitski, Phys. Rev. B 91, 235115 (2015); P. Rodriguez-Lopez, A. Popescu, I. Fialkovsky, N. Khusnutdinov, and L. M. Woods, Commun. Mater. 1, 14 (2020); M. B. Farias, A. A. Zyuzin, and T. L. Schmidt, Phys. Rev. B 101, 235446 (2020).

[7] B.-S. Lu, Universe 7, 237 (2021).

[8] V. P. Gusynin and S. G. Sharapov, Phys. Rev. B 73, 245411 (2006).
[9] M. Bordag, I. V. Fialkovsky, D. M. Gitman, and D. V. Vassilevich, Phys. Rev. B **80**, 245406 (2009); I. V. Fialkovsky, V. N. Marachevsky, and D. V. Vassilevich, Phys. Rev. B **84**, 035446 (2011).

[10] M. Bordag, G. Petrov, and D. Robaschik, Sov. J. Nucl. Phys. **39**, 828 (1984); M. Bordag, F. M. Dittes, and D. Robaschik, **43**, 1034 (1986).

[11] V. Dodonov, Physics **2**, 67–104 (2020).

[12] S. I. Maslovski, Phys. Rev. A **84**, 022506 (2011).

[13] J. B. Pendry, J. Phys. Condens. Matter **9**, 10301–10320 (1997); T. G. Philbin and U. Leonhardt, New J. Phys. **11**, 033035 (2009); I. Brevik, B. Shapiro, and M. G. Silveirinha, Int. J. Mod. Phys. A **37**, 22410123 (2022).

[14] A. I. Volokitin and B. N. J. Persson, Rev. Mod. Phys. **79**, 1291–1329 (2007).

[15] M. B. Farias, C. D. Fosco, F. C. Lombardo, and F. D. Mazzitelli, Phys. Rev. D **95**, 065012 (2017); M. B. Farías, C. D. Fosco, F. C. Lombardo, F. D. Mazzitelli, and A. E. R. López, Phys. Rev. D **91**, 105020 (2015).

[16] M. F. Maghrebi, R. Golestanian, and M. Kardar, Phys. Rev. A **88**, 042509 (2013).

[17] I. Fialkovsky, N. Khusnutdinov, and D. Vassilevich, Phys. Rev. B **97**, 165432 (2018).

[18] N. V. Kampen, B. Nijboer, and K. Schram, Phys. Lett. A **26**, 307–308 (1968); E. Gerlach, Phys. Rev. B **4**, 393–396 (1971).

[19] C. Henkel, K. Joulain, J.-P. Mulet, and J.-J. Greffet, Phys. Rev. A **69**, 023808 (2004).

[20] L. D. Landau and E. M. Lifshitz, *Electrodynamics of Continuous Media*, 2nd ed., Vol. Volume 8 (Pergamon, 1984).

[21] K. Joulain, J.-P. Mulet, F. Marquier, R. Carminati, and J.-J. Greffet, Surf. Sci. Rep. **57**, 59–112 (2005).

[22] M. Bordag, Phys. Rev. D **85**, 025005 (2012).

[23] M. G. Silveirinha, Phys. Rev. X **4**, 031013 (2014); M. F. Maghrebi, R. Golestanian, and M. Kardar, Phys. Rev. D **87**, 025016 (2013).

[24] N. Khusnutdinov, D. Drosdoff, and L. M. Woods, Phys. Rev. D **89**, 85033 (2014).

[25] C. Itzykson and J.-B. Zuber, *Quantum Field Theory* (Dover Publications Inc., 2006).