# Casimir effect in a Lorentz-violating tensor extension of a scalar field
theory

M. C. Araújo<EMAIL_ADDRESS>Universidade Federal do Ceará (UFC),
Departamento de Física,
Campus do Pici, Fortaleza - CE, C.P. 6030, 60455-760 - Brazil.  J. Furtado
<EMAIL_ADDRESS>Universidade Federal do Cariri (UFCA), Av. Tenente
Raimundo Rocha,
Cidade Universitária, Juazeiro do Norte, Ceará, CEP 63048-080, Brasil  R. V.
Maluf<EMAIL_ADDRESS>Universidade Federal do Ceará (UFC),
Departamento de Física,
Campus do Pici, Fortaleza - CE, C.P. 6030, 60455-760 - Brazil.  Departamento
de Física Teórica and IFIC, Centro Mixto Universidad de Valencia - CSIC.
Universidad de Valencia, Burjassot-46100, Valencia, Spain.

(August 28, 2024)

###### Abstract

This paper investigates the Casimir Energy modifications due to the Lorentz-
violating CPT-even contribution in an extension of the scalar QED. We have
considered the complex scalar field satisfying Dirichlet boundary conditions
between two parallel plates separated by a small distance. An appropriate
tensor parametrization allowed us to study the Casimir effect in three setups:
isotropic, anisotropic parity-odd, and anisotropic parity-even. We have shown
that the Lorentz-violating contributions promote increased Casimir energy for
both the isotropic and anisotropic parity-odd configurations. However, in the
parity-even case, the Lorentz-violating terms can promote either an increase
or a decrease in the Casimir energy. We have shown that both the increased and
decreased amounts in the Casimir energy depend on the momentum projection over
the Lorentz-violating vectors.

Casimir effect, Lorentz symmetry breaking, CPT violation, Lorentz-violating
scalar QED.

##  I Introduction

The Casimir effect is a purely quantum phenomenon and presents as one of the
most direct manifestations of the existence of vacuum quantum fluctuations
Bordag:2001qi ; Bordag:2009zz . Predicted by H. Casimir in 1948 Casimir:1948dh
, the Casimir effect is characterized by the force of attraction taking place
between two parallel, electrically neutral conducting plates separated by a
short distance in the quantum vacuum, and it emerges as a consequence of the
boundary conditions imposed on the electromagnetic field by the plates. These
boundaries force the frequencies of the field to be in a discrete spectrum
where only specific, well-defining values are allowed. Consequently, the zero-
point energy (vacuum energy) also goes through changes. After some
regularization procedure, we obtain a finite amount that can be interpreted as
the needed energy to keep the plates in the desired configuration.
Experimentally, M. J. Sparnaay confirmed the Casimir effect at the micrometer
scale by Al plates in 1958 Sparnaay:1958wg . About thirty-nine years later,
the experiment was again carried out with higher accuracy using layers of Cu
and Au by Lamoreaux Lamoreaux:1996wh  and for a metallic sphere and flat plate
by Mohideen and Roy Mohideen:1998iz . A modern review of the experimental
methods and realistic measurements can be found in the references
Mohideen:1998iz ; Klimchitskaya:2009cw ; Milton:2019izg .

Although the Casimir effect was initially studied for the electromagnetic
field, it is known to occur for any quantum field under certain boundary
conditions. These boundaries can be material means, interfaces between two
phases of the vacuum, or even space-time topologies Bordag:2001qi ;
Plunien:1986ca . In fact, the Casimir force can depend on many parameters: as
the geometry of the confining materials and the type of boundary conditions
(Dirichlet, Neumann, or mixed) Cruz:2017kfo ; Muniz:2015jba , or even the
presence of extra-dimensions Ponton:2001hq . In this sense, the Casimir effect
was widely evaluated in a number of different scenarios: as black holes
Muniz:2015jba ; Alencar:2018avb ; Sorge:2019ecb ; Christodoulakis:2001ps ;
Vagenas:2003tv , modified gravity Lambiase:2016bjy ; Buoninfante:2018bkc ,
Hořava-Lifshitz-like field theory Muniz:2013uva ; Ferrari:2010dj ;
MoralesUlion:2015tve ; Muniz:2014dga ; daSilva:2019iwn , and Lorentz-violation
contexts Cruz:2017kfo ; Blasone:2018nfy ; Escobar-Ruiz:2021dxi ; Cruz:2018thz
.

As is known, Lorentz covariance and CPT symmetry are preserved by the Standard
Model (SM) of particle physics. However, although experimentally successful,
the SM still leaves some open questions, such as the neutrino oscillation
Barger:2012pxa , or the lack of a quantum description of gravity, making the
search for a physics beyond the SM a quite natural choice. In these new
scenarios, Lorentz and CPT violation could arise from an underlying theory
combining gravity with quantum mechanics such as string theory
Kostelecky:1988zi ; Kostelecky:1991ak , Horava-Lifshitz gravity Horava:2009uw
; Cognola:2016gjy , loop quantum gravity Gambini:1998it ; Bojowald:2004bb ,
among others. Alternatively, in order to investigate signs of such violations
at low energies, Colladay and Kostelecký proposed the Standard Model Extension
(SME) Kostelecky:1994rn ; Colladay:1996iz ; Colladay:1998fq ;
Kostelecky:2003fs . As a comprehensive framework for treating CPT violation at
the level of effective field theory, the SME has led to several experimental
investigations Link:2002fg ; BaBar:2007run ; KLOE-2:2013ozx ; KLOE:2010yad ;
D0:2015ycz ; LHCb:2016vdl . Moreover, it has been extensively studied in
numerous contexts such as neutrino oscillation MINOS:2008fnv ; MINOS:2010kat ;
Katori:2012pe ; MiniBooNE:2011pix , radiative corrections Jackiw:1999yp ;
Kostelecky:2001jc ; Ferrari:2021eam ; Furtado:2014cja , supersymmetric models
Lehum:2018tpi ; Ferrari:2017rwk ; Belich:2015qxa ; Colladay:2010tx , and
noncommutative quantum field theories Carroll:2001ws ; Hayakawa:1999yt .

More recently, Edwards and Kostelecký proposed an extension of the SM scalar
sector where small CPT deviations in neutral-meson oscillations could be
evaluated Edwards:2019lfb . Structurally, the extension consists of a general
effective scalar field theory in any spacetime dimension that contains
explicit-perturbative spin-independent Lorentz-violating operators of
arbitrary mass dimension Edwards:2018lsn . The significance of spin-
independency lies on the fact that even for a particle with nonzero spin,
Lorentz-violating effects could be handled as though the particle had zero
spin. This is extremely interesting since the Higgs boson is the only example
of a fundamental spinless particle in the SM.

In this paper, the four-dimensional tensor sector of the Lorentz violating
extension of the scalar electrodynamics proposed in Ref. Edwards:2018lsn  will
be considered in the calculation of the Casimir energy for a massive complex
scalar field satisfying Dirichlet boundary conditions between two large
parallel plates separated by a small distance. As will be shown, by an
appropriate tensor parametrization in terms of vectors comprising the Lorentz-
violating coefficients, we will be able to analyze the Casimir effect in three
different scenarios: isotropic, anisotropic parity-odd, and anisotropic
parity-even scenarios. Consequently, effects of Lorentz violation on the
Casimir energy may be direction-dependent. In all the cases, we apply the
dimensional regularization process to regularize the Casimir energy.

This paper is organized as follows. In Sec. II, the model is presented and the
tensor parameterization is discussed. In Sec. III, the isotropic case is
considered. Sections IV and V are respectively dedicated to the anisotropic
parity-even and anisotropic parity-odd sectors. In Sec. VI, the main results
of the paper are summarized.

##  II The model

In this section, we present some features concerning the model we are going to
investigate the Casimir effect. As proposed by Edwards and Kostelecký in
reference Edwards:2018lsn , the Lagrange density describing the propagation of
a complex scalar field $\phi$ of mass $m$ in the presence of arbitrary
Lorentz-violating effects can be written in the form

$\displaystyle\mathcal{L}=\partial_{\mu}\phi^{\ast}\partial^{\mu}\phi+\hat{K}_{c}^{\mu\nu}\partial_{\mu}\phi^{\ast}\partial_{\nu}\phi-\frac{i}{2}\left[\hat{K}^{\mu}_{a}(\phi^{\ast}\partial_{\mu}\phi-\partial_{\mu}\phi^{\ast}\phi)\right]-m^{2}\phi^{\ast}\phi,$
(1)

where $\hat{K}_{c}^{\mu\nu}$ and $\hat{K}^{\mu}_{a}$ are respectively a null
trace constant tensor and a constant vector assumed to introduce only small
Lorentz violation to conventional physics. In this work, our focus will be on
the tensor sector of the model (1). General vector contributions will be left
as perspective. Moreover, Lorentz violation effects on the Casimir energy due
to similar vector terms in a massive real scalar field model were investigated
in Ref. Cruz:2017kfo .

The corresponding modified Klein-Gordon equation is

$\displaystyle(\Box+m^{2}+\hat{K}_{c}^{\mu\nu}\partial_{\mu}\partial_{\nu})\phi=0,$
(2)

with similar equation for $\phi^{\ast}$.

Now, since $\hat{K}_{c}^{\mu\nu}$ is a symmetric traceless tensor, we can use
the general parametrization

$\displaystyle\hat{K}_{c}^{\mu\nu}=\frac{1}{2}(U^{\mu}V^{\nu}+U^{\nu}V^{\mu})-\frac{1}{4}\eta^{\mu\nu}(U\cdot
V),$ (3)

where $U$ and $V$ are two arbitrary 4-vectors which comprise the Lorentz-
violating coefficients Casana:2010nd . This parametrization is useful to
investigate the different configurations of this tensor individually
dosSantos:2018rns : the anisotropic parity-even sector, for example, is
parameterized by two pure spacelike 4-vector $U=(0,\mathbf{u})$ and
$V=(0,\mathbf{v})$, with $\mathbf{u}\cdot\mathbf{v}=0$, whereas the
anisotropic parity-odd sector is represented by $U=(0,\mathbf{u})$ and
$V=(v_{0},\mathbf{0})$. The isotropic sector is in turn recovered by two
timelike 4-vector $U=(u_{0},\mathbf{0})$ and $V=(v_{0},\mathbf{0})$. In the
following sections, we will calculate the vacuum energy for each sector above
considering the complex scalar field satisfying Dirichlet boundary conditions,

$\displaystyle\phi(x,y,0,t)=\phi(x,y,a,t)=0,$ (4)

between two large parallel plates separated by a small distance $a$, see Fig.
(1).

Figure 1: Two parallel plates with area $L^{2}$ separated by a small distance
$a$. We mean here by “small” that $a\ll L$.

##  III The isotropic sector

As discussed before, the tensor $\hat{K}^{\mu\nu}_{c}$ in this case must be
parameterized by two timelike 4-vectors $U=(u_{0},\mathbf{0})$ and
$V=(v_{0},\mathbf{0})$ so that Eq. (2) can be rearranged as

$\displaystyle(\Box+m^{2}+\frac{3}{4}\,u_{0}\,v_{0}\,\partial_{0}^{2}+\frac{1}{4}\,u_{0}\,v_{0}\,\nabla^{2})\phi=0.$
(5)

Adopting the standard procedure Ryder:1985wq , the quantum field can be
expanded as

$\displaystyle\phi(x)=\sqrt{\frac{2}{a}}\sum_{n=1}^{\infty}\int\frac{d^{2}p}{(2\pi)^{2}}\frac{\sin{(n\pi
z/a)}}{2\omega_{n}(\mathbf{p})}\left[a_{n}({\mathbf{p}})e^{-ip\cdot
x}+b^{\dagger}_{n}({\mathbf{p}})e^{ip\cdot x}\right],$ (6)

where $n$ is an integer and the angular frequency written as

$\displaystyle\omega_{n}(\mathbf{p})=\frac{1}{\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\sqrt{\left(1-\frac{1}{4}u_{0}v_{0}\right)\left[\mathbf{p}^{2}+\left(\frac{n\pi}{a}\right)^{2}\right]+m^{2}},$
(7)

with $p\cdot x\equiv\omega_{n}(\mathbf{p})t-\mathbf{p}\cdot\mathbf{x}$, and
$\mathbf{p}=(p_{x},p_{y})$. Note that the conditions (4) are directly
verified.

Now, as a result of canonical quantization, the Hamiltonian can be settled as

$\displaystyle
H=\sum_{n=1}^{\infty}\int\frac{d^{2}p}{(2\pi)^{2}}\frac{1}{2}\left(1+\frac{3}{4}u_{0}v_{0}\right)\left[a^{\dagger}_{n}({\mathbf{p}})a_{n}({\mathbf{p}})+b^{\dagger}_{n}({\mathbf{p}})b_{n}({\mathbf{p}})+\frac{L^{2}}{1+\frac{3}{4}u_{0}v_{0}}2\omega_{n}(\mathbf{p})\right],$
(8)

where $L^{2}$ is the area of each plate and the annihilation and creation
operators $a_{n}({\mathbf{p}})$ and $a^{\dagger}_{n}({\mathbf{p}})$,
respectively, satisfy the algebra

$\displaystyle[a_{n}(\mathbf{p}),a^{\dagger}_{m}(\mathbf{q})]=\frac{1}{1+\frac{3}{4}u_{0}v_{0}}\,(2\pi)^{2}\,2\omega_{n}(\mathbf{p})\,\delta^{2}(\mathbf{p}-\mathbf{q})\,\delta_{nm},$
(9)

with all the other commutation relations involving them being zero. Operators
$b_{n}({\mathbf{p}})$ and $b^{\dagger}_{n}({\mathbf{p}})$ obey a similar
algebra. Therefore, taking the vacuum expectation value of the Hamiltonian
(8), we find

$\displaystyle E_{0\,(iso)}=\frac{L^{2}}{(2\pi)^{2}}\sum_{n=1}^{\infty}\int
d^{2}p\,\frac{1}{\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\sqrt{\left(1-\frac{1}{4}u_{0}v_{0}\right)\left[\mathbf{p}^{2}+\left(\frac{n\pi}{a}\right)^{2}\right]+m^{2}},$
(10)

as the vacuum energy for the massive complex scalar field in the isotropic
sector.

In order to reach a finite quantity of the equation (10), we will invoke the
dimensional regularization procedure milton2001casimir  and write

$\displaystyle
E^{Reg}_{0\,(iso)}=\frac{L^{d}}{(2\pi)^{d}}\sum_{n=1}^{\infty}\int
d^{d}p\,\frac{1}{\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\sqrt{\left(1-\frac{1}{4}u_{0}v_{0}\right)\left[\mathbf{p}^{2}+\left(\frac{n\pi}{a}\right)^{2}\right]+m^{2}},$
(11)

where $d$ is the transverse dimension assumed as a continuous, complex
variable. In this way we can perform the momentum integral by using the
Schwinger proper time representation:

$\displaystyle\frac{1}{a^{z}}=\frac{1}{\Gamma(z)}\int_{0}^{\infty}dt\,t^{z-1}e^{-at}.$
(12)

After performing the moment integration, the energy density becomes

$\displaystyle
E_{0\,(iso)}^{Reg}=\frac{L^{d}}{(2\pi)^{d}\,\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\frac{\pi^{d/2}\,\Gamma\left(-\frac{d+1}{2}\right)}{\Gamma\left(-\frac{1}{2}\right)\left(1-\frac{1}{4}u_{0}v_{0}\right)^{d/2}}\sum_{n=1}^{\infty}\left[m^{2}+\left(1-\frac{1}{4}u_{0}v_{0}\right)\left(\frac{n\pi}{a}\right)^{2}\right]^{\frac{d+1}{2}}.$

One put the sum over $n$ into a more helpful form using an Epstein-Hurwitz
Zeta function type Zhai:2010mr

$\displaystyle\sum_{n=-\infty}^{+\infty}(bn^{2}+\mu^{2})^{-s}=\sqrt{\frac{\pi}{b}}\,\frac{\Gamma\left(s-\frac{1}{2}\right)}{\Gamma(s)}\mu^{1-2s}+\frac{2\pi^{s}}{\sqrt{b}\Gamma(s)}\sum_{n=-\infty}^{+\infty\,\prime}\mu^{\frac{1}{2}-s}\left(\frac{n}{\sqrt{b}}\right)^{s-\frac{1}{2}}K_{\frac{1}{2}-s}\left(\frac{2\pi\mu
n}{\sqrt{b}}\right),$

where $K_{\nu}(z)$ is the modified Bessel function, and the prime in the sum
means that the term $n=0$ has to be excluded. After some algebra, we can write
the dimensionally regularized energy in Eq. (III) as

$\displaystyle E_{0\,(iso)}^{Reg}$ $\displaystyle=$
$\displaystyle-\frac{L^{d}}{(4\pi)^{\frac{d+2}{2}}\left(1-\frac{1}{4}u_{0}v_{0}\right)^{\frac{d+1}{2}}\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\left\\{-m^{d+1}\Gamma\left(-\frac{d+1}{2}\right)\sqrt{\left(1-\frac{1}{4}u_{0}v_{0}\right)\pi}\right.$
(15) $\displaystyle+$
$\displaystyle\left.a\,m^{d+2}\Gamma\left(-\frac{d+2}{2}\right)+\frac{4m^{\frac{d+2}{2}}\left(1-\frac{1}{4}u_{0}v_{0}\right)^{\frac{d+2}{4}}}{a^{d/2}}\sum_{n=1}^{\infty}\frac{K_{\frac{d+2}{2}}\left(\frac{2amn}{\sqrt{1-\frac{1}{4}u_{0}v_{0}}}\right)}{n^{(d+2)/2}}\right\\}.$

Note that only the third term in brackets is relevant for the Casimir energy.
The other two terms do not contribute because they do not depend on the
distance $a$ or because they are linearly related to it. Hence, the Casimir
energy for an arbitrary $d$ can be expressed as

$\displaystyle
E^{Cas}_{(iso)}(d)=-\frac{4\,L^{d}\,a^{-d/2}}{\left(1-\frac{1}{4}u_{0}v_{0}\right)^{-\frac{d}{4}}\sqrt{1+\frac{3}{4}u_{0}v_{0}}}\left(\frac{m}{4\pi}\right)^{\frac{d+2}{2}}\sum_{n=1}^{\infty}\frac{K_{\frac{d+2}{2}}\left(\frac{2amn}{\sqrt{1-\frac{1}{4}u_{0}v_{0}}}\right)}{n^{(d+2)/2}}.$
(16)

Finally, taking $d=2$ in the expression above, some interesting limits can be
analyzed. For example, to the case where
$ma\ll\left(1-\frac{1}{4}u_{0}v_{0}\right)^{1/2}$ we get

$\displaystyle
E_{(iso)}^{Cas}=-\sqrt{\frac{4-u_{0}v_{0}}{4+3u_{0}v_{0}}}\frac{L^{2}\,\pi^{2}}{720\,a^{3}}+\frac{4}{\sqrt{16-8u_{0}v_{0}-3u_{0}^{2}v_{0}^{2}}}\frac{L^{2}\,m^{2}}{48\,a}+\cdots,$
(17)

where the first term can be recognized as the usual Casimir energy for the
massless complex scalar field and the second one as the usual mass term, both
fixed by multiplicative Lorentz-violating factors. A quick inspection of the
result (17) allows us to conclude that the Lorentz-violating coefficient
decreases the magnitude of the Casimir energy in the first term. In contrast,
the opposite occurs in the second term, where the Casimir energy increase in
magnitude. In addition, since the product $u_{0}v_{0}$ is supposed to be very
small, the Eq. (17) can be expanded, up to first order, as

$\displaystyle
E^{Cas}_{(iso)}=-\frac{\pi^{2}L^{2}}{720a^{3}}+\frac{L^{2}m^{2}}{48a}+\frac{\pi^{2}L^{2}u_{0}v_{0}}{1440a^{3}}+\frac{L^{2}m^{2}u_{0}v_{0}}{192a}.$
(18)

The above expression allow us to see that the net effect of the isotropic
Lorentz-violating contribution, expressed by the product $u_{0}v_{0}$, is in
order to increase the Casimir energy. The behaviour of the Casimir energy for
this case is graphically depicted in Fig. (2). As we can see, when the value
of the product $u_{0}v_{0}$ increases, the Casimir energy increases as well.
Therefore, for too small values of $u_{0}v_{0}$, like the current bounds on it
($u_{0}v_{0}\approx 10^{-14}GeV$), the curves are practically overlapping.

Figure 2: Casimir energy for some values of the product $u_{0}v_{0}$. For this
plot we have considered $L=1$, $m=1$ and $a=(10^{-7},10^{-6})$.

On the other extreme, i.e. $ma\gg\left(1-\frac{1}{4}u_{0}v_{0}\right)^{1/2}$,
we have

$\displaystyle
E^{Cas}_{(iso)}=-\frac{1}{\sqrt{\left(1+\frac{3}{4}u_{0}v_{0}\right)\sqrt{1-\frac{1}{4}u_{0}v_{0}}}}\frac{L^{2}\,m^{2}}{8\,\pi^{2}\,a}\left(\frac{\pi}{m\,a}\right)^{1/2}\exp{\left(-\frac{2\,m\,a}{\sqrt{1-\frac{1}{4}u_{0}v_{0}}}\right)},$
(19)

such that the expected exponential decay with respect to the mass is observed,
however shaped by Lorentz-violating factors as well. Note that the Lorentz-
violating coefficient into the exponential factor causes the magnitude of the
Casimir energy to decrease faster when compared to the usual case.
Accordingly, the Casimir energy must again lead to a small force at the non-
relativistic limit. Note that both results (17) and (19) carry an extra factor
of $2$ that takes into account the degrees of freedom of the complex scalar
field. Also, the standard results are recovered when Lorentz violation is
“turned off”. Accordingly, up to first order in $u_{0}v_{0}$, Eq. (19) reads

$\displaystyle
E^{Cas}_{(iso)}=-\frac{L^{2}m^{2}e^{-2am}\sqrt{\frac{1}{am}}}{8\pi^{3/2}a}+\frac{L^{2}m^{2}e^{-2am}\sqrt{\frac{1}{am}}(4am+5)}{128\pi^{3/2}a}u_{0}v_{0},$
(20)

and also here, the Lorentz-violating contribution is strictly positive, thus
promoting an increase in the Casimir energy. As we can see from Fig. (3), the
behaviour of the Casimir energy in this regime is very similar to that
presented in Fig. (2).

Figure 3: Casimir energy for some values of the product $u_{0}v_{0}$. For this
plot we have considered $L=1$, $m=10^{7}$ and $a=(10^{-7},10^{-6})$.

##  IV The anisotropic parity-even sector

The anisotropic parity-even sector is in turn parameterized by two pure
spacelike 4-vector $U=(0,\mathbf{u})$ and $V=(0,\mathbf{v})$, with
$\mathbf{u}\cdot\mathbf{v}=0$. Thus, the modified Klein-Gordon equation takes
the form

$\displaystyle\left[\Box+m^{2}+\left(\mathbf{u}\cdot\nabla\right)\left(\mathbf{v}\cdot\nabla\right)\right]\phi=0.$
(21)

In order to apply the method of separation of variables for solving Eq. (21),
two types of solutions will be considered here. The first one will correspond
to the case where there is no Lorentz violation in the z-direction, and the
second one to the case where only Lorentz violation in the z-direction is
taken into account. Let us treat each case separately.

###  IV.1 No Lorentz violation in the z-direction

The general situation for this case is reached when
$\mathbf{u}=(u_{x},u_{y},0)$ and $\mathbf{v}=(v_{x},v_{y},0)$, so that the
solution to Eq. (21) can again be written like Eq. (6), but now with a new
dispersion relation:

$\displaystyle\omega_{n}(\mathbf{p})=\sqrt{\mathbf{p}^{2}+m^{2}-\left(\mathbf{u}\cdot\mathbf{p}\right)\left(\mathbf{v}\cdot\mathbf{p}\right)+\left(\frac{n\pi}{a}\right)^{2}},$
(22)

where $\mathbf{p}=(p_{x},p_{y})$ as well.

This leads to a new Hamiltonian, that now reads

$\displaystyle
H=\frac{1}{2}\sum_{n=1}^{\infty}\int\frac{d^{2}p}{(2\pi)^{2}}\left[a^{\dagger}_{n}({\mathbf{p}})a_{n}({\mathbf{p}})+b^{\dagger}_{n}({\mathbf{p}})b_{n}({\mathbf{p}})+L^{2}2\omega_{n}(\mathbf{p})\right],$
(23)

where the annihilation and creation operators $a_{n}(\mathbf{p})$ and
$a^{\dagger}_{n}(\mathbf{p})$, respectively, satisfy

$\displaystyle[a_{n}(\mathbf{p}),a^{\dagger}_{m}(\mathbf{q})]=(2\pi)^{2}\,2\omega_{n}(\mathbf{p})\,\delta^{2}(\mathbf{p}-q)\,\delta_{nm}$
(24)

with all the other commutation relations involving them being zero and with
$b_{n}({\mathbf{p}})$ and $b^{\dagger}_{n}({\mathbf{p}})$ obeying a similar
algebra.

Consequently, the vacuum energy for the massive complex scalar field in the
anisotropic parity-even sector with no Lorentz violation in the z-direction is

$\displaystyle E_{(APE)}=\frac{L^{2}}{(2\pi)^{2}}\sum_{n=1}^{\infty}\int
d^{2}p\,\sqrt{\mathbf{p}^{2}+m^{2}-\left(\mathbf{u}\cdot\mathbf{p}\right)\left(\mathbf{v}\cdot\mathbf{p}\right)+\left(\frac{n\pi}{a}\right)^{2}}.$
(25)

Proceeding as before, we can then write the corresponding Casimir energy for
an arbitrary dimension $d$:

$\displaystyle
E^{Cas}_{(APE)}(d)=-\frac{4\,L^{d}\,a^{-d/2}}{\left(1-|\mathbf{u}||\mathbf{v}|\cos{\theta_{u}}\,\cos{\theta_{v}}\right)^{d/2}}\left(\frac{m}{4\pi}\right)^{\frac{d+2}{2}}\sum_{n=1}^{\infty}\frac{K_{\frac{d+2}{2}}\left(2amn\right)}{n^{(d+2)/2}}.$
(26)

Here $\theta_{u}$ and $\theta_{v}$ are respectively the angles that the
vectors $\mathbf{u}$ and $\mathbf{v}$ have with respect to the momentum
$\mathbf{p}$.

For $d=2$, the asymptotic behavior $ma\ll 1$ reads

$\displaystyle
E^{Cas}_{(APE)}=\frac{1}{1-|\mathbf{u}||\mathbf{v}|\cos{\theta_{u}}\,\cos{\theta_{v}}}\left[-\frac{L^{2}\pi^{2}}{720\,a^{3}}+\frac{L^{2}m^{2}}{48\,a}+\cdots\right],$
(27)

where the terms in brackets can be recognized as the usual Casimir energy and
its correction due to the mass. Unlike the isotropic case, where the usual and
mass terms are each corrected by different Lorentz-violating factors, see Eq.
(17), we see in Eq. (27) that a single factor breaking the Lorentz symmetry
fixes both terms in the anisotropic parity-even sector. Note that no Lorentz
violation effect will be observed if one of the vectors $\mathbf{u}$ and
$\mathbf{v}$, or both, are orthogonal to the momentum $\mathbf{p}$ as well. In
Fig. (4), we depict the behaviour of the Casimir energy in terms of the angles
between the momentum and the vectors $\mathbf{u}$ and $\mathbf{v}$. As it can
be noted, the effect of the Lorentz symmetry violation can be either in order
to increase or decrease the Casimir energy. For instance, when
$\theta_{u}=\theta_{v}=0$ we have a non-vanishing Lorentz contribution and the
Casimir energy reaches a minimum value. On the other hand, when $\theta_{u}=0$
and $\theta_{v}=\pi$ we have also a non-vanishing Lorentz violating
contribution, but the Casimir energy reaches a maximum value instead.

Figure 4: Casimir energy in terms of the angles between the momentum and the
vectors $\mathbf{u}$ and $\mathbf{v}$. For this plot we have considered $L=1$,
$m=1$, $|\mathbf{u}|=|\mathbf{v}|=10^{-1}$ and $a=10^{-7}$.

Now, in the limit $ma\gg 1$, although the Eq. (27) has to be replaced by

$\displaystyle
E^{Cas}_{(APE)}=-\frac{1}{1-|\mathbf{u}||\mathbf{v}|\cos{\theta_{u}}\,\cos{\theta_{v}}}\left[\frac{L^{2}\,m^{2}}{8\pi^{2}a}\left(\frac{\pi}{ma}\right)^{1/2}e^{-2ma}\right],$
(28)

we note that some features are preserved. For example, the usual mass
exponential decay is fixed by a Lorentz-violating factor that disappear when
$\mathbf{u}$ and $\mathbf{v}$ are orthogonal to the momentum $\mathbf{p}$. The
general behaviour of the Casimir energy in this particular regime is similar
to the one presented in Fig. (4), being the energy scale the only
modification.

###  IV.2 Lorentz violation only in the z-direction

In this case, let us assume $\mathbf{u}=u_{z}\,\mathbf{\hat{z}}$ and
$\mathbf{v}=v_{z}\,\mathbf{\hat{z}}$. As before, the solution to Eq. (21) is
given by Eq. (6) with a new dispersion relation, namely

$\displaystyle\omega_{n}(\mathbf{p})=\sqrt{\mathbf{p}^{2}+m^{2}+\left(1-u_{z}v_{z}\right)\left(\frac{n\pi}{a}\right)^{2}}.$
(29)

Accordingly, the Hamiltonian is once again given by Eq. (23) with
$\omega_{n}(\mathbf{p})$ given by Eq. (29). Therefore, the vacuum energy for
the massive complex scalar field in the anisotropic parity-even sector with
Lorentz violation only occurring in the z-direction is

$\displaystyle E_{(APE)}=\frac{L^{2}}{(2\pi)^{2}}\sum_{n=1}^{\infty}\int
d^{2}p\,\sqrt{\mathbf{p}^{2}+m^{2}+\left(1-u_{z}v_{z}\right)\left(\frac{n\pi}{a}\right)^{2}}.$
(30)

By using dimensional regularization, we can then write

$\displaystyle
E^{Cas}_{(APE)}(d)=-(1-u_{z}v_{z})^{d/4}\,4L^{d}\,a^{-d/2}\left(\frac{m}{4\pi}\right)^{\frac{d+2}{2}}\sum_{n=1}^{\infty}\frac{K_{\frac{d+2}{2}}\left(\frac{2amn}{\sqrt{1-u_{z}v_{z}}}\right)}{n^{(d+2)/2}}$
(31)

and by taking $d=2$ we can analyze the limit
$ma\ll\left(1-u_{z}v_{z}\right)^{1/2}$. The result is

$\displaystyle
E^{Cas}_{(APE)}=\left(1-u_{z}v_{z}\right)^{3/2}\left(-\frac{L^{2}\pi^{2}}{720a^{3}}\right)+\left(1-u_{z}v_{z}\right)^{1/2}\left(\frac{L^{2}m^{2}}{48a}\right)+\cdots.$
(32)

Similar to the isotropic case, the usual Casimir energy and its mass
correction are both fixed by different factors violating the Lorentz symmetry.
Also here we can expand the previous expression in terms of $u_{z}v_{z}$ up to
first order as a way of better understanding the Lorentz violating
contribution. Such expansion yields

$E^{Cas}_{(APE)}=-\frac{\pi^{2}L^{2}}{720a^{3}}+\frac{L^{2}m^{2}}{48a}+\frac{\pi^{2}L^{2}}{480a^{3}}u_{z}v_{z}-\frac{L^{2}m^{2}}{96a}u_{z}v_{z}.$
(33)

From the above expression we can clearly see that the Lorentz violating
contributions are not necessarily in order to increase or decrease the Casimir
energy, but the increase or decrease depends on some conditions on the
parameters of the system. If $m^{2}<\pi^{2}/5a$ ($m^{2}>\pi^{2}/5a$) the
Lorentz violating contribution will increase (decrease) the Casimir energy.
Also, for a certain defined range of parameters, the increase or decrease in
the Casimir energy may depend on the signal of the product $u_{z}v_{z}$. If
the vectors $\mathbf{u}$ and $\mathbf{v}$ are parallel, then $u_{z}v_{z}>0$.
Otherwise, if the vectors $\mathbf{u}$ and $\mathbf{v}$ are anti-parallel,
then $u_{z}v_{z}<0$. The effects of parallel and anti-parallel $\mathbf{u}$
and $\mathbf{v}$ in the Casimir energy are depicted in Fig. (5).

Figure 5: Casimir energy for parallel and anti-parallel $\mathbf{u}$ and
$\mathbf{v}$ vectors. For this plot we have considered $L=1$ and $m=1$.

Accordingly, in the limit $ma\gg\left(1-u_{z}v_{z}\right)^{1/2}$, we have

$\displaystyle
E^{Cas}_{(APE)}=-\left(1-u_{z}v_{z}\right)^{3/4}\frac{L^{2}m^{2}}{8\pi^{2}a}\left(\frac{\pi}{ma}\right)^{1/2}\exp\left(-\frac{2ma}{\sqrt{1-u_{z}v_{z}}}\right),$
(34)

or

$E^{Cas}_{(APE)}=-\frac{L^{2}m^{2}e^{-2am}\sqrt{\frac{1}{am}}}{8\pi^{3/2}a}+\frac{L^{2}m^{3}e^{-2am}(4am+3)}{32\pi^{3/2}(am)^{3/2}}u_{z}v_{z},$
(35)

up to first order in $u_{z}v_{z}$. In this regime, the Casimir energy
similarly behaves to that one depicted in Fig. (5), but with a different
energy scale.

##  V The anisotropic parity-odd sector

In the anisotropic parity-odd sector, the tensor in Eq. (3) is parameterized
by a spacelike and a timelike $4$-vector so that the modified Klein-Gordon
equation take the form

$\displaystyle[\Box+m^{2}+v_{0}(\mathbf{u}\cdot\nabla)\partial_{0}]\phi=0.$
(36)

Similarly to the previous section, in order to apply the method of separation
of variables for solving Eq. (36), we will be interested in analyzing the case
where there is no Lorentz violation in the z-direction, that is
$\mathbf{u}=(u_{x},u_{y},0)$. The field solution for this situation is again
given by Eq. (6) with

$\displaystyle\omega_{n}(\mathbf{p})=\sqrt{\frac{v_{0}^{2}}{4}(\mathbf{u\cdot\mathbf{p}})^{2}+\mathbf{p}^{2}+m^{2}+\left(\frac{n\pi}{a}\right)^{2}}.$
(37)

The Hamiltonian can then be read as

$\displaystyle
H=\sum_{n=1}^{\infty}\int\frac{d^{2}p}{(2\pi)^{2}}\frac{1}{4\omega_{n}(\mathbf{p})^{2}}\left[2\omega_{n}(\mathbf{p})^{2}+v_{0}\,\omega_{n}(\mathbf{p})\,(\mathbf{u}\cdot\mathbf{p})\right]\left[a^{\dagger}_{n}({\mathbf{p}})a_{n}({\mathbf{p}})+b^{\dagger}_{n}({\mathbf{p}})b_{n}({\mathbf{p}})+L^{2}2\omega_{n}(\mathbf{p})\right],$
(38)

where all the annihilation and creation operators satisfy the same algebra as
in Eq. (24). Therefore, the vacuum energy for the massive complex scalar field
in the anisotropic parity-odd sector with no Lorentz violation effects
occurring in the z direction is

$\displaystyle E_{(APO)}=\frac{L^{2}}{(2\pi)^{2}}\sum_{n=1}^{\infty}\int
d^{2}p\left[\frac{v_{0}}{2}(\mathbf{u}\cdot\mathbf{p})+\sqrt{\frac{v_{0}^{2}}{4}(\mathbf{u\cdot\mathbf{p}})^{2}+\mathbf{p}^{2}+m^{2}+\left(\frac{n\pi}{a}\right)^{2}}\,\right]$
(39)

By using dimensional regularization and proceeding as before, we can then
write the corresponding Casimir energy for an arbitrary dimension d:

$\displaystyle
E^{Cas}_{(APO)}(d)=-4\,L^{d}\,a^{-d/2}\left[\frac{1}{1+\left(\frac{v_{0}|\mathbf{u}|\cos\theta_{u}}{2}\right)^{2}}\right]^{\frac{d}{2}}\left(\frac{m}{4\pi}\right)^{\frac{d+2}{2}}\sum_{n=1}^{\infty}\frac{K_{\frac{d+2}{2}}\left(2amn\right)}{n^{(d+2)/2}}.$
(40)

where $\theta_{u}$ is again the angle that the vector $\mathbf{u}$ has with
the momentum $\mathbf{p}$. For $d=2$, this expression in the limit $ma\ll 1$
takes the form

$\displaystyle
E^{Cas}_{(APO)}=\frac{1}{1+\left(\frac{v_{0}|\mathbf{u}|\cos\theta_{u}}{2}\right)^{2}}\left[-\frac{L^{2}\pi^{2}}{720\,a^{3}}+\frac{L^{2}m^{2}}{48\,a}+\cdots\right].$
(41)

On the other hand, when $ma\gg\left(1-\frac{1}{4}u_{0}v_{0}\right)^{1/2}$, we
get

$\displaystyle
E^{Cas}_{(APO)}=-\frac{L^{2}m^{2}}{8\pi^{2}a}\left(\frac{\pi}{ma}\right)^{1/2}\frac{e^{-2ma}}{1+\left(\frac{v_{0}|\mathbf{u}|\cos\theta_{u}}{2}\right)^{2}}$
(42)

As we can see, similarly to the anisotropic parity-even case with no Lorentz
violation in the z-direction, the anisotropic parity-odd Casimir energy is
just the usual result for the massive complex scalar field fixed by a Lorentz-
violating factor depending on the angle between the vectors $\mathbf{u}$ and
$\mathbf{p}$. Furthermore, no Lorentz violation effect will be observed if
those vectors are orthogonal to each other. The behaviour of the Casimir
energy for the regime $ma\ll 1$ is depicted in Fig. (6), and an entirely
similar behaviour is observed when $ma\gg 1$, being the energy scale the only
difference.

Figure 6: Casimir energy when $ma\ll 1$. For this plot we have considered
$L=1$, $m=1$ and $a=10^{-7}$.

##  VI Conclusion

In this paper we investigate the influence of the Lorentz-violating CPT-even
contributions of the scalar electrodynamics proposed by Kostelecky and Edwards
in the Casimir energy. We have considered the complex scalar field satisfying
Dirichlet boundary conditions between two parallel plates separated by a small
distance. An appropriate tensor parametrization allwed us to study the Casimir
effect in three configurations: isotropic, anisotropic parity-odd and
anisotropic parity-even. In each case we have employed dimensional
regularization in order to regularize the Casimir energy.

For all three configurations we have showed that the Lorentz violating
contributions affect the Casimir energy by promoting an increase or decrease
on it, depending on the configuration of the system itself. Such an increase
or decrease are indeed very small since the current bounds on the Lorentz
violating tensor $\hat{K}_{c}^{\mu\nu}$ are no greater than $10^{-14}GeV$. The
increase (decrease) in the Casimir energy consequently yields to a decrease
(increase) in the Casimir force. The Casimir force from the Lorentz violating
contributions must be at least fourteen magnitude order smaller than the usual
Casimir force.

For the isotropic configuration we have found a general solution which was
studied in two regimes, namely
$ma\ll\left(1-\frac{1}{4}u_{0}v_{0}\right)^{1/2}$ and
$ma\gg\left(1-\frac{1}{4}u_{0}v_{0}\right)^{1/2}$. Both regimes exhibit a
similar behaviour, i.e., the influence of the Lorentz symmetry violation is in
order to promote an increase in the Casimir energy, being the energy scale the
only difference between the two regimes. Analogously, for the anisotropic
parity-odd configuration we have a similar behaviour of the Casimir energy in
the regimes where $ma\ll 1$ and $ma\gg 1$. In both cases the Lorentz-violating
terms contribute by increasing the Casimir energy, but the amount by which the
Casimir is increased depends on the projection of the momentum on the Lorentz-
violating vectors.

Finally, in the parity-even case, the Lorentz symmetry violation can increase
or decrease the Casimir energy. And in this case as well, the amount by which
the Casimir energy is increased or decreased depends on the projection of the
momentum on the Lorentz-violating vectors. Such a direction dependency on the
Casimir energy can be used as a signature to detect small violation in the
Lorentz symmetry. It is important to highlight that our results are a direct
and important generalization of previous results found in the literature in
regard the influence of Lorentz symmetry violation in the Casimir effect.

## Acknowledgments

The authors thank the Fundação Cearense de Apoio ao Desenvolvimento Científico
e Tecnológico (FUNCAP), Grant no. PNE0112-00085.01.00/16 (JF), and the
Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq), Grant
no. 200879/2022-7 (RVM) for financial support. R. V. Maluf acknowledges the
Departament de Física Teòrica de la Universitat de València for the kind
hospitality.

## References

  * (1) M. Bordag, U. Mohideen and V. M. Mostepanenko, Phys. Rept. 353, 1-205 (2001) doi:10.1016/S0370-1573(01)00015-1 [arXiv:quant-ph/0106045 [quant-ph]].
  * (2) M. Bordag, G. L. Klimchitskaya, U. Mohideen and V. M. Mostepanenko, Int. Ser. Monogr. Phys. 145, 1-768 (2009) Oxford University Press, 2009,
  * (3) H. B. G. Casimir, Indag. Math. 10, 261-263 (1948)
  * (4) M. J. Sparnaay, Physica 24, 751-764 (1958) doi:10.1016/S0031-8914(58)80090-7
  * (5) S. K. Lamoreaux, Phys. Rev. Lett. 78, 5-8 (1997) [erratum: Phys. Rev. Lett. 81, 5475-5476 (1998)] doi:10.1103/PhysRevLett.78.5
  * (6) U. Mohideen and A. Roy, Phys. Rev. Lett. 81, 4549-4552 (1998) doi:10.1103/PhysRevLett.81.4549 [arXiv:physics/9805038 [physics]].
  * (7) G. L. Klimchitskaya, U. Mohideen and V. M. Mostepanenko, Rev. Mod. Phys. 81, 1827-1885 (2009) doi:10.1103/RevModPhys.81.1827 [arXiv:0902.4022 [cond-mat.other]].
  * (8) K. Milton and I. Brevik, Symmetry 11, no.2, 201 (2019) doi:10.3390/sym11020201
  * (9) G. Plunien, B. Muller and W. Greiner, Phys. Rept. 134, 87-193 (1986) doi:10.1016/0370-1573(86)90020-7
  * (10) M. B. Cruz, E. R. Bezerra de Mello and A. Y. Petrov, Phys. Rev. D 96, no.4, 045019 (2017) doi:10.1103/PhysRevD.96.045019 [arXiv:1705.03331 [hep-th]].
  * (11) C. R. Muniz, M. O. Tahim, M. S. Cunha and H. S. Vieira, JCAP 01, 006 (2018) doi:10.1088/1475-7516/2018/01/006 [arXiv:1510.06655 [gr-qc]].
  * (12) E. Ponton and E. Poppitz, JHEP 06, 019 (2001) doi:10.1088/1126-6708/2001/06/019 [arXiv:hep-ph/0105021 [hep-ph]].
  * (13) G. Alencar, R. N. Costa Filho, M. S. Cunha and C. R. Muniz, Eur. Phys. J. Plus 135, no.1, 18 (2020) doi:10.1140/epjp/s13360-019-00044-9 [arXiv:1805.01791 [hep-th]].
  * (14) F. Sorge and J. H. Wilson, Phys. Rev. D 100, no.10, 105007 (2019) doi:10.1103/PhysRevD.100.105007 [arXiv:1909.07357 [gr-qc]].
  * (15) T. Christodoulakis, G. A. Diamandis, B. C. Georgalas and E. C. Vagenas, Phys. Rev. D 64, 124022 (2001) doi:10.1103/PhysRevD.64.124022 [arXiv:hep-th/0107049 [hep-th]].
  * (16) E. C. Vagenas, Phys. Rev. D 68, 024015 (2003) doi:10.1103/PhysRevD.68.024015 [arXiv:hep-th/0301149 [hep-th]].
  * (17) G. Lambiase, A. Stabile and A. Stabile, Phys. Rev. D 95, no.8, 084019 (2017) doi:10.1103/PhysRevD.95.084019 [arXiv:1611.06494 [gr-qc]].
  * (18) L. Buoninfante, G. Lambiase, L. Petruzziello and A. Stabile, Eur. Phys. J. C 79, no.1, 41 (2019) doi:10.1140/epjc/s10052-019-6574-1 [arXiv:1811.12261 [gr-qc]].
  * (19) C. R. Muniz, V. B. Bezerra and M. S. Cunha, Phys. Rev. D 88, 104035 (2013) doi:10.1103/PhysRevD.88.104035 [arXiv:1311.2570 [gr-qc]].
  * (20) A. F. Ferrari, H. O. Girotti, M. Gomes, A. Y. Petrov and A. J. da Silva, Mod. Phys. Lett. A 28, 1350052 (2013) doi:10.1142/S0217732313500521 [arXiv:1006.1635 [hep-th]].
  * (21) I. J. Morales Ulion, E. R. Bezerra de Mello and A. Y. Petrov, Int. J. Mod. Phys. A 30, no.36, 1550220 (2015) doi:10.1142/S0217751X15502206 [arXiv:1511.00489 [hep-th]].
  * (22) C. R. Muniz, V. B. Bezerra and M. S. Cunha, Annals Phys. 359, 55-63 (2015) doi:10.1016/j.aop.2015.04.014 [arXiv:1405.5424 [hep-th]].
  * (23) D. R. da Silva, M. B. Cruz and E. R. Bezerra de Mello, Int. J. Mod. Phys. A 34, no.20, 1950107 (2019) doi:10.1142/S0217751X19501070 [arXiv:1905.01295 [hep-th]].
  * (24) M. Blasone, G. Lambiase, L. Petruzziello and A. Stabile, Eur. Phys. J. C 78, no.11, 976 (2018) doi:10.1140/epjc/s10052-018-6464-y [arXiv:1808.04425 [hep-th]].
  * (25) A. M. Escobar-Ruiz, A. Martín-Ruiz, E. C. A. and R. Linares, Int. J. Mod. Phys. A 36, no.23, 2150168 (2021) doi:10.1142/S0217751X21501682 [arXiv:2105.12953 [hep-th]].
  * (26) M. B. Cruz, E. R. Bezerra de Mello and A. Y. Petrov, Phys. Rev. D 99, no.8, 085012 (2019) doi:10.1103/PhysRevD.99.085012 [arXiv:1812.05428 [hep-th]].
  * (27) V. Barger, D. Marfatia and K. Whisnant, Princeton Univ. Pr., 2012, ISBN 978-0-691-12853-5
  * (28) V. A. Kostelecky and S. Samuel, Phys. Rev. D 39, 683 (1989) doi:10.1103/PhysRevD.39.683
  * (29) V. A. Kostelecky and R. Potting, Nucl. Phys. B 359, 545-570 (1991) doi:10.1016/0550-3213(91)90071-5
  * (30) P. Horava, Phys. Rev. D 79, 084008 (2009) doi:10.1103/PhysRevD.79.084008 [arXiv:0901.3775 [hep-th]].
  * (31) G. Cognola, R. Myrzakulov, L. Sebastiani, S. Vagnozzi and S. Zerbini, Class. Quant. Grav. 33, no.22, 225014 (2016) doi:10.1088/0264-9381/33/22/225014 [arXiv:1601.00102 [gr-qc]].
  * (32) R. Gambini and J. Pullin, Phys. Rev. D 59, 124021 (1999) doi:10.1103/PhysRevD.59.124021 [arXiv:gr-qc/9809038 [gr-qc]].
  * (33) M. Bojowald, H. A. Morales-Tecotl and H. Sahlmann, Phys. Rev. D 71, 084012 (2005) doi:10.1103/PhysRevD.71.084012 [arXiv:gr-qc/0411101 [gr-qc]].
  * (34) V. A. Kostelecky and R. Potting, Phys. Rev. D 51, 3923-3935 (1995) doi:10.1103/PhysRevD.51.3923 [arXiv:hep-ph/9501341 [hep-ph]].
  * (35) D. Colladay and V. A. Kostelecky, Phys. Rev. D 55, 6760-6774 (1997) doi:10.1103/PhysRevD.55.6760 [arXiv:hep-ph/9703464 [hep-ph]].
  * (36) D. Colladay and V. A. Kostelecky, Phys. Rev. D 58, 116002 (1998) doi:10.1103/PhysRevD.58.116002 [arXiv:hep-ph/9809521 [hep-ph]].
  * (37) V. A. Kostelecky, Phys. Rev. D 69, 105009 (2004) doi:10.1103/PhysRevD.69.105009 [arXiv:hep-th/0312310 [hep-th]].
  * (38) J. M. Link, M. Reyes, P. M. Yager, J. C. Anjos, I. Bediaga, C. Gobel, J. Magnin, A. Massafferri, J. M. de Miranda and I. M. Pepe, et al. Phys. Lett. B 556, 7-13 (2003) doi:10.1016/S0370-2693(03)00103-5 [arXiv:hep-ex/0208034 [hep-ex]].
  * (39) B. Aubert et al. [BaBar], Phys. Rev. Lett. 100, 131802 (2008) doi:10.1103/PhysRevLett.100.131802 [arXiv:0711.2713 [hep-ex]].
  * (40) D. Babusci et al. [KLOE-2], Phys. Lett. B 730, 89-94 (2014) doi:10.1016/j.physletb.2014.01.026 [arXiv:1312.6818 [hep-ex]].
  * (41) A. Di Domenico et al. [KLOE], Found. Phys. 40, 852-866 (2010) doi:10.1007/s10701-009-9366-x
  * (42) V. M. Abazov et al. [D0], Phys. Rev. Lett. 115, no.16, 161601 (2015) doi:10.1103/PhysRevLett.115.161601 [arXiv:1506.04123 [hep-ex]].
  * (43) R. Aaij et al. [LHCb], Phys. Rev. Lett. 116, no.24, 241601 (2016) doi:10.1103/PhysRevLett.116.241601 [arXiv:1603.04804 [hep-ex]].
  * (44) P. Adamson et al. [MINOS], Phys. Rev. Lett. 101, 151601 (2008) doi:10.1103/PhysRevLett.101.151601 [arXiv:0806.4945 [hep-ex]].
  * (45) P. Adamson et al. [MINOS], Phys. Rev. Lett. 105, 151601 (2010) doi:10.1103/PhysRevLett.105.151601 [arXiv:1007.2791 [hep-ex]].
  * (46) T. Katori [MiniBooNE], Mod. Phys. Lett. A 27, 1230024 (2012) doi:10.1142/S0217732312300248 [arXiv:1206.6915 [hep-ex]].
  * (47) A. A. Aguilar-Arevalo et al. [MiniBooNE], Phys. Lett. B 718, 1303-1308 (2013) doi:10.1016/j.physletb.2012.12.020 [arXiv:1109.3480 [hep-ex]].
  * (48) R. Jackiw and V. A. Kostelecky, Phys. Rev. Lett. 82, 3572-3575 (1999) doi:10.1103/PhysRevLett.82.3572 [arXiv:hep-ph/9901358 [hep-ph]].
  * (49) V. A. Kostelecky, C. D. Lane and A. G. M. Pickering, Phys. Rev. D 65, 056006 (2002) doi:10.1103/PhysRevD.65.056006 [arXiv:hep-th/0111123 [hep-th]].
  * (50) A. F. Ferrari, J. Furtado, J. F. Assunção, T. Mariz and A. Y. Petrov, EPL 136, no.2, 21002 (2021) doi:10.1209/0295-5075/ac419a [arXiv:2109.11901 [hep-th]].
  * (51) J. Furtado and T. Mariz, Phys. Rev. D 89, no.2, 025021 (2014) doi:10.1103/PhysRevD.89.025021 [arXiv:1401.0492 [hep-ph]].
  * (52) A. C. Lehum, T. Mariz, J. R. Nascimento and A. Y. Petrov, Phys. Lett. B 790, 129-133 (2019) doi:10.1016/j.physletb.2019.01.014 [arXiv:1802.01368 [hep-th]].
  * (53) A. F. Ferrari and A. C. Lehum, EPL 122, no.3, 31001 (2018) doi:10.1209/0295-5075/122/31001 [arXiv:1712.02650 [hep-th]].
  * (54) H. Belich, L. D. Bernald, P. Gaete, J. A. Helayël-Neto and F. J. L. Leal, Eur. Phys. J. C 75, no.6, 291 (2015) doi:10.1140/epjc/s10052-015-3447-0 [arXiv:1502.06126 [hep-th]].
  * (55) D. Colladay and P. McDonald, Phys. Rev. D 83, 025021 (2011) doi:10.1103/PhysRevD.83.025021 [arXiv:1010.1781 [hep-ph]].
  * (56) S. M. Carroll, J. A. Harvey, V. A. Kostelecky, C. D. Lane and T. Okamoto, Phys. Rev. Lett. 87, 141601 (2001) doi:10.1103/PhysRevLett.87.141601 [arXiv:hep-th/0105082 [hep-th]].
  * (57) M. Hayakawa, Phys. Lett. B 478, 394-400 (2000) doi:10.1016/S0370-2693(00)00242-2 [arXiv:hep-th/9912094 [hep-th]].
  * (58) B. R. Edwards and V. A. Kostelecký, Phys. Lett. B 795, 620-626 (2019) doi:10.1016/j.physletb.2019.07.012 [arXiv:1907.05206 [hep-ph]].
  * (59) B. R. Edwards and V. A. Kostelecky, Phys. Lett. B 786, 319-326 (2018) doi:10.1016/j.physletb.2018.10.011 [arXiv:1809.05535 [hep-th]].
  * (60) R. Casana, M. M. Ferreira, A. R. Gomes and F. E. P. dos Santos, Phys. Rev. D 82, 125006 (2010) doi:10.1103/PhysRevD.82.125006 [arXiv:1010.2776 [hep-th]].
  * (61) F. E. P. dos Santos and M. M. Ferreira, Symmetry 10, no.8, 302 (2018) doi:10.3390/sym10080302
  * (62) L. H. Ryder, Cambridge University Press, 1996, ISBN 978-0-521-47814-4, 978-1-139-63239-3, 978-0-521-23764-2 doi:10.1017/CBO9780511813900
  * (63) K. A. Milton, The Casimir effect: Physical manifestations of zeropoint energy (World Scientific, River Edge, 2001).
  * (64) X. h. Zhai, X. z. Li and C. J. Feng, Mod. Phys. Lett. A 26, 669-679 (2011) doi:10.1142/S0217732311035110 [arXiv:1008.3020 [hep-th]].


# Extending third quantization with commuting observables: a dissipative spin-
boson model

Luka Medic1,2, Anton Ramšak1,2, and Tomaž Prosen2 $\phantom{}{}^{1}$ Jožef
Stefan Institute, Jamova 39, Ljubljana, Slovenia  $\phantom{}{}^{2}$
University of Ljubljana, Faculty for Mathematics and Physics, Jadranska 19,
Ljubljana, Slovenia

###### Abstract

We consider the spectral and initial value problem for the Lindblad-Gorini-
Kossakowski-Sudarshan master equation describing an open quantum system of
bosons and spins, where the bosonic parts of the Hamiltonian and Lindblad jump
operators are quadratic and linear respectively, while the spins couple to
bosons via mutually commuting spin operators. Needless to say, the spin
degrees of freedom can be replaced by any set of finite-level quantum systems.
A simple, yet non-trivial example of a single open spin-boson model is worked
out in some detail.

##  1 Introduction

Exact solutions of simple but nontrivial models describing characteristic
physical phenomena are of paramount importance for understanding statistical
physics in nonequilibrium. While there is an abundance of such exactly
solvable nonequilibrium models in classical statistical mechanics (see, e.g.,
a review paper [1]), very few explicit analytic approaches are known in the
quantum realm (e.g., review papers [2, 3]).

A broad route to quantum non-equilibrium physics leads via the theory of open
quantum systems, especially in the many-body realm. Considering a large (many-
body) quantum system, one can often describe the dynamics of its (so-called
central) parts within the framework of the so-called Markov approximation,
neglecting the back information flow from the rest of the system (the so-
called environment) to its central part. The differential equation describing
the density matrix of the central system within the Markov approximation and
the rotating-wave approximation [4] is the unique mathematical evolution law
that preserves the hermeticity, positivity, and trace of the density matrix,
called the Lindblad-Gorini-Kossakowski-Sudarshan [5, 6] equation, or Lindblad
equation for short. We note that the many-body Lindblad equation provides a
perfect mathematical platform for preparing engineered quantum states or
quantum phases of matter within cold atom and ion trap setups [7, 8].

Some time ago, one of the authors developed canonical formalism of
quantization over operator spaces for the diagonalization of the generator of
the many-body Lindblad equation – the so-called Liouvillian superoperator, or
completely solving the Lindbladian initial value problem, for a general
quadratic Hamiltonian and a set of Lindblad jump operators which are linear in
canonical creation/annihilation operators. The original proposal for fermionic
systems [9, 10] was later extended to bosonic systems [11] (see also Ref. [12]
for a more abstract discussion of quantization over operator spaces), and
further developed by other authors [13, 14, 15]. Specifically, the latter
reference extended the technique to include quadratic Hermitian jump
operators, allowing the analytical treatment of nonequilibrium phase
transitions [16] and crossovers between ballistic and diffusive as well as
quantum and classical transport [17, 18]. Note, however, that even within the
class of linear jump operators, one can discuss non-trivial critical phenomena
in translationally invariant [19, 20] and so-called boundary-driven systems
[3, 21, 22] (in the latter, the jump processes are confined to the boundaries
of the system).

Experience has shown that the kinds of systems that can be treated efficiently
under the closed system (2nd quantization) formulation can also be treated
efficiently under the open system (3rd quantization) formulation. However,
there is one type of systems with very important applications, e.g., in
quantum optics, that has somehow been left out so far, namely spin-boson
systems. Some aspects of integrability and exact solvability of these systems
in the closed system framework have been extensively discussed in the
literature (see, e.g. [23] for the Rabi model and [24] for the Jaynes-Cummings
and Dicke models). A recent case is the exactly solvable model of an electron
in a driven harmonic oscillator with Rashba coupling [25, 26]. Coupling to a
thermal bath can be treated numerically [27], but the analytical approach to
such problems is precisely the goal of the present work. Here we essentially
provide a small extension of the third quantization method that allows us to
include additional degrees of freedom with finite-dimensional phase space,
provided that these degrees of freedom enter the Hamiltonian and jump
operators only through commuting operators.

##  2 Formal solution

The aim of this paper is to solve the following Lindblad-Gorini-Kossakowski-
Sudarshan master equation for $n$ particles:

$\frac{\textrm{d}\rho}{\textrm{d}t}=\hat{\mathcal{L}}\rho:=-i\left[H,\rho\right]+\sum_{\mu}\left(2L_{\mu}\rho
L_{\mu}^{\dagger}-\left\\{L_{\mu}^{\dagger}L_{\mu},\rho\right\\}\right),$ (1)

where $H$ and $L_{\mu}$ are the Hamiltonian and Lindblad operators,
respectively. $\rho$ is the density operator describing the state of the $n$
particles. The operators $H$ and $L_{\mu}$ are given by:

$\displaystyle H$ $\displaystyle=$
$\displaystyle\underline{{a}}^{\dagger}\cdot\mathbf{H}\underline{{a}}+\underline{{a}}\cdot\mathbf{K}\underline{{a}}+\underline{{a}}^{\dagger}\cdot\bar{\mathbf{K}}\underline{{a}}^{\dagger}+\underline{\sigma}\cdot\mathbf{\Omega}\underline{{a}}+\underline{{a}}^{\dagger}\cdot\mathbf{\Omega}^{\dagger}\underline{\sigma},$
(2) $\displaystyle L_{\mu}$ $\displaystyle=$
$\displaystyle\underline{l}_{\mu}\cdot\underline{{a}}+\underline{k}_{\mu}\cdot\underline{{a}}^{\dagger}+\underline{w}_{\mu}\cdot\underline{\sigma},$
(3)

where $\underline{{a}}$ and $\underline{{a}}^{\dagger}$ are $n$-dimensional
vectors of canonical bosonic annihilation/creation operators, and
$\underline{\sigma}$ is an $n$-dimensional vector of mutually commuting
Hermitian operators, with a finite discrete spectrum. We also assume that
$[a_{j},\sigma_{j^{\prime}}]=0$. The matrix $\mathbf{H}$ is Hermitian
($\mathbf{H}^{\dagger}=\mathbf{H}$), the matrix $\mathbf{K}$ is symmetric
($\mathbf{K}=\mathbf{K}^{T}$), and the matrix $\mathbf{\Omega}$ is arbitrary.
The vectors $\underline{l}_{\mu}$, $\underline{k}_{\mu}$, and
$\underline{w}_{\mu}$ are $n$-dimensional vectors of constants.

Below we follow the notation and formalism developed in Ref. [11] extending
$n$ bosonic degrees of freedom with additional $n$ finite level quantum
systems (e.g. spins where $\sigma_{j}$ can be considered as their
$z-$projections). We introduce operators $\hat{a}_{\nu,j}$ and
$\hat{a}^{\prime}_{\nu,j}$, where $\nu=0,1$ and $j=1,...,n$:

$\displaystyle\hat{a}_{0,j}=\hat{a}_{j}^{L},\qquad\hat{a}_{0,j}^{\prime}=\hat{a}_{j}^{\dagger}{}^{L}-\hat{a}_{j}^{\dagger}{}^{R},$
$\displaystyle\hat{a}_{1,j}=\hat{a}_{j}^{\dagger}{}^{R},\qquad\hat{a}_{1,j}^{\prime}=\hat{a}_{j}^{R}-\hat{a}_{j}^{L}.$
(4)

Here, the superscripts $L$ and $R$ indicate the left- and right-multiplication
maps acting on a vector space of operators

$\hat{b}^{L}|\rho\rangle=|b\rho\rangle,\quad\hat{b}^{R}|\rho\rangle=|\rho
b\rangle,\quad\textrm{where}\quad b\equiv
a_{j},a^{\dagger}_{j},\;\textrm{or}\;\sigma_{j}\,.$ (5)

The operators (4) satisfy the almost-canonical commutation relations:

$\left[\hat{a}_{\nu,j},\hat{a}^{\prime}_{\mu,k}\right]=\delta_{\nu,\mu}\delta_{j,k},\quad\left[\hat{a}_{\nu,j},\hat{a}_{\mu,k}\right]=\left[\hat{a}^{\prime}_{\nu,j},\hat{a}^{\prime}_{\mu,k}\right]=0.$
(6)

In terms of these operators, we rewrite the Liouvillean as follows:

$\displaystyle\kern-40.00006pt\hat{\mathcal{L}}=$ $\displaystyle-$
$\displaystyle
i\hat{H}^{L}+i\hat{H}^{R}+\sum_{\mu}2\hat{L_{\mu}}^{L}\hat{L_{\mu}^{\dagger}}^{R}-\hat{L_{\mu}^{\dagger}}^{L}\hat{L_{\mu}}^{L}-\hat{L_{\mu}}^{R}\hat{L_{\mu}^{\dagger}}^{R}$
(7) $\displaystyle\kern-30.00005pt=$ $\displaystyle-$ $\displaystyle
i\underline{\hat{a}}_{0}^{\prime}\cdot\mathbf{H}\underline{\hat{a}}_{0}+i\underline{\hat{a}}_{1}^{\prime}\cdot\mathbf{H}\underline{\hat{a}}_{1}+i\underline{\hat{a}}_{1}^{\prime}\cdot\mathbf{K}\left(2\underline{\hat{a}}_{0}+\underline{\hat{a}}_{1}^{\prime}\right)-i\underline{\hat{a}}_{0}^{\prime}\cdot\bar{\mathbf{K}}\left(2\underline{\hat{a}}_{1}+\underline{\hat{a}}_{0}^{\prime}\right)$
(8) $\displaystyle-$ $\displaystyle
i\underline{\hat{a}}_{0}^{\prime}\cdot\mathbf{\Omega}^{\dagger}\underline{\sigma}^{L}+i\underline{\sigma}^{R}\cdot\mathbf{\Omega}\underline{\hat{a}}_{1}^{\prime}-i\left(\underline{\sigma}^{L}-\underline{\sigma}^{R}\right)\cdot\mathbf{\Omega}\underline{\hat{a}}_{0}+i\underline{\hat{a}}_{1}\cdot\mathbf{\Omega}^{\dagger}\left(\underline{\sigma}^{R}-\underline{\sigma}^{L}\right)$
$\displaystyle+$
$\displaystyle\underline{\hat{a}}_{0}^{\prime}\cdot\left(\mathbf{N}-\bar{\mathbf{M}}\right)\underline{\hat{a}}_{0}+\underline{\hat{a}}_{1}^{\prime}\cdot\left(\bar{\mathbf{N}}-\mathbf{M}\right)\underline{\hat{a}}_{1}+\underline{\hat{a}}_{0}^{\prime}\cdot\left(\bar{\mathbf{L}}^{T}-\bar{\mathbf{L}}\right)\underline{\hat{a}}_{1}+\underline{\hat{a}}_{1}^{\prime}\cdot\left(\mathbf{L}^{T}-\mathbf{L}\right)\underline{\hat{a}}_{0}$
$\displaystyle-$
$\displaystyle\underline{\hat{a}}_{0}^{\prime}\cdot\bar{\mathbf{L}}\underline{\hat{a}}_{0}-\underline{\hat{a}}_{1}^{\prime}\cdot\mathbf{L}\underline{\hat{a}}_{1}+2\underline{\hat{a}}_{0}^{\prime}\cdot\mathbf{N}\underline{\hat{a}}_{1}^{\prime}+2\underline{\sigma}^{L}\cdot\mathbf{W}\underline{\sigma}^{R}-\underline{\sigma}^{L}\cdot\mathbf{W}\underline{\sigma}^{L}-\underline{\sigma}^{R}\cdot\mathbf{W}\underline{\sigma}^{R}$
$\displaystyle+$
$\displaystyle\underline{\hat{a}}_{0}\cdot\left(\mathbf{\bar{\mathbf{F}}-\mathbf{E}}\right)\left(\underline{\sigma}^{L}-\underline{\sigma}^{R}\right)+\underline{\hat{a}}_{1}\cdot\left(\mathbf{\mathbf{F}-\bar{\mathbf{E}}}\right)\left(\underline{\sigma}^{R}-\underline{\sigma}^{L}\right)-\underline{\hat{a}}_{0}^{\prime}\cdot\bar{\mathbf{E}}\underline{\sigma}^{L}-\underline{\hat{a}}_{1}^{\prime}\cdot\mathbf{E}\underline{\sigma}^{R}$
$\displaystyle-$
$\displaystyle\underline{\hat{a}}_{0}^{\prime}\cdot\mathbf{F}\left(\underline{\sigma}^{L}-2\underline{\sigma}^{R}\right)-\underline{\hat{a}}_{1}^{\prime}\cdot\bar{\mathbf{F}}\left(\underline{\sigma}^{R}-2\underline{\sigma}^{L}\right),$

where we define $n\times n$ matrices:

$\displaystyle\mathbf{M}$ $\displaystyle:=$
$\displaystyle\sum_{\mu}\underline{l}_{\mu}\otimes\underline{\bar{l}}_{\mu}=\mathbf{M}^{\dagger},\quad\mathbf{N}:=\sum_{\mu}\underline{k}_{\mu}\otimes\underline{\bar{k}}_{\mu}=\mathbf{N}^{\dagger},\quad\mathbf{L}:=\sum_{\mu}\underline{l}_{\mu}\otimes\underline{\bar{k}}_{\mu},$
$\displaystyle\mathbf{W}$ $\displaystyle:=$
$\displaystyle\sum_{\mu}\underline{w}_{\mu}\otimes\underline{\bar{w}}_{\mu}=\mathbf{W}^{\dagger},\quad\mathbf{E}:=\sum_{\mu}\underline{l}_{\mu}\otimes\underline{\bar{w}}_{\mu},\quad\mathbf{F}:=\sum_{\mu}\underline{k}_{\mu}\otimes\underline{\bar{w}}_{\mu}.$
(9)

To get rid of terms linear in $\underline{\hat{a}}$ and
$\underline{\hat{a}}^{\dagger}$ we use a map
$\underline{\hat{a}}\to\underline{\hat{a}}+\underline{s},\quad\underline{\hat{a}}^{\prime}\to\underline{\hat{a}}^{\prime}+\underline{s}^{\prime}$
where $\underline{s},\underline{s}^{\prime}$ are constant $n$-vectors to be
determined. Let us define
$\hat{\underline{b}}=\left(\underline{\hat{a}},\underline{\hat{a}}^{\prime}\right)^{T}=\left(\underline{\hat{a}}_{0},\underline{\hat{a}}_{1},\underline{\hat{a}}^{\prime}_{0},\underline{\hat{a}}^{\prime}_{1}\right)^{T}$,
$\underline{d}=\left(\underline{s},\underline{s}^{\prime}\right)^{T}=\left(\underline{s}_{0},\underline{s}_{1},\underline{s}^{\prime}_{0},\underline{s}^{\prime}_{1}\right)^{T}$
and
$\underline{\tilde{\sigma}}=\left(\underline{\sigma}^{L},\underline{\sigma}^{R}\right)^{T}$.
The Liouvillean can be compactly written in a matrix form

$\hat{\mathcal{L}}=\left(\hat{\underline{b}}-\underline{d}\right)\cdot\mathbf{S}\left(\hat{\underline{b}}-\underline{d}\right)-\hat{S}_{0}$
(10)

where

$\mathbf{S}=\pmatrix{\mathbf{0}&-\mathbf{X}\cr-\mathbf{X}^{T}&\mathbf{Y}},$
(11)
$\mathbf{X}:=\frac{1}{2}\pmatrix{i\bar{\mathbf{H}}-\bar{\mathbf{N}}+\mathbf{M}&-2i\mathbf{K}-\mathbf{L}+\mathbf{L}^{T}\cr
2i\bar{\mathbf{K}}-\bar{\mathbf{L}}+\bar{\mathbf{L}}^{T}&-i\bar{\mathbf{H}}-\mathbf{N}+\bar{\mathbf{M}}},$
(12)
$\mathbf{Y}:=\frac{1}{2}\pmatrix{-2i\bar{\mathbf{K}}-\bar{\mathbf{L}}-\bar{\mathbf{L}}^{T}&2\mathbf{N}\cr
2\mathbf{N}^{T}&2i\mathbf{K}-\mathbf{L}-\mathbf{L}^{T}}=\mathbf{Y}^{T}.$ (13)

In order to eliminate the linear terms, we must fulfill the following
condition for the vector $\underline{d}$:

$2\mathbf{S}\underline{d}=-\mathbf{G}\underline{\tilde{\sigma}},\quad\textnormal{where}\quad\mathbf{G}:=\pmatrix{\bar{\mathbf{F}}-\mathbf{E}-i\mathbf{\Omega}^{T}&-\bar{\mathbf{F}}+\mathbf{E}+i\mathbf{\Omega}^{T}\cr-\mathbf{F}+\bar{\mathbf{E}}-i\bar{\mathbf{\Omega}}^{T}&\mathbf{F}-\bar{\mathbf{E}}\cr-\mathbf{F}-\bar{\mathbf{E}}-i\bar{\mathbf{\Omega}}^{T}&2\mathbf{F}\cr
2\bar{\mathbf{F}}&-\bar{\mathbf{F}}-\mathbf{E}+i\mathbf{\Omega}^{T}}.$ (14)

If $\mathbf{S}$ is nonsingular, then we can compute the vector
$\underline{d}=-\frac{1}{2}\mathbf{S}^{-1}\mathbf{G}\underline{\tilde{\sigma}}$.
Further, we can calculate the last term in Eq. (10):

$\displaystyle\hat{S}_{0}$ $\displaystyle=$
$\displaystyle\tr\mathbf{X}-\underline{d}\cdot\mathbf{S}\underline{d}-\underline{d}\cdot\mathbf{G}\underline{\tilde{\sigma}}-\underline{\sigma}^{L}\cdot\mathbf{W}\left(\underline{\sigma}^{R}-\underline{\sigma}^{L}\right)-\left(\underline{\sigma}^{L}-\underline{\sigma}^{R}\right)\cdot\mathbf{W}\underline{\sigma}^{R}$
(15) $\displaystyle=$
$\displaystyle\tr\mathbf{M}-\tr\mathbf{N}+\frac{1}{4}\underline{\tilde{\sigma}}\cdot\mathbf{G}^{T}\mathbf{S}^{-1}\mathbf{G}\underline{\tilde{\sigma}}-\left(\underline{\sigma}^{L}-\underline{\sigma}^{R}\right)\cdot\left(\mathbf{W}\underline{\sigma}^{R}-\bar{\mathbf{W}}\underline{\sigma}^{L}\right).$

By literally following the procedure introduced in [11], we successfully
transform the Liouvillean into its final form:

$\hat{\mathcal{L}}=-2\sum_{r=1}^{2n}\beta_{r}\hat{\zeta}_{r}^{\prime}\hat{\zeta}_{r}-\hat{S}_{0}$
(16)

where

$\underline{\hat{\zeta}}=\mathbf{P}\left(\left(\underline{\hat{a}}-\underline{s}\right)-\mathbf{Z}\left(\underline{\hat{a}}^{\prime}-\underline{s}^{\prime}\right)\right),\qquad\underline{\hat{\zeta}}^{\prime}=\mathbf{P}^{-1}\left(\underline{\hat{a}}^{\prime}-\underline{s}^{\prime}\right).$
(17)

The rapidities $\beta_{r}$ and the matrix $\mathbf{P}$ are obtained from the
diagonalization of the matrix $\mathbf{X}$, namely:

$\mathbf{X}=\mathbf{P}\mathbf{\Delta}\mathbf{P}^{-1},\qquad\mathbf{\Delta}=\textrm{diag}\left\\{\beta_{1},\dots,\beta_{2n}\right\\},$
(18)

and the matrix $\mathbf{Z}$ is a solution of the continuous Lyapunov equation

$\mathbf{X}^{T}\mathbf{Z}+\mathbf{Z}\mathbf{X}=\mathbf{Y}.$ (19)

If all the rapidities are non-negative, $\forall\beta_{r}\geq 0$, one can
obtain non-equilibrium steady states $|\textrm{NESS}^{\underline{s}}\rangle$
and all the corresponding decay modes
$|\Theta_{\underline{m}}^{\underline{s}}\rangle$
($\underline{m}\in\mathbb{Z}_{+}^{2n}$ is a multi-index) of the Liouvillean
for each combination of joint eigenvalues
$(\underline{s}^{L},\underline{s}^{R})\equiv\underline{s}$ of the hermitian
operators $(\underline{\sigma}^{L},\underline{\sigma}^{R})$, such that

$\hat{\mathcal{L}}|\textrm{NESS}^{\underline{s}}\rangle=0,$ (20)

or more precisely

$\hat{\zeta}_{r}|\textrm{NESS}^{\underline{s}}\rangle=0,$ (21)

and

$\hat{\mathcal{L}}|\Theta_{\underline{m}}^{\underline{s}}\rangle=\lambda_{\underline{m}}|\Theta_{\underline{m}}^{\underline{s}}\rangle,$
(22)

where

$|\Theta_{\underline{m}}^{\underline{s}}\rangle=\prod_{r}\frac{\left(\zeta_{r}^{\prime}\right)^{m_{r}}}{\sqrt{m_{r}!}}|\textrm{NESS}^{\underline{s}}\rangle,\qquad\lambda_{\underline{m}}=-2\sum_{r}m_{r}\beta_{r}.$
(23)

##  3 Examples

Let us consider two simple examples using the above formalism. First, let us
take $H=0$ and $L_{1}=a+z_{1}\sigma^{z}$, where $z_{1}\in\mathbb{C}$ and
$\sigma^{z}$ is the Pauli matrix for the spin degree of freedom. Thus we have
$\mathbf{M}=1$, $\mathbf{W}=|z_{1}|^{2}$, and $\mathbf{E}=\bar{z}_{1}$
($\mathbf{H}=\mathbf{K}=\mathbf{\Omega}=\mathbf{N}=\mathbf{L}=\mathbf{F}=0$).
We compute the vector $\underline{d}$ which equals to

$\underline{d}=\pmatrix{-z_{1}\sigma^{L}\cr-\bar{z}_{1}\sigma^{R}\cr-\bar{z}_{1}\sigma^{L}+\bar{z}_{1}\sigma^{R}\cr
z_{1}\sigma^{L}-z_{1}\sigma^{R}},$ (24)

and from here we get $\hat{S}_{0}=\tr\mathbf{M}=1$. Combining the
$\underline{b}$ and $\underline{d}$ vectors, we find that the result can be
written in a simple form as
$\hat{\underline{b}}-\underline{d}=\left(\hat{A}^{L},\,\hat{A}^{\dagger}{}^{R},\,\hat{A}^{\dagger}{}^{L}-\hat{A}^{\dagger}{}^{R},\,\hat{A}^{R}-\hat{A}^{L}\right)^{T}:=\left(\hat{A}_{0},\,\hat{A}_{1},\,\hat{A}_{0}^{\prime},\,\hat{A}_{1}^{\prime}\right)^{T}$,
where we have defined $\hat{A}^{L/R}=\hat{a}^{L/R}+z_{1}\sigma^{L/R}$, which
are, in a sense, bosonic operators shifted by a $z_{1}\sigma^{z}$ operator.

In this case $\mathbf{X}=\mathbf{\Delta}$ is already diagonal
($\mathbf{P}=I$), and the rapidities are $\beta_{1}=\beta_{2}=1/2$. Moreover,
since $\mathbf{Y}=0$, we also find that $\mathbf{Z}=0$. Consequently, we can
conclude that
$(\underline{\zeta},\underline{\zeta}^{\prime})^{T}=\underline{\hat{b}}-\underline{d}$.

For each sector $\underline{s}=(s^{L},s^{R})\in\left\\{-1,1\right\\}^{2}$, we
can compute the NESS:

$|\textrm{NESS}^{\underline{s}}\rangle=|\alpha_{1}\rangle\langle\alpha_{2}|\otimes|s^{L}\rangle\langle
s^{R}|,$ (25)

where $|\alpha_{1,2}\rangle$ with $\alpha_{1}=-z_{1}s^{L}$ and
$\alpha_{2}=-z_{1}s^{R}$ are coherent states, and $|s^{L}\rangle\langle
s^{R}|$ represents the spin degree of freedom. One can also compute all the
decay modes $|\Theta_{\underline{m}}^{\underline{s}}\rangle$ for each
$\underline{s}$ and $\underline{m}$.

In the second example, we consider the case where $H=0$,
$L_{1}=a+z_{1}\sigma^{z}$ and $L_{2}=z_{2}\sigma^{z}$, with
$z_{1},z_{2}\in\mathbb{C}$. This modification leads to a change in
$\mathbf{W}$ compared to the previous example, and we have
$\mathbf{W}=|z_{1}|^{2}+|z_{2}|^{2}$. Additionally, an extra term appears in
$\hat{S}_{0}$, which is now given by:

$\hat{S}_{0}=\tr\mathbf{M}+|z_{2}|^{2}\left(\sigma^{L}-\sigma^{R}\right)^{2},$
(26)

however, the vector $\underline{d}$ remains unchanged. The extra term
introduced by $L_{2}$ (i.e., the term involving $z_{2}$) in $\hat{S}_{0}$
contributes only to a spin-dependent dephasing effect.

###  3.1 Initial condition and expectation values of spin

Here, we present a calculation of the time evolution of the expectation value
of spin,
$\langle\sigma^{\mu}(t)\rangle=\textrm{Tr}\left(\sigma^{\mu}\rho(t)\right)$
($\mu=x,y,z$), given an initial condition. For this calculation, we will focus
on the first example, where $L_{1}=a+z_{1}\sigma^{z}$.

We start by writing a general solution for the density matrix:

$|\rho(t)\rangle=\sum_{\underline{s},\underline{m}}c_{\underline{m}}^{\underline{s}}e^{-\lambda_{\underline{m}}t}|\Theta_{\underline{m}}^{\underline{s}}\rangle$
(27)

where $c_{\underline{m}}^{\underline{s}}$ are coefficients related to the
initial condition
$|\rho_{0}\rangle=\sum_{\underline{s}}|\tilde{\rho}_{0}^{\underline{s}}\rangle\otimes|s^{L}\rangle\langle
s^{R}|$, with $|\tilde{\rho}_{0}^{\underline{s}}\rangle$ being the bosonic
part. To compute the coefficients $c_{\underline{m}}^{\underline{s}}$, we need
to solve the following system of linear equations for
$c_{\underline{m}}^{\underline{s}}$ obtained by combining Eqs. (23), (25), and
(27):

$|\tilde{\rho}_{0}^{\underline{s}}\rangle=\sum_{\underline{m}}c_{\underline{m}}^{\underline{s}}\,\frac{\left(\zeta_{1}^{\prime}\right)^{m_{1}}}{\sqrt{m_{1}!}}\frac{\left(\zeta_{2}^{\prime}\right)^{m_{2}}}{\sqrt{m_{2}!}}\,|-z_{1}s^{L}\rangle\langle-
z_{1}s^{R}|,\qquad\forall\underline{s}.$ (28)

We consider a scenario where the system is prepared in the initial state
$|\rho_{0}\rangle=|0\rangle\langle 0|\otimes\frac{1}{2}(I+\sigma_{x})$. In
Fig. 1, we illustrate the numerical convergence of the coefficients
$c_{\underline{m}}^{\underline{s}}$ for different truncations of the linear
system Eq. (28), and two choices of the parameter $z_{1}$.

(a)

(b)

Figure 1: Convergence of the coefficients $c_{\underline{m}}^{\underline{s}}$
for different truncations (trunc) of the linear system Eq. (28) with respect
to the truncation $\textrm{trunc}_{0}=30$ for (a) $z_{1}=0.2$, and (b)
$z_{1}=1$.

Once the coefficients $c_{\underline{m}}^{\underline{s}}$ are calculated, one
can compute the expectation values $\langle\sigma^{\mu}(t)\rangle$. For the
example discussed above, the expectation values at $t=0$ are
$\langle\sigma^{x}(0)\rangle=1$,
$\langle\sigma^{y}(0)\rangle=\langle\sigma^{z}(0)\rangle=0$. It is evident
that only the non-zero components will exhibit time-dependence, which, in this
example, is $\langle\sigma^{x}(t)\rangle$. To test our results, we also
computed $\langle\sigma^{x}(t)\rangle$ using the Python library QuTiP [28,
29], which allows for the computation of the time evolution of the density
matrix following the Lindblad master equation. Numerical calculations for the
time dependence of $\langle\sigma^{x}(t)\rangle$ (Fig. 2a) are within
numerical precision fully consistent with the conjectured analytical
expression:

$\langle\sigma^{x}(t)\rangle=\langle\sigma^{x}(0)\rangle\exp{\left[-4|z_{1}|^{2}\left(1-e^{-t}\right)\right]}.$
(29)

In Fig. 2b, we illustrate the discrepancies between $\langle\sigma^{x}\rangle$
computed using the coefficients $c_{\underline{m}}^{\underline{s}}$ obtained
from solving the linear system (28) for different truncations (trunc) and the
analytical expression (29). As expected, the discrepancy decreases as we
increase trunc. However, beyond $\textrm{trunc}\approx 18$, the solution does
not improve further due to limitations in numerical precision when inverting
the linear system (28) to obtain the coefficients.

(a)

(b)

Figure 2: Numerical and analytical solutions of the expectation value
$\langle\sigma^{x}(t)\rangle$: (a) the numerical solution obtained by QuTiP
(crosses), the analytical expression (29) (solid line), and the approximation
(33) (dashed line); (b) discrepancies between
$\langle\sigma^{x}\rangle_{coef.}$ – i.e. expectation value computed using the
coefficients $c_{\underline{m}}^{\underline{s}}$ for different truncations
(trunc) – and the analytical expression $\langle\sigma^{x}(t)\rangle_{a.exp.}$
(29) at $t=0.1$, $1$ and $10$.

We can analytically justify the expression in the limit $|z_{1}|\to 0$ (29).
To obtain the coefficients $c_{\underline{m}}^{\underline{s}}$, we perform a
Taylor expansion of the expression (28) up to quadratic terms in $z_{1}$. This
expansion results in a $4\times 4$ system of equations to solve for each
$\underline{s}=(s^{L},s^{R})$:

$\kern-63.00012pt\pmatrix{1/2\cr 0\cr 0\cr
0}\kern-5.0pt=\kern-3.00003pte^{-|z_{1}|^{2}}\kern-3.99994pt\pmatrix{1&z_{1}s^{R}&z_{1}^{*}s^{L}&-1+|z_{1}|^{2}s^{L}s^{R}\cr-
z_{1}^{*}s^{R}&1-|z_{1}|^{2}&-(z_{1}^{*})^{2}s^{L}s^{R}&z_{1}^{*}(s^{L}+s^{R})\cr-
z_{1}s^{L}&-z_{1}^{2}s^{L}s^{R}&1-|z_{1}|^{2}&z_{1}(s^{L}+s^{R})\cr|z_{1}|^{2}s^{L}s^{R}&-z_{1}s^{L}&-z_{1}^{*}s^{R}&1-|z_{1}|^{2}(2+s^{L}s^{R})}\kern-6.99997pt\pmatrix{c_{(0,0)}^{\underline{s}}\cr
c_{(0,1)}^{\underline{s}}\cr c_{(1,0)}^{\underline{s}}\cr
c_{(1,1)}^{\underline{s}}}\kern-3.00003pt,$ (30)

where the solution for the coefficients is given by:

$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\pmatrix{c_{(0,0)}^{\underline{s}}\cr
c_{(0,1)}^{\underline{s}}\cr c_{(1,0)}^{\underline{s}}\cr
c_{(1,1)}^{\underline{s}}}$ $\displaystyle=$
$\displaystyle\frac{e^{|z_{1}|^{2}}}{2(1+3|z_{1}|^{4}+2|z_{1}|^{6})}\pmatrix{1-|z_{1}|^{2}(2-s^{L}s^{R})+|z_{1}|^{4}(4+2s^{L}s^{R})\cr
z_{1}^{*}(1-|z_{1}|^{2})s^{R}\cr
z_{1}\,\,\,(1-|z_{1}|^{2})s^{L}\cr|z_{1}|^{2}(1+2|z_{1}|^{2})s^{L}s^{R}}$ (31)
$\displaystyle\approx$
$\displaystyle\frac{1}{2}\pmatrix{1-|z_{1}|^{2}(1-s^{L}s^{R})\cr
z_{1}^{*}s^{R}\cr z_{1}\,\,\,s^{L}\cr|z_{1}|^{2}s^{L}s^{R}}.$

We also need the expectation values of $\sigma^{x}$ for the slowest decay
modes to obtain the desired result. Following the notation in [11], i.e.,
$(A|\rho\rangle=\textrm{Tr}(A\rho)$, we have:

$\displaystyle(\sigma^{x}|\Theta_{(0,0)}^{\underline{s}}\rangle$
$\displaystyle=$ $\displaystyle\delta_{(-s^{L},s^{R})}e^{-2|z_{1}|^{2}},$
$\displaystyle(\sigma^{x}|\Theta_{(0,1)}^{\underline{s}}\rangle$
$\displaystyle=$
$\displaystyle\delta_{(-s^{L},s^{R})}e^{-2|z_{1}|^{2}}2z_{1}s^{R},$
$\displaystyle(\sigma^{x}|\Theta_{(1,0)}^{\underline{s}}\rangle$
$\displaystyle=$
$\displaystyle\delta_{(-s^{L},s^{R})}e^{-2|z_{1}|^{2}}2z_{1}^{*}s^{L},\qquad\forall\underline{s},$
(32)

where $\delta$ is the Kronecker delta. The final result in the limit
$|z_{1}|\to 0$ for the slowest decay modes
($\lambda_{(0,0)}=0,\lambda_{(0,1)}=\lambda_{(1,0)}=-1$) is then:

$\langle\sigma^{x}(t)\rangle=\sum_{\underline{s},\underline{m}}c_{\underline{m}}^{\underline{s}}e^{-\lambda_{\underline{m}}t}(\sigma^{x}|\Theta_{\underline{m}}^{\underline{s}}\rangle\approx
e^{-2|z_{1}|^{2}}\left(1-2|z_{1}|^{2}+4|z_{1}|^{2}e^{-t}\right).$ (33)

If we take the logarithm and approximate it for $|z_{1}|\to 0$, we arrive at
the analytical expression (29). Additionally, from Fig. 2a, we observe that
Eq. (33) approximates the analytical expression (29) and the QuTiP solution
well even for $z_{1}=0.2$ and $t\sim 1$.

It is essential to emphasize that the expressions (29) and (33) are valid only
for the specific initial condition
$|\tilde{\rho}_{0}^{\underline{s}}\rangle\sim|0\rangle\langle 0|$. For
different initial conditions, such as
$|\tilde{\rho}_{0}^{\underline{s}}\rangle\sim|n\rangle\langle n|$ with $n\neq
0$, one must include higher orders of $\underline{m}$ in the computation of
the coefficients $c_{\underline{m}}^{\underline{s}}$ to obtain the correct
result.

For completeness, let us provide the solution for the second example with
$L_{1}=a+z_{1}\sigma^{z}$ and $L_{2}=z_{2}\sigma^{z}$. As mentioned before,
the additional term $|z_{2}|^{2}\left(\sigma^{L}-\sigma^{R}\right)^{2}$ from
Eq. (26) just adds trivially to the spin-dependent dephasing. Thus, combined
with the conjectured expression (29), we obtain:

$\langle\sigma^{x}(t)\rangle=\langle\sigma^{x}(0)\rangle\exp{\left[-4|z_{1}|^{2}\left(1-e^{-t}\right)-4|z_{2}|^{2}\,t\right]}.$
(34)

In Fig. 3, we present the time evolution (34) for various values of $z_{2}$,
expressed in terms of $z_{1}$. Evidently, for large $t$, the solutions exhibit
exponential decay with a decay rate proportional to $|z_{2}|^{2}$.

Figure 3: Time evolution of the expectation value
$\langle\sigma^{x}(t)\rangle$ (see Eq. (34)) for different parameters $z_{2}$
and initial value $\langle\sigma^{x}(0)\rangle=1$.

###  3.2 Limiting cases

Lastly, let us consider two limiting cases. In the first case suppose again
$|z_{1}|\to 0$. Then we have $\langle\sigma^{x}(t)\rangle\approx
1-4|z_{1}|^{2}(1-e^{-t})$ and in the limit $t\to\infty$ this reduces to
$1-4|z_{1}|^{2}$.

In the second case, suppose slightly modified
$L_{1}=wa+\sigma^{z}=w(a+\sigma^{z}/w)$ and $|w|\to 0$ such that $|w|^{2}t\ll
1$ for some value of $t$. Now we get

$\langle\sigma^{x}(t)\rangle=\exp\left[-\frac{4}{|w|^{2}}(1-e^{-|w|^{2}t})\right]\mathchoice{\vtop{\halign{#\cr$\hfil\phantom{}{}^{\phantom{}{}^{|w|^{2}t\ll
1}}\hfil$\crcr\nointerlineskip\kern-3.60004pt\cr$\displaystyle\relbar\mkern-7.0mu\cleaders\hbox{$\displaystyle\mkern-2.0mu\relbar\mkern-2.0mu$}\hfill\mkern-7.0mu\rightarrow$\crcr}}}{\vtop{\halign{#\cr$\hfil\phantom{}{}^{\phantom{}{}^{|w|^{2}t\ll
1}}\hfil$\crcr\nointerlineskip\kern-3.60004pt\cr$\textstyle\relbar\mkern-7.0mu\cleaders\hbox{$\textstyle\mkern-2.0mu\relbar\mkern-2.0mu$}\hfill\mkern-7.0mu\rightarrow$\crcr}}}{\vtop{\halign{#\cr$\hfil\phantom{}{}^{\phantom{}{}^{|w|^{2}t\ll
1}}\hfil$\crcr\nointerlineskip\kern-3.60004pt\cr$\scriptstyle\relbar\mkern-7.0mu\cleaders\hbox{$\scriptstyle\mkern-2.0mu\relbar\mkern-2.0mu$}\hfill\mkern-7.0mu\rightarrow$\crcr}}}{\vtop{\halign{#\cr$\hfil\phantom{}{}^{\phantom{}{}^{|w|^{2}t\ll
1}}\hfil$\crcr\nointerlineskip\kern-3.60004pt\cr$\scriptscriptstyle\relbar\mkern-7.0mu\cleaders\hbox{$\scriptscriptstyle\mkern-2.0mu\relbar\mkern-2.0mu$}\hfill\mkern-7.0mu\rightarrow$\crcr}}}\,e^{-4t}$
(35)

and we reconstructed the exponential decoherence of the Lindbladian
$L=\sigma^{z}$ with the rate $\gamma=2$ (note the factor 2 in the definition
of the Liouvillean (1)).

##  4 Conclusion

In conclusion, we have studied the dynamics of a quantum system – which is
quadratic in bosonic degrees of freedom and linear in additional commuting
(classical) degrees of freedom – in the framework of the formalism of the
third quantization. In particular, we succeeded in diagonalizing the Liouville
superoperator and obtaining nonequilibrium states and the corresponding decay
modes for various combinations of eigenvalues of the mutually commuting
Hermitian operators. A special case is the harmonic oscillator with spin
degrees of freedom. Analyzing two specific examples solved exactly, we
observed intriguing behavior with exponential decay rates of the spin
expectation values.

For the first example, we analyzed the time evolution of the expectation value
of spin $\langle\sigma^{x}(t)\rangle$ for $L_{1}=a+z_{1}\sigma^{z}$. We have
verified numerically the dependence on the parameter $z_{1}$, but the analytic
proof of full generality remains an open question. In the second example, we
considered a modified system with the additional term $L_{2}=z_{2}\sigma^{z}$
in the Lindbladian. The analysis of this case showed an exponential decay with
a decay rate proportional to $|z_{2}|^{2}$ for large $t$. The results
emphasize the influence of the parameter $z_{2}$ on the long-term relaxation
behavior of the system. Moreover, we have investigated two limiting cases and
recovered the known time dependence for the Lindbladians $L=a$ and
$L=\sigma_{z}$.

The results presented provide valuable contributions to understanding the
dynamics of open quantum systems and provide a basis for further studies of
related physical phenomena in which spin degrees of freedom are coupled with
bosonic ones. For example, these findings could be particularly important for
understanding decoherence in spin systems coupled to a bosonic thermal bath
and confined in quadratic potential traps [27]. Such studies could provide
important insights into the behavior of quantum systems under realistic
experimental conditions and pave the way for advances in quantum technology
and quantum information processing.

## Acknowledgments

LM and AR acknowledge support from the Slovenian Research and Innovation
Agency (ARIS) under Contract No. P1-0044; AR was also supported by Grant
J2-2514. TP was supported by Programme P1-0402 and Grants N1-0219 and N1-0233
(ARIS).

## References

## References

  * [1] Schütz G M 2001 Exactly solvable models for many-body systems far from equilibrium
  * [2] Prosen T 2015 Journal of Physics A: Mathematical and Theoretical 48 373001
  * [3] Landi G T, Poletti D and Schaller G 2022 Rev. Mod. Phys. 94(4) 045006
  * [4] Breuer H P and Petruccione F 2007 The Theory of Open Quantum Systems (Oxford University Press) ISBN 9780199213900
  * [5] Lindblad G 1976 Communications in Mathematical Physics 48 119–130
  * [6] Gorini V, Kossakowski A and Sudarshan E C G 2008 Journal of Mathematical Physics 17 821–825 ISSN 0022-2488
  * [7] Kraus B, Büchler H P, Diehl S, Kantian A, Micheli A and Zoller P 2008 Phys. Rev. A 78(4) 042307
  * [8] Diehl S, Micheli A, Kantian A, Kraus B, Büchler H P and Zoller P 2008 Nature Physics 4 878–883
  * [9] Prosen T 2008 New Journal of Physics 10 043026
  * [10] Prosen T 2010 Journal of Statistical Mechanics: Theory and Experiment 2010 P07020
  * [11] Prosen T and Seligman T H 2010 Journal of Physics A: Mathematical and Theoretical 43 392004
  * [12] Prosen T, Martignon L and Seligman T H 2015 Physica Scripta 90 074036
  * [13] Dzhioev A A and Kosov D S 2011 The Journal of Chemical Physics 134 044121 ISSN 0021-9606
  * [14] Guo C and Poletti D 2017 Phys. Rev. A 95(5) 052107
  * [15] B T and Z Y 2021 Journal of Statistical Mechanics: Theory and Experiment 2022
  * [16] Zhang Y and Barthel T 2022 Phys. Rev. Lett. 129(12) 120401
  * [17] Eisler V 2011 Journal of Statistical Mechanics: Theory and Experiment 2011 P06007
  * [18] Temme K, Wolf M M and Verstraete F 2012 New Journal of Physics 14 075004
  * [19] Eisert J and Prosen T 2010 arXiv: Quantum Physics
  * [20] Höning M, Moos M and Fleischhauer M 2012 Phys. Rev. A 86(1) 013606
  * [21] Prosen T and Pižorn I 2008 Phys. Rev. Lett. 101(10) 105701
  * [22] Prosen T and Žunkovič B 2010 New Journal of Physics 12 025016
  * [23] Braak D 2011 Phys. Rev. Lett. 107(10) 100401
  * [24] Skrypnyk T 2009 Journal of Mathematical Physics 50 103523 ISSN 0022-2488
  * [25] Čadež T, Jefferson J H and Ramšak A 2013 New Journal of Physics 15 013029
  * [26] Čadež T, Jefferson J H and Ramšak A 2014 Phys. Rev. Lett. 112(15) 150402
  * [27] Donvil B, Ulčakar L, Rejec T and Ramšak A 2020 Phys. Rev. B 101(20) 205427
  * [28] Johansson J R, Nation P D and Nori F 2012 Computer Physics Communications 183 1760–1772 ISSN 0010-4655
  * [29] Johansson J R, Nation P D and Nori F 2013 Computer Physics Communications 184 1234–1240 ISSN 0010-4655