Bosonic Casimir Effect in an Aether-like Lorentz-Violating Scenario with Higher Order Derivatives

Abstract

In this paper, we investigate the bosonic Casimir effect in a Lorentz-violating symmetry scenario. The theoretical model adopted consists of a real massive scalar quantum field confined in a region between two large parallel plates, having its dynamics governed by a modified Klein–Gordon equation that presents a Lorentz symmetry-breaking term. In this context, we admit that the quantum field obeys specific boundary conditions on the plates. The Lorentz-violating symmetry is implemented by the presence of an arbitrary constant space-like vector in a CPT-even aether-like approach, considering a direct coupling between this vector with the derivative of the field in higher order. The modification of the Klein–Gordon equation produces important corrections on the Casimir energy and pressure. Thus, we show that these corrections strongly depend on the order of the higher derivative term and the specific direction of the constant vector, as well as the boundary conditions considered.

1. Introduction

The Casimir effect is one of the most important macroscopic consequences of the existence of the quantum vacuum. Although being theoretically proposed in 1948 by H. B. Casimir [1], the theoretical result was only verified to be compatible with experiments 10 years later by M. J. Sparnaay [2]. In the 1990s, experiments confirmed the Casimir effect with a high degree of accuracy [3,4]. In their original work, Casimir predicted that due to the quantum vacuum fluctuations associated with the electromagnetic field, two parallel flat neutral (grounded) plates, separated by a distance a, attract each other with a force per unit area given by:

$$\frac{F}{A}=-\frac{{\pi }^2\hslash c}{240a^4},$$

(1)

where A is the area of the plates.

In general, the Casimir effect is defined as being a force per unit area, when boundary conditions are imposed on quantum fields. The simplest theoretical device to study the Casimir effect is constituted by two neutral parallel plates placed in the (classical) vacuum. As the quantum vacuum consists of an infinite set of waves that contemplates all possible wavelengths, when the plates are considered, only a few wavelengths are allowed between them.

The Lorentz invariance, which is a cornerstone of Quantum Field Theory, was questioned in a work by V. A. Kostelecky and S. Samuel [5], that describes a mechanism in string theory that allows the violation of Lorentz symmetry at the Planck energy scale. According to this mechanism, the violation of the Lorentz symmetry is introduced by the emergence of non-vanishing vacuum expectation values of some vector and tensor components, which imply preferential directions, providing in this way a space–time anisotropy. In the quantum gravity context, Hořava–Lifshitz (HL) proposed a theory [6] with the objective to implement the possible convergence of quantum corrections by imposing different properties of scales in which the space and time coordinates are set. The HL approach clearly provides an anisotropy between space and time and can be applied not only to gravity, but also to other field theory models, including scalar, spinor, and gauge theories. Among the most important results achieved in these scenarios, one can emphasize the calculation of the one-loop effective potential in the HL-like QED and HL-like Yukawa model [7,8,9] and the study of different issues related to the renormalization of these theories [10,11,12,13]. If there is a violation of the Lorentz symmetry at the Planck energy scale in a more fundamental theory, the effects of this breakdown must manifest themselves in other energy scales in different QFT models. Other mechanisms of violation of Lorentz symmetry are possible, such as space–time noncommutativity [14,15,16,17,18].

The violation of the Lorentz symmetry became of great experimental interest. For instance, in the search for vestiges left by the violation, the high-accuracy experimental measurements of the Casimir pressure became a great ally in the study of the Lorentz symmetry breakdown in theoretical models within field theory.

The first analyses of Casimir energy in the Lorentz-violating (LV) theories were developed in [19,20,21], considering different Lorentz-breaking extensions of the QED. In addition, the studies of Casimir effects associated with massless scalar and fermionic quantum fields confined in the region between two large parallel plates, taking into account HL formalism, have been investigated in Refs. [22,23], respectively. More recently, the Casimir effect associated with a massive real scalar field was developed in Ref. [24].

Considering direct coupling between the derivative of the field with an arbitrary constant four-vector in an aether-like CPT-even Lorentz symmetry breaking, the analysis of Casimir effects associated with real scalar and fermionic massive fields has been investigated in Refs. [25,26], respectively. Moreover, local Casimir densities in an LV scenario have been analyzed in [27,28]. In Ref. [29], the influence of a constant magnetic field on the Casimir effect in the Lorentz-violating scalar field was considered. The thermal effect on the Casimir energy and pressure caused by the Lorentz-violating scalar field was investigated in [30]. The analyses of the Casimir energy and topological mass associated with a massive scalar field in the LV scenario were considered in Refs. [31,32].

In this paper, we intend to continue in the same line of investigation, i.e., analyzing the Casimir effect associated with a massive scalar quantum field in an LV scenario; however, at this time we shall take into consideration that the LV is implemented by a new term that involves higher order derivatives of the field, coupled to a space-like constant vector. In this way, we may understand this term as a combination of the HL methodology with the aether-like CPT-even Lorentz symmetry breaking.

This paper is organized as follows: In Section 2, we briefly introduce the theoretical model that governs the dynamics of the real scalar field. We present the LV bosonic action and the corresponding modified Klein–Gordon equation. In Section 3, we develop the calculation of the Casimir energies in cases where the constant vector is parallel and orthogonal to the plates. In order to confine the bosonic field between the two parallel plates, we should impose that the flux of virtual particles crossing the plates is zero. This can be conducted by imposing Dirichlet, Neumann, or mixed boundary conditions on the field at the plates. Finally, in Section 4 are conclusions of our most relevant remarks found in this paper. Here, units are assumed to be $\hslash =c=1$, and the metric signature will be taken as $(+,-,-,-)$.

2. Klein–Gordon Equation in Aether-like Lorentz Symmetry Violation Scenario with Higher-Order Derivatives

In this section, we introduce the theoretical model that we want to investigate. It is composed of a massive scalar quantum field in a Lorentz-violating symmetry scenario introduced by the presence of a constant space-like vector in an aether-like approach, considering its direct coupling with the derivative of the field in higher order. In this sense, the Lorentz violation symmetry is caused by the presence of a constant background vector and by an anisotropy between space and time coordinates due to scaling transformations. This model is formally given by the Lagrangian density below:

$$\mathcal{L}=\frac{1}{2}\left[\left({\partial }_{\mu }\varphi \right)\left({\partial }^{\mu }\varphi \right)-l^{2(ϵ-1)}{(-1)}^{ϵ}{\left[{\left(u^{\mu }{\partial }_{\mu }\right)}^{ϵ}\varphi \right]}^2-m^2{\varphi }^2\right].$$

(2)

In the above Lagrangian, the parameter l is of the order of the inverse of the energy scale where the Lorentz symmetry is broken. The dimensionless constant vector $u^{\mu }$ that is associated with a preferential direction couples to the scalar field through its derivative, as explained above, and the parameter $ϵ$ is an integer number.

In this formalism, the modified Klein–Gordon Equation (KG) reads,

$$\left[\square +l^{2(ϵ-1)}{\left(u^{\mu }{\partial }_{\mu }\right)}^{2ϵ}+m^2\right]\varphi =0.$$

(3)

For $ϵ=1$, the above equation coincides with the one presented in [25]. In the latter, the analyses of the Casimir energy and pressure were considered, admitting that the constant four-vector is both time-like and space-like, separately. Because we are interested in investigating the behavior of the scalar field in higher order derivatives theory, we will consider $ϵ\ge 2$. Moreover, to avoid unitarity problems, we will also assume that the vector $u^{\mu }$ is only space-like. Another way to implement a higher-order derivative in an LV scenario is by considering a higher-order time derivative of the field; however, the presence of this term may violate unitarity of the theory.

Imposing that the action associated with the Lagrangian (2) is invariant under the infinitesimal translation, $x^{\mu }\to x^{\mu }+\delta a^{\mu }$, the obtained energy–momentum tensor (EM) reads (our derivation of EM was developed by an induction procedure, i.e, we first assumed $ϵ=2$, followed by $ϵ=3$, and $ϵ=4$),

$$\begin{array}{ccc}\hfill T^{\mu \nu }& =& \frac{\partial \mathcal{L}}{\partial \left({\partial }_{\mu }\varphi \right)}{\partial }^{\nu }\varphi +\frac{\partial \mathcal{L}}{\partial ({\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }_{\mu }\varphi )}{\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }^{\nu }\varphi -{\partial }_{{\mu }_1}\frac{\partial \mathcal{L}}{\partial ({\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }_{\mu }\varphi )}{\partial }_{{\mu }_2}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }^{\nu }\varphi \hfill \\ & +& {\partial }_{{\mu }_1}{\partial }_{{\mu }_2}\frac{\partial \mathcal{L}}{\partial ({\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }_{\mu }\varphi )}{\partial }_{{\mu }_3}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }^{\nu }\varphi +\dots +{(-1)}^{ϵ-1}{\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}\frac{\partial \mathcal{L}}{\partial ({\partial }_{{\mu }_1}\dots {\partial }_{{\mu }_{ϵ-1}}{\partial }_{\mu }\varphi )}{\partial }^{\nu }\varphi -{\eta }^{\mu \nu }\mathcal{L}.\hfill \end{array}$$

(4)

Substituting the Lagrangian (2) into (4), we obtain:

$$\begin{array}{ccc}\hfill T^{\mu \nu }& =& \left({\partial }^{\mu }\varphi \right)\left({\partial }^{\nu }\varphi \right)+ϵ!l^{2(ϵ-1)}{u}^{\mu }\{\left[{(u·\partial )}^{2ϵ-1}\varphi \right]\left({\partial }^{\nu }\varphi \right)-\left[{(u·\partial )}^{2ϵ-2}\varphi \right](u·\partial )\left({\partial }^{\nu }\varphi \right)\hfill \\ & +& \left[{(u·\partial )}^{2ϵ-3}\varphi \right]{(u·\partial )}^2\left({\partial }^{\nu }\varphi \right)+\dots -{(-1)}^{ϵ}\left[{(u·\partial )}^{ϵ}\varphi \right]{(u·\partial )}^{ϵ-1}\left({\partial }^{\nu }\varphi \right)\}-{\eta }^{\mu \nu }\mathcal{L}.\hfill \end{array}$$

(5)

Although by construction, the general expression for the energy–momentum tensor, Equation (4), satisfies the condition,

$$\begin{array}{c}\hfill {\partial }_{\mu }{T}^{\mu \nu }=0,\end{array}$$

(6)

we have explicitly shown, for different values of $ϵ$, that the above equation is obeyed by the Lagrangian (2).

Furthermore, we can see that the energy–momentum tensor is not symmetric; its anti-symmetric part is given by

$$\begin{array}{ccc}\hfill T^{\mu \nu }-T^{\nu \mu }& =& ϵ!l^{2(ϵ-1)}\{\left[{(u·\partial )}^{2ϵ-1}\varphi \right]-\left[{(u·\partial )}^{2ϵ-2}\varphi \right](u·\partial )\hfill \\ & +& \left[{(u·\partial )}^{2ϵ-3}\varphi \right]{(u·\partial )}^2+\dots -{(-1)}^{ϵ}\left[{(u·\partial )}^{ϵ}\varphi \right]{(u·\partial )}^{ϵ-1}\}\left[u^{\mu }\left({\partial }^{\nu }\varphi \right)-u^{\nu }\left({\partial }^{\mu }\varphi \right)\right].\hfill \end{array}$$

(7)

This anti-symmetry is typical for Lorentz symmetry violation formalism. Moreover, we have explicitly checked that ${\partial }_{\nu }{T}^{\mu \nu }\ne 0$.

3. The Casimir Effect in the Context of Higher-Order Derivatives Lorentz Symmetry Violation

The main objective of this section is to analyze how the LV symmetry represented by the presence of a higher-order derivative term of the field, along a specific direction, modifies the dispersion relations responsible for the deviation of the Casimir energy and pressure when compared with the scenario preserving the Lorentz symmetry. As we have already mentioned, we will assume that the scalar field operator, $\widehat{\varphi }\left(x\right)$, satisfies specific boundary conditions on the plates exhibited in Figure 1.

In order to obtain the Hamiltonian operator, $\widehat{H}$, we first calculate the normalized set of positive/negative energy solutions of Equation (3) that obey specific boundary conditions on the plates. Performing this, we can calculate the total vacuum energy of the system and then determine the Casimir energy for each case considered.

3.1. Dirichlet Condition

It has been shown in [25] that the quantum field operator, $\widehat{\varphi }\left(x\right)$ is expressed in terms of the normalized positive/negative energy solution of the Klein–Gordon equation that satifies Dirichlet boundary condition on the plates, $z=0$ and $z=a$; that is,

$$\varphi {\left(x\right)}_{z=0}=\varphi {\left(x\right)}_{z=a}=0,$$

(8)

has the general form given below,

$$\begin{array}{c}\hfill \widehat{\varphi }\left(x\right)=\int {\mathrm{d}}^2\mathbf{k}\sum _{n=1}^{\infty }\frac{1}{{\left[{\left(2\pi \right)}^{2}a{\omega }_{\mathbf{k},n}\right]}^{1/2}}sin\left({\displaystyle \frac{n\pi }{a}}z\right)[{\widehat{a}}_{\mathbf{k},n}{e}^{-ikx}+{\widehat{a}}_{\mathbf{k},n}^{†}{e}^{ikx}],\end{array}$$

(9)

where ${\widehat{a}}_{\mathbf{k},n}$ and ${\widehat{a}}_{\mathbf{k},n}^{†}$ correspond to the annihilation and creation operators, respectively. These operators satisfy the following commutation relations

$$\begin{array}{cc}\hfill [{\widehat{a}}_{\mathbf{k},n},{\widehat{a}}_{{\mathbf{k}}^{\prime },n^{\prime }}^{†}]& ={\delta }_{n,n^{\prime }}{\delta }^2(\mathbf{k}-{\mathbf{k}}^{\prime }),\hfill \\ \hfill [{\widehat{a}}_{\mathbf{k},n},{\widehat{a}}_{{\mathbf{k}}^{\prime },n^{\prime }}]=& [{\widehat{a}}_{\mathbf{k},n}^{†},{\widehat{a}}_{{\mathbf{k}}^{\prime },n^{\prime }}^{†}]=0,\hfill \end{array}$$

(10)

being $kx\equiv {\omega }_{\mathbf{k},n}t-k_{x}x-k_{y}y-k_{z}z$. The explicit form of ${\omega }_{\mathbf{k},n}$ will depend on the specific model adopted for the LV term. This subject will be explored in the following subsections.

3.1.1. Vector Parallel to the Plates

The space-like four-vector $u^{\mu }$ can be in three different directions, parallel to the plates, $u^{\mu }=(0,1,0,0)$ and $u^{\mu }=(0,0,1,0)$, and perpendicular to the plates, $u^{\mu }=(0,0,0,1)$. The dispersion relations associated with the first two vectors are the same. Thereby, here in this subsection, we will consider a parallel vector as being,

$$\begin{array}{c}\hfill u^{\mu }=(0,1,0,0).\end{array}$$

(11)

The corresponding dispersion relation is

$$\begin{array}{c}\hfill {\omega }_{\mathbf{k},n}^2=k_x^2+k_y^2+l^{2(ϵ-1)}{(-1)}^{ϵ}{k}_x^{2ϵ}+{\left(\frac{n\pi }{a}\right)}^2+m^2.\end{array}$$

(12)

Hence, the Hamiltonian operator $\widehat{H}$, resulting from the canonical quantization, reads

$$\begin{array}{c}\hfill \widehat{H}=\frac{1}{2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=1}^{\infty }{\omega }_{\mathbf{k},n}\left[2{\widehat{a}}_{\mathbf{k},\mathbf{n}}^{†}{\widehat{a}}_{\mathbf{k},\mathbf{n}}+\frac{L^2}{{\left(2\pi \right)}^2}\right].\end{array}$$

(13)

The vacuum energy is obtained by taking the vacuum expectation value of $\widehat{H}$:

$$\begin{array}{c}\hfill E_0=\langle 0|\widehat{H}|0\rangle =\frac{L^2}{8{\pi }^2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=1}^{\infty }{\omega }_{\mathbf{k},n}.\end{array}$$

(14)

Performing a change of coordinates from Cartesian coordinate $(k_x,k_y)$ to polar one, $(k,\theta )$, and making a change in variable $u=ak$, we obtain

$$\begin{array}{c}\hfill E_0=\frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }d\theta {\int }_0^{\infty }udu\sum _{n=1}^{\infty }{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2+{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{u}^{2ϵ}{cos}^{2ϵ}\theta \right]}^{\frac{1}{2}}.\end{array}$$

(15)

Because it is our interest to investigate the LV correction on the Casimir energy due to higher-order space derivative term, we will consider $ϵ\ge 2$. Although for this case the integral over the variable u can be evaluated, the result is not very enlightening. Furthermore, we have not found in the literature the integral over $\theta$ for general values of $ϵ$, even for $ϵ=2$. Thus, in order to provide a quantitative result for the correction of the Casimir energy caused by the Lorentz-violating term, we develop an expansion in the parameter associated with the Lorentz violation. By performing an expansion up to the first order in the parameter $\frac{l}{a}\ll 1$, expression (15) can be written as

$$\begin{array}{ccc}\hfill E_0& \approx & \frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }d\theta {\int }_0^{\infty }udu\sum _{n=1}^{\infty }\{{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\hfill \\ & +& \frac{1}{2}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{u}^{2ϵ}{cos}^{2ϵ}\theta {\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}\},\hfill \end{array}$$

(16)

where the first term is associated with the vacuum energy without Lorentz violation. Thus, after integration over the angular coordinate, the LV term becomes

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\sum _{n=1}^{\infty }{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(17)

The Casimir energy by unit area, associated with a massive scalar quantum field confined between two large and parallel plates of area $L^2$, separated by a distance a, that obeys the Dirichlet boundary condition, was obtained in [25] in an integral representation by,

$$\begin{array}{c}\hfill \frac{E_C}{L^2}=-\frac{a{m}^4}{6{\pi }^2}{\int }_1^{\infty }\frac{{(v^2-1)}^{\frac{3}{2}}dv}{e^{2amv}-1}.\end{array}$$

(18)

Because our main interest in this research is to investigate the contribution of the Lorentz symmetry breaking in the Casimir energy, we will focus our analysis on Equation (17). As our first step to evaluate this contribution, we will use the Abel–Plana summation formula below [33] to develop the summation over the quantum number n, i.e.,

$$\begin{array}{c}\hfill \sum _{n=0}^{\infty }F\left(n\right)=\frac{1}{2}F\left(0\right)+{\int }_0^{\infty }F\left(t\right)\mathrm{d}t+i{\int }_0^{\infty }\frac{\mathrm{d}t}{e^{2\pi t}-1}[F\left(it\right)-F(-it)].\end{array}$$

(19)

Then, expression (17) becomes

$$\begin{array}{ccc}\hfill \widetilde{E_0}& =& \frac{L^2}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\{-\frac{1}{2}F\left(0\right)+{\int }_0^{\infty }F\left(t\right)\mathrm{d}t\hfill \\ & +& i{\int }_0^{\infty }\frac{F\left(it\right)-F(-it)}{e^{2\pi t}-1}\mathrm{d}t\},\hfill \end{array}$$

(20)

where

$$\begin{array}{c}\hfill F\left(n\right)={\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(21)

Note that the first term on the right side of (20) refers to the vacuum energy in the presence of only one plate, whereas the second term refers to the vacuum energy without boundary. Both terms are divergent and do not contribute to the Casimir energy. As a result, the LV contribution to the Casimir energy per unit area of the plates is given by

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& \frac{i}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}t\frac{{[u^2+{\left(it\pi \right)}^2+{\left(ma\right)}^2]}^{-1/2}-{[u^2+{(-it\pi )}^2+{\left(ma\right)}^2]}^{-1/2}}{e^{2\pi t}-1}.\hfill \end{array}$$

(22)

Performing a change of variable, with $t\pi =v$, we obtain

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& \frac{i}{8{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}v\frac{{[u^2+{\left(ma\right)}^2+{\left(iv\right)}^2]}^{-1/2}-{[u^2+{\left(ma\right)}^2+{(-iv)}^2]}^{-1/2}}{e^{2v}-1}.\hfill \end{array}$$

(23)

The integral over the variable v must be considered in two cases, for ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}>v$ and ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}<v$, taking into account that we have

• For the case ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}>v$:

  $$\begin{array}{c}\hfill {[u^2+{\left(ma\right)}^2+{(\pm iv)}^2]}^{-1/2}={[u^2+{\left(ma\right)}^2-v^2]}^{-1/2};\end{array}$$

  (24)
• For the case ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}<v$:

  $$\begin{array}{c}\hfill {[u^2+{\left(ma\right)}^2+{(\pm iv)}^2]}^{-1/2}=\mp i{\left[v^2-(u^2+{\left(ma\right)}^2)\right]}^{-\frac{1}{2}}.\end{array}$$

  (25)

Consequently, the integral in u over the interval $[0,{(u^2+{\left(ma\right)}^2)}^{1/2}]$ vanishes. Thus, it remains

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& \frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_{{\left[u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}}^{\infty }\mathrm{d}v\frac{{[v^2-(u^2+{\left(ma\right)}^2)]}^{-1/2}}{e^{2v}-1}.\hfill \end{array}$$

(26)

Furthermore, performing the new change of variable ${\rho }^2=v^2-(u^2+{\left(ma\right)}^2)$, we find

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& \frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\frac{\mathrm{d}\rho }{{\left[{\rho }^2+u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\left[e^{2{({\rho }^2+u^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}-1\right]}.\hfill \end{array}$$

(27)

Now, by making a change of coordinates from the plane $(u,\rho )$ to the polar one, we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_0^{\infty }\frac{{\sigma }^{2(ϵ+1)}\mathrm{d}\sigma }{{\left[{\sigma }^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\left[e^{2{({\sigma }^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}-1\right]}.\end{array}$$

(28)

For massless field, the above integral becomes,

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_0^{\infty }\frac{{\sigma }^{2ϵ+1}\mathrm{d}\sigma }{e^{2\sigma }-1}.\end{array}$$

(29)

We can further make use of the integral given by [34]

$$\begin{array}{c}\hfill {\int }_0^{\infty }\frac{x^{\nu -1}\mathrm{d}x}{e^{\mu x}-1}=\frac{1}{{\mu }^{\nu }}\Gamma \left(\nu \right)\zeta \left(\nu \right),\end{array}$$

(30)

where $\Gamma \left(\nu \right)$ and $\zeta \left(\nu \right)$ correspond to the Gamma and Riemann zeta functions [35], respectively. Thus, we find

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}\frac{\Gamma (2ϵ+2)\zeta (2ϵ+2)}{2^{(2ϵ+2)}}.\end{array}$$

(31)

In this case, for instance, we have

• For $ϵ=2$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\pi }^4}{10,080a^3}{\left(\frac{l}{a}\right)}^2;\end{array}$$

  (32)
• For $ϵ=3$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\pi }^6}{13,440a^3}{\left(\frac{l}{a}\right)}^4.\end{array}$$

  (33)

At this point, we would like to emphasize that Equation (31) has been derived considering $ϵ\ge 2$. For this reason, our results for LV Casimir energies cannot reduce to the corresponding ones obtained in [26,27,28,31] by taking its limit $ϵ=1$.

Unfortunately for the massive case, Equation (28) can only be given in terms of an infinite sum of modified Bessel functions as shown below; let us evaluate its asymptotic limits for small and large values of the dimensionless parameter $ma$. In order to perform this, we make the following change of variables ${\xi }^2={\sigma }^2+{\left(ma\right)}^2$ and $\xi =mav$. Thus, it gives

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{ϵ+\frac{1}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

(34)

Knowing that the geometric series can be represented as

$$\begin{array}{c}\hfill \frac{1}{e^{2amv}-1}=\sum _{j=1}^{\infty }{e}^{-2amvj},\end{array}$$

(35)

we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}\sum _{j=1}^{\infty }{\int }_1^{\infty }{\left(v^2-1\right)}^{ϵ+\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(36)

The integral representation of the modified Bessel function, $K_{\mu }\left(z\right)$ [35]

$$\begin{array}{c}\hfill K_{\nu }\left(x\right)=\frac{{\left(\frac{x}{2}\right)}^{\nu }\Gamma \left(\frac{1}{2}\right)}{\Gamma \left(\nu +\frac{1}{2}\right)}{\int }_1^{\infty }\mathrm{d}t{\left(t^2-1\right)}^{\nu -\frac{1}{2}}{e}^{-xt},\end{array}$$

(37)

allows us to put expression (42) in the form

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{ϵ+1}\Gamma \left(ϵ+\frac{3}{2}\right){(-1)}^{ϵ}}{4{\left(\pi \right)}^{\frac{5}{2}}{a}^3(2ϵ+1)}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\sum _{j=1}^{\infty }\frac{K_{ϵ+1}\left(2amj\right)}{j^{ϵ+1}}.\end{array}$$

(38)

Let us now consider the expression above in two asymptotic regime cases:

(i)

The LV Casimir energy, $\widetilde{E_C}$, for large values of $am>>1$, can be obtained by using the asymptotic expression for the modified Bessel function for large arguments [35]:

$$\begin{array}{c}\hfill K_{\nu }\left(z\right)\approx \sqrt{\frac{\pi }{2z}}{e}^{-z}.\end{array}$$

(39)

The dominant contribution is for $j=1$. Thus, we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx \frac{{\left(am\right)}^{ϵ+\frac{1}{2}}\Gamma \left(ϵ+\frac{3}{2}\right){(-1)}^{ϵ}}{8{\pi }^2{a}^3(2ϵ+1)}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{e}^{-2am}.\end{array}$$

(40)

We can observe that the Casimir energy decays exponentially.

(ii)

For $am<<1$, it is better to consider Equation (34):

• In the case $ϵ=2$, expression (34) becomes

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^6}{20{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{\frac{5}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

  (41)

We can approximate the integrand as shown below and obtain a series expansion, i.e.,

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & \frac{{\left(am\right)}^6}{20{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{\left(v^5-\frac{5}{2}{v}^3+\frac{15}{8}v\right)\mathrm{d}v}{e^{2amv}-1}\hfill \\ & \approx & \frac{1}{10,080{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2\left[{\pi }^6-\frac{21}{4}{\pi }^4{\left(am\right)}^2+\frac{315}{8}{\pi }^2{\left(am\right)}^4\right].\hfill \end{array}$$

(42)

The LV Casimir pressure can be obtained through the standard procedure:

$$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=-\frac{1}{L^2}\frac{\partial \widetilde{E_C}}{\partial a}.\end{array}$$

(43)

Taking the approximated expression (42), we obtain:

$$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=\frac{1}{80,640a^4}{\left(\frac{l}{a}\right)}^2\left[40{\pi }^4-126{\pi }^2{\left(am\right)}^2+315{\left(am\right)}^4\right].\end{array}$$

(44)

• In the case $ϵ=3$, expression (34) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^8}{28{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{\frac{7}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

  (45)

Adopting the same procedure to approximate the integrand as above, the series expansion is

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & -\frac{{\left(am\right)}^8}{28{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{\left(v^7-\frac{7}{2}{v}^5+\frac{35}{8}{v}^3-\frac{35}{16}v\right)\mathrm{d}v}{e^{2amv}-1}\hfill \\ & \approx & -\frac{1}{13,440{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4\left[{\pi }^8-\frac{10}{3}{\pi }^6{\left(am\right)}^2+\frac{35}{4}{\pi }^4{\left(am\right)}^4\right].\hfill \end{array}$$

(46)

Taking (43), the corresponding LV Casimir pressure is

$$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=-\frac{1}{161,280a^4}{\left(\frac{l}{a}\right)}^4\left[84{\pi }^6-200{\pi }^4{\left(am\right)}^2+315{\pi }^2{\left(am\right)}^4\right].\end{array}$$

(47)

In Figure 2, we present the behavior of the Casimir energy per unit area multiplied by $a^3$, ${\epsilon }_c=\frac{\widetilde{E_C}}{L^2}{a}^3$, as a function of $ma$, considering as only an illustrative example $\frac{l}{a}=0.01$, for two distinct values of $ϵ$. In the plot on the left we consider $ϵ=2$, whereas on the right we consider $ϵ=3$.

3.1.2. Vector Perpendicular to the Plates

Let us now consider that the four-vector $u^{\mu }$ is perpendicular to the plates, i.e.,

$$\begin{array}{c}\hfill u^{\mu }=(0,0,0,1).\end{array}$$

(48)

The corresponding dispersion relation is

$$\begin{array}{c}\hfill {\omega }_{\mathbf{k},n}^2=k_x^2+k_y^2+{\left(\frac{n\pi }{a}\right)}^2+l^{2(ϵ-1)}{(-1)}^{ϵ}{\left(\frac{n\pi }{a}\right)}^{2ϵ}+m^2.\end{array}$$

(49)

For this case, the Hamiltonian operator, $\widehat{H}$, has the same structure as (13), and consequently the vacuum energy is given by

$$\begin{array}{c}\hfill E_0=\langle 0|\widehat{H}|0\rangle =\frac{L^2}{8{\pi }^2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=1}^{\infty }{\omega }_{\mathbf{k},n}.\end{array}$$

(50)

Performing a change of coordinates $(k_x,k_y)$ to polar ones, $(k,\theta )$, and also a change of variable $u=ak$, we obtain

$$\begin{array}{c}\hfill E_0=\frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }d\theta {\int }_0^{\infty }udu\sum _{n=1}^{\infty }{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2+{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{\left(n\pi \right)}^{2ϵ}\right]}^{\frac{1}{2}}.\end{array}$$

(51)

For this case, the integral over the angular variable is trivial. However, to obtain the Casimir energy, we must develop the summation over n. In this sense, by using the Abel–Plana formula, Equation (19), we have not found in the literature a very enlightening result for the integral over the variable t, for any value of $ϵ$. Thus, by adopting an analogous procedure, as in the last subsection, we perform an expansion in the parameter $\frac{l}{a}\ll 1$ in the integrand of (51). Performing this, the leading term in the approximated expression for $E_0$ is given by

$$\begin{array}{ccc}\hfill E_0& \approx & \frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }d\theta {\int }_0^{\infty }udu\sum _{n=1}^{\infty }\{{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\hfill \\ & +& \frac{1}{2}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{\left(n\pi \right)}^{2ϵ}{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}\}.\hfill \end{array}$$

(52)

Again, the first term is associated with the Casimir energy without Lorentz violation. Thus, after integration over the angular coordinate, the second term becomes

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2{(-1)}^{ϵ}}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\sum _{n=1}^{\infty }{\left(n\pi \right)}^{2ϵ}{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(53)

Using the summation Formula (19), we find

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2{(-1)}^{ϵ}}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\left\{-\frac{1}{2}F\left(0\right)+{\int }_0^{\infty }F\left(t\right)\mathrm{d}t+i{\int }_0^{\infty }\frac{F\left(it\right)-F(-it)}{e^{2\pi t}-1}\mathrm{d}t\right\},\end{array}$$

(54)

where

$$\begin{array}{c}\hfill F\left(n\right)={\left(n\pi \right)}^{2ϵ}{\left[u^2+{\left(n\pi \right)}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(55)

Discarding the divergent contributions coming from the first two terms, the LV Casimir energy per unit area of the plates is given by

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& \frac{i}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}v\frac{v^{2ϵ}\left\{{[u^2+{\left(iv\right)}^2+{\left(ma\right)}^2]}^{-1/2}-{[u^2+{(-iv)}^2+{\left(ma\right)}^2]}^{-1/2}\right\}}{e^{2v}-1},\hfill \end{array}$$

(56)

where we have performed a change of variable $t\pi =v$.

Again, analyzing the integral in the variable v over the two intervals $v<{\left[u^2+{\left(ma\right)}^2\right]}^{1/2}$ and $v>{\left[u^2+{\left(ma\right)}^2\right]}^{1/2}$, we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u{\int }_{{\left[u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}}^{\infty }\mathrm{d}v\frac{v^{2ϵ}{[v^2-(u^2+{\left(ma\right)}^2)]}^{-1/2}}{e^{2v}-1}.\end{array}$$

(57)

Next, performing a changing of variable ${\rho }^2=v^2-(u^2+{\left(ma\right)}^2)$ and also a change of coordinates in the plane $(u,\rho )$ to polar ones, we are able to rewrite Equation (57) in the form

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }\frac{{\left[{\sigma }^2+{\left(ma\right)}^2\right]}^{ϵ-\frac{1}{2}}{\sigma }^2\mathrm{d}\sigma }{e^{2{({\sigma }^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}-1}.\end{array}$$

(58)

For the massless scalar field case we have

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }\frac{{\sigma }^{2ϵ+1}\mathrm{d}\sigma }{e^{2\sigma }-1}.\end{array}$$

(59)

Consequently, by using (30), we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{\Gamma (2ϵ+2)\zeta (2ϵ+2)}{2^{(2ϵ+2)}}.\end{array}$$

(60)

Thus, we can analyze two cases:

• For $ϵ=2$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\pi }^4}{2016a^3}{\left(\frac{l}{a}\right)}^2;\end{array}$$

  (61)
• For $ϵ=3$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\pi }^6}{1920a^3}{\left(\frac{l}{a}\right)}^4.\end{array}$$

  (62)

The integral in (58) can only be expressed in terms of an infinite series in modified Bessel functions for $m\ne 0$, so let us evaluate its asymptotic limits. To perform this, we make the following change of variables ${\xi }^2={\sigma }^2+{\left(ma\right)}^2$ and $\xi =mav$. Thus, we have

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_1^{\infty }\frac{v^{2ϵ}{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

(63)

Expressing the denominator in terms of a geometric series as shown in (35), we can rewrite (63) as

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\sum _{j=1}^{\infty }{\int }_1^{\infty }{v}^{2ϵ}{\left(v^2-1\right)}^{\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(64)

The use of the identity,

$$\begin{array}{c}\hfill \frac{1}{{\left(2mj\right)}^{2ϵ}}\frac{{\mathrm{d}}^{2ϵ}\left(e^{-2amvj}\right)}{\mathrm{d}{a}^{2ϵ}}=v^{2ϵ}{e}^{-2amvj},\end{array}$$

(65)

allows us to put Equation (64) in the form

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{1}{{\left(2m\right)}^{2ϵ}}\sum _{j=1}^{\infty }\frac{1}{j^{2ϵ}}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}{\int }_1^{\infty }{\left(v^2-1\right)}^{\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(66)

Again, by making use of the integral representation for the modified Bessel function, the above expression reads

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{1}{{\left(2m\right)}^{2ϵ+1}}\sum _{j=1}^{\infty }\frac{1}{j^{2ϵ+1}}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}\left(\frac{K_1\left(2amj\right)}{a}\right).\end{array}$$

(67)

Let us now consider the asymptotic limits of the above result for $am>>1$ and for $am<<1$:

(i)

For large values of $am>>1$, and using the asymptotic expression (39) for the modified Bessel function, the dominant term provides,

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx \frac{{\left(am\right)}^{2(ϵ+1)}}{8{\left(\pi \right)}^{\frac{3}{2}}{a}^3{m}^{\frac{1}{2}}{\left(2m\right)}^{2ϵ+1}}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}\left(\frac{e^{-2am}}{a^{\frac{3}{2}}}\right).\end{array}$$

(68)

In this sense, we have

• For case $ϵ=2$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx \frac{{\left(am\right)}^{\frac{9}{2}}}{16{\left(\pi \right)}^{\frac{3}{2}}{a}^3}{\left(\frac{l}{a}\right)}^2{e}^{-2am};\end{array}$$

  (69)
• For case $ϵ=3$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx \frac{{\left(am\right)}^{\frac{13}{2}}}{16{\left(\pi \right)}^{\frac{3}{2}}{a}^3}{\left(\frac{l}{a}\right)}^4{e}^{-2am}.\end{array}$$

  (70)

(ii)

For $am<<1$, we have to take the integral representation (63)

• In the case $ϵ=2$, Equation (63) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^6}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{v^4{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

  (71)

  By approximating the integrand as shown below, we obtain

  $$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & \frac{{\left(am\right)}^6}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{\left(v^5-\frac{1}{2}{v}^3-\frac{1}{8}v\right)\mathrm{d}v}{e^{2amv}-1}\hfill \\ & \approx & \frac{1}{2016{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2\left[{\pi }^6-\frac{21}{20}{\pi }^4{\left(am\right)}^2-\frac{21}{8}{\pi }^2{\left(am\right)}^4\right].\hfill \end{array}$$

  (72)

  For this case, the pressure is given by

  $$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=\frac{1}{80,640a^4}{\left(\frac{l}{a}\right)}^2\left[200{\pi }^4-126{\pi }^2{\left(am\right)}^2+105{\left(am\right)}^4\right].\end{array}$$

  (73)
• In the case $ϵ=3$, expression (63) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^8}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{v^6{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}-1}.\end{array}$$

  (74)

Thus, the series expansion is now given by

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & \frac{{\left(am\right)}^8}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{\left(v^7-\frac{1}{2}{v}^5-\frac{1}{8}{v}^3-\frac{1}{16}v\right)\mathrm{d}v}{e^{2amv}-1}\hfill \\ & \approx & \frac{1}{1920{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4\left[{\pi }^8-\frac{10}{21}{\pi }^6{\left(am\right)}^2-\frac{1}{4}{\pi }^4{\left(am\right)}^4\right].\hfill \end{array}$$

(75)

Consequently, the pressure is found to be

$$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=\frac{1}{161,280a^4}{\left(\frac{l}{a}\right)}^4\left[588{\pi }^6-200{\pi }^4{\left(am\right)}^2+63{\pi }^2{\left(am\right)}^4\right].\end{array}$$

(76)

In Figure 3, we present the behavior of the Casimir energy per unit area multiplied by $a^3$, ${\epsilon }_c=\frac{\widetilde{E_C}}{L^2}{a}^3$, as function of $ma$, considering as an illustrative example $\frac{l}{a}=0.01$, for two distinct values of $ϵ$. In the left panel we adopted $ϵ=2$, and in the right panel $ϵ=3$.

3.2. Neumann Condition

We turn now to the Neumann boundary condition. Thus, the field operator solution of Equation (3), which obeys the boundary condition below,

$$\begin{array}{c}\hfill \frac{\partial \varphi \left(\mathbf{x}\right)}{\partial z}{|}_{z=0}=\frac{\partial \varphi \left(\mathbf{x}\right)}{\partial z}{|}_{z=a}=0,\end{array}$$

(77)

has been given in [25]. It reads,

$$\begin{array}{c}\hfill \widehat{\varphi }\left(x\right)=\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }{c}_{n}cos\left({\displaystyle \frac{n\pi }{a}}z\right)[{\widehat{a}}_{\mathbf{k},n}{e}^{-ikx}+{\widehat{a}}_{\mathbf{k},n}^{†}{e}^{ikx}],\end{array}$$

(78)

where

$$c_n=\left\{\begin{array}{cc}{\displaystyle \frac{1}{\sqrt{2{\left(2\pi \right)}^{2}a{\omega }_{\mathbf{k},n}}}}\hfill & \mathrm{for}n=0,\hfill \\ {\displaystyle \frac{1}{\sqrt{{\left(2\pi \right)}^{2}a{\omega }_{\mathbf{k},n}}}}\hfill & \mathrm{for}n\ge 0.\hfill \end{array}\right.$$

(79)

Note that, although the field operator is different from the one presented in the previous subsection, (9), the Hamiltonian operator and the dispersion relations remain the same as to the Dirichlet boundary condition for each choice of the space-like four-vector. Thus, for this reason we decided to not present all the calculations for this case.

3.3. Mixed Boundary Condition

Now, let us consider that the scalar quantum field solution of Equation (3) obeys the Dirichlet boundary condition on one plate and Neumann on the other. In this case, we have two different configurations for the scalar quantum field:

• First configuration,

  $$\begin{array}{c}\hfill \varphi (z=0)=\frac{\partial \varphi \left(\mathbf{x}\right)}{\partial z}{|}_{z=a}=0;\end{array}$$

  (80)
• Second configuration,

  $$\begin{array}{c}\hfill \frac{\partial \varphi \left(\mathbf{x}\right)}{\partial z}{|}_{z=0}=\varphi (z=a)=0.\end{array}$$

  (81)

For this case, the field operators are given by

$$\begin{array}{c}\hfill {\widehat{\varphi }}_{\left(i\right)}\left(x\right)=\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }\frac{1}{\sqrt{{\left(2\pi \right)}^{2}a{\omega }_{\mathbf{k},n}}}sin\left((n+1/2)\frac{\pi }{a}z\right)[{\widehat{a}}_{\mathbf{k},n}{e}^{-ikx}+{\widehat{a}}_{\mathbf{k},n}^{†}{e}^{ikx}]\end{array}$$

(82)

for the first configuration and

$$\begin{array}{c}\hfill {\widehat{\varphi }}_{\left(ii\right)}\left(x\right)=\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }\frac{1}{\sqrt{{\left(2\pi \right)}^{2}a{\omega }_{\mathbf{k},n}}}cos\left((n+1/2)\frac{\pi }{a}z\right)[{\widehat{a}}_{\mathbf{k},n}{e}^{-ikx}+{\widehat{a}}_{\mathbf{k},n}^{†}{e}^{ikx}],\end{array}$$

(83)

for the second configuration. However, the field operators ${\widehat{\varphi }}_1$ and ${\widehat{\varphi }}_2$ provide the same Hamiltonian operator and dispersion relations.

3.3.1. Vector Parallel to the Plates

Considering $u^{\mu }=(0,1,0,0)$, the field presents the following dispersion relation,

$$\begin{array}{c}\hfill {\omega }_{\mathbf{k},n}^2=k_x^2+k_y^2+l^{2(ϵ-1)}{(-1)}^{ϵ}{k}_x^{2ϵ}+{\left[(n+1/2)\frac{\pi }{a}\right]}^2+m^2.\end{array}$$

(84)

Note that both field operators, ${\widehat{\varphi }}_{\left(i\right)}\left(x\right)$ and ${\widehat{\varphi }}_{\left(ii\right)}\left(x\right)$, provide the same Hamiltonian operator, i.e.,

$$\begin{array}{c}\hfill \widehat{H}=\frac{1}{2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }{\omega }_{\mathbf{k},n}\left[2{\widehat{a}}_{\mathbf{k},\mathbf{n}}^{†}{\widehat{a}}_{\mathbf{k},\mathbf{n}}+\frac{L^2}{{\left(2\pi \right)}^2}\right].\end{array}$$

(85)

Consequently, the corresponding vacuum energy is expressed as

$$\begin{array}{c}\hfill E_0=\langle 0|\widehat{H}|0\rangle =\frac{L^2}{8{\pi }^2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }{\omega }_{\mathbf{k},n}.\end{array}$$

(86)

By performing a change of coordinates in the plane $(k_x,k_y)$ to the polar ones, making the change of variable $u=ak$ and performing an expansion in the dimensionless parameter $\frac{l}{a}<<1$, up to the first order in this parameter, we find

$$\begin{array}{ccc}\hfill E_0& \approx & \frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }\mathrm{d}\theta {\int }_0^{\infty }u\mathrm{d}u\sum _{n=0}^{\infty }\{{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\hfill \\ & +& \frac{1}{2}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{u}^{2ϵ}{cos}^{2ϵ}\theta {\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}\}.\hfill \end{array}$$

(87)

Note that the first term on the right side of the above expression is associated with the vacuum energy without Lorentz violation. An integral representation for the corresponding Casmir energy by unity area for this case [25], reads

$$\begin{array}{c}\hfill \frac{E_C}{L^2}=\frac{a{m}^4}{6{\pi }^2}{\int }_1^{\infty }\frac{{(v^2-1)}^{\frac{3}{2}}}{e^{2amv}+1}dv.\end{array}$$

(88)

As to the second term, after integration over the angular coordinate, we have

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\sum _{n=0}^{\infty }{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(89)

In order to develop the summation over the half-integer number, we will use the Abel–Plana formula below [33]:

$$\begin{array}{c}\hfill \sum _{n=0}^{\infty }F\left(n+\frac{1}{2}\right)={\int }_0^{\infty }F\left(t\right)\mathrm{d}t-i{\int }_0^{\infty }\frac{\mathrm{d}t}{e^{2\pi t}+1}[F\left(it\right)-F(-it)].\end{array}$$

(90)

Then, the expression (89) becomes

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\left\{{\int }_0^{\infty }F\left(t\right)\mathrm{d}t-i{\int }_0^{\infty }\frac{F\left(it\right)-F(-it)}{e^{2\pi t}+1}\mathrm{d}t\right\},\end{array}$$

(91)

where

$$\begin{array}{c}\hfill F\left(n+\frac{1}{2}\right)={\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(92)

The first term on the right side of (91) refers to the free vacuum energy, so it is discarded in the renormalization process. Then, the LV Casimir energy will be given by

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{i}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}t\frac{{[u^2+{\left(it\pi \right)}^2+{\left(ma\right)}^2]}^{-1/2}-{[u^2+{(-it\pi )}^2+{\left(ma\right)}^2]}^{-1/2}}{e^{2\pi t}+1}.\hfill \end{array}$$

(93)

Performing a change of variable, where $t\pi =v$, we can rewrite Equation (93) in the form

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{i}{8{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}v\frac{{[u^2+{\left(ma\right)}^2+{\left(iv\right)}^2]}^{-1/2}-{[u^2+{\left(ma\right)}^2+{(-iv)}^2]}^{-1/2}}{e^{2v}+1}.\hfill \end{array}$$

(94)

Now, carrying out the same analysis as in the previous cases, for the interval of the integral in v, we obtain

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_{{\left[u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}}^{\infty }\mathrm{d}v\frac{{[v^2-(u^2+{\left(ma\right)}^2)]}^{-1/2}}{e^{2v}+1}.\hfill \end{array}$$

(95)

Additionally, the new change of variable ${\rho }^2=v^2-(u^2+{\left(ma\right)}^2)$ provides

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}(2ϵ-1)!!}{\left(2ϵ\right)!!}{\int }_0^{\infty }{u}^{{}^{(2ϵ+1)}}\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\frac{\mathrm{d}\rho }{{\left[{\rho }^2+u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\left[e^{2{({\rho }^2+u^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}+1\right]}.\hfill \end{array}$$

(96)

Finally, by changing the coordinates in the plane $(u,\rho )$ to polar ones, we find

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_0^{\infty }\frac{{\sigma }^{2(ϵ+1)}\mathrm{d}\sigma }{{\left[{\sigma }^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\left[e^{2{({\sigma }^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}+1\right]}.\end{array}$$

(97)

The vacuum energy above provides the massless scalar field case by setting $m=0$. This gives

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_0^{\infty }\frac{{\sigma }^{2ϵ+1}\mathrm{d}\sigma }{e^{2\sigma }+1}.\end{array}$$

(98)

Using the result of the integral below [34]

$$\begin{array}{c}\hfill {\int }_0^{\infty }\frac{x^{\nu -1}\mathrm{d}x}{e^{\mu x}+1}=\frac{1}{{\mu }^{\nu }}\left(1-2^{1-\nu }\right)\Gamma \left(\nu \right)\zeta \left(\nu \right),\end{array}$$

(99)

where $\Gamma \left(\nu \right)$ corresponds to the Gamma function and $\zeta \left(\nu \right)$ is the Riemann zeta function, we arrive at

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}\frac{\left(1-2^{-(2ϵ+1)}\right)\Gamma (2ϵ+2)\zeta (2ϵ+2)}{2^{(2ϵ+2)}}.\end{array}$$

(100)

Consequently, we can exhibit two cases:

• For $ϵ=2$:

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{31{\pi }^4}{322,560a^3}{\left(\frac{l}{a}\right)}^2;\end{array}$$

(101)

• For $ϵ=3$:

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{127{\pi }^6}{1,720,320a^3}{\left(\frac{l}{a}\right)}^4.\end{array}$$

(102)

Because the integral in (97) cannot be expressed in terms of elementary functions, let us evaluate its asymptotic limits. To perform this, we make the following changes of variables ${\xi }^2={\sigma }^2+{\left(ma\right)}^2$ and $\xi =mav$; thus, we have

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{ϵ+\frac{1}{2}}\mathrm{d}v}{e^{2amv}+1}.\end{array}$$

(103)

Knowing that the geometric series can be represented as

$$\begin{array}{c}\hfill \frac{1}{e^{2amv}+1}=\sum _{j=1}^{\infty }{\left(-1\right)}^{j+1}{e}^{-2amvj},\end{array}$$

(104)

we obtain,

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{(-1)}^{ϵ}}{(2ϵ+1)}\sum _{j=1}^{\infty }{\left(-1\right)}^{j+1}{\int }_1^{\infty }{\left(v^2-1\right)}^{ϵ+\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(105)

Finally, by using the integral representation of the modified Bessel function, the expression above becomes

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{ϵ+1}\Gamma \left(ϵ+\frac{3}{2}\right){(-1)}^{ϵ}}{4{\left(\pi \right)}^{\frac{5}{2}}{a}^3(2ϵ+1)}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\sum _{j=1}^{\infty }\frac{{\left(-1\right)}^{j+1}}{j^{ϵ+1}}{K}_{ϵ+1}\left(2amj\right).\end{array}$$

(106)

We can now examine the asymptotic limits $ma>>1$ and $ma<<1$. They are:

(i)

For large arguments, $am>>1$, the modified Bessel can be expressed in an exponential form as shown in (39); consequently the dominant term, $j=1$, provides

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx -\frac{{\left(am\right)}^{ϵ+\frac{1}{2}}\Gamma \left(ϵ+\frac{3}{2}\right){(-1)}^{ϵ}}{8{\pi }^2{a}^3(2ϵ+1)}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{e}^{-2am};\end{array}$$

(107)

(ii)

In order to analyze the case for $am<<1$, let us consider $ϵ=2,3$:

• In the case $ϵ=2$, expression (103) becomes

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^6}{20{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{\frac{5}{2}}\mathrm{d}v}{e^{2amv}+1}.\end{array}$$

  (108)

  By expanding the integrand in a series of positive powers of v, we can obtain an approximated expression for the LV Casimir energy per unity area, as shown below:

  $$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & -\frac{{\left(am\right)}^6}{20{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{\left(v^5-\frac{5}{2}{v}^3+\frac{15}{8}v\right)\mathrm{d}v}{e^{2amv}+1}\hfill \\ & \approx & -\frac{31}{322,560{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2\left[{\pi }^6-\frac{147}{31}{\pi }^4{\left(am\right)}^2+\frac{630}{31}{\pi }^2{\left(am\right)}^4\right].\hfill \end{array}$$

  (109)

  As for the pressure, we have,

  $$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=-\frac{1}{322,560a^4}{\left(\frac{l}{a}\right)}^2\left[155{\pi }^4-441{\pi }^2{\left(am\right)}^2+630{\left(am\right)}^4\right].\end{array}$$

  (110)
• In the case $ϵ=3$, expression (103) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=\frac{{\left(am\right)}^8}{28{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{{\left(v^2-1\right)}^{\frac{7}{2}}\mathrm{d}v}{e^{2amv}+1}.\end{array}$$

  (111)

  Adopting a similar procedure as above, we have:

  $$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & \frac{{\left(am\right)}^8}{28{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{\left(v^7-\frac{7}{2}{v}^5+\frac{35}{8}{v}^3-\frac{35}{16}v\right)\mathrm{d}v}{e^{2amv}+1}\hfill \\ & \approx & \frac{127}{1,720,320{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4\left[{\pi }^8-\frac{1240}{381}{\pi }^6{\left(am\right)}^2+\frac{980}{127}{\pi }^4{\left(am\right)}^4\right].\hfill \end{array}$$

  (112)

  In this case, the pressure is given by

  $$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=\frac{1}{5,160,960a^4}{\left(\frac{l}{a}\right)}^4\left[2667{\pi }^6-6200{\pi }^4{\left(am\right)}^2+8820{\pi }^2{\left(am\right)}^4\right].\end{array}$$

  (113)

  In Figure 4 we present the behavior of the Casimir energy per unit area as function of $ma$, considering again as an illustrative example $\frac{l}{a}=0.01$, for two distinct values of $ϵ$. We consider $ϵ=2$ in the left plot and $ϵ=3$ in the right plot.

3.3.2. Vector Perpendicular to the Plates

Now, let us consider the four-vector $u^{\mu }$ as being perpendicular to the plates by choosing

$$\begin{array}{c}\hfill u^{\mu }=(0,0,0,1).\end{array}$$

(114)

For this case, the corresponding dispersion relation is

$$\begin{array}{c}\hfill {\omega }_{\mathbf{k},n}^2=k_x^2+k_y^2+{\left[(n+1/2)\frac{\pi }{a}\right]}^2+l^{2(ϵ-1)}{(-1)}^{ϵ}{\left[(n+1/2)\frac{\pi }{a}\right]}^{2ϵ}+m^2.\end{array}$$

(115)

Consequently, the Hamiltonian operator, $\widehat{H}$, is now given by

$$\begin{array}{c}\hfill \widehat{H}=\frac{1}{2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }{\omega }_{\mathbf{k},n}\left[2{\widehat{a}}_{\mathbf{k},\mathbf{n}}^{†}{\widehat{a}}_{\mathbf{k},\mathbf{n}}+\frac{L^2}{{\left(2\pi \right)}^2}\right].\end{array}$$

(116)

The corresponding vacuum energy in this case is given by

$$\begin{array}{c}\hfill E_0=\langle 0|\widehat{H}|0\rangle =\frac{L^2}{8{\pi }^2}\int {\mathrm{d}}^2\mathbf{k}\sum _{n=0}^{\infty }{\omega }_{\mathbf{k},n}.\end{array}$$

(117)

Performing a change of coordinates in the plane $(k_x,k_y)$ to polar ones, a change of variable $u=ak$ and performing an expansion in $\frac{l}{a}<<1$, we obtain

$$\begin{array}{ccc}\hfill E_0& \approx & \frac{L^2}{8{\pi }^2{a}^3}{\int }_0^{2\pi }\mathrm{d}\theta {\int }_0^{\infty }u\mathrm{d}u\sum _{n=0}^{\infty }\{{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}\hfill \\ & +& \frac{1}{2}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{(-1)}^{ϵ}{\left[\left(n+1/2\right)\pi \right]}^{2ϵ}{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}\}.\hfill \end{array}$$

(118)

The first term on the right side is associated with energy without Lorentz violation. Thus, after integration over the angular coordinate, the second term becomes

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2{(-1)}^{ϵ}}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\sum _{n=0}^{\infty }{\left[(n+1/2)\pi \right]}^{2ϵ}{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(119)

By using Formula (90) again, we find

$$\begin{array}{c}\hfill \widetilde{E_0}=\frac{L^2{(-1)}^{ϵ}}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\left\{{\int }_0^{\infty }F\left(t\right)\mathrm{d}t-i{\int }_0^{\infty }\frac{F\left(it\right)-F(-it)}{e^{2\pi t}+1}\mathrm{d}t\right\},\end{array}$$

(120)

where

$$\begin{array}{c}\hfill F\left(n+\frac{1}{2}\right)={\left[(n+1/2)\pi \right]}^{2ϵ}{\left[u^2+{\left[(n+1/2)\pi \right]}^2+{\left(ma\right)}^2\right]}^{-\frac{1}{2}}.\end{array}$$

(121)

Note that the integral in the first term on the right side of Equation (120) is the free vacuum contribution, whereas the integral in the second term gives the renormalized vacuum energy per unit area, i.e.,

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{i}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}t\frac{{\left(t\pi \right)}^{2ϵ}{[u^2+{\left(it\pi \right)}^2+{\left(ma\right)}^2]}^{-1/2}-{\left(t\pi \right)}^{2ϵ}{[u^2+{(-it\pi )}^2+{\left(ma\right)}^2]}^{-1/2}}{e^{2\pi t}+1}.\hfill \end{array}$$

(122)

By making the change of variable $t\pi =v$, we find

$$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& =& -\frac{i}{8\pi a^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u\hfill \\ & \times & {\int }_0^{\infty }\mathrm{d}t\frac{v^{2ϵ}\left\{{[u^2+{\left(iv\right)}^2+{\left(ma\right)}^2]}^{-1/2}-{[u^2+{(-iv)}^2+{\left(ma\right)}^2]}^{-1/2}\right\}}{e^{2v}+1}.\hfill \end{array}$$

(123)

Again, analyzing the integral on the variable v over the two intervals: ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}>v$ and ${\left[u^2+{\left(ma\right)}^2\right]}^{1/2}<v$, we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }u\mathrm{d}u{\int }_{{\left[u^2+{\left(ma\right)}^2\right]}^{\frac{1}{2}}}^{\infty }\mathrm{d}v\frac{v^{2ϵ}{[v^2-(u^2+{\left(ma\right)}^2)]}^{-1/2}}{e^{2v}+1}.\end{array}$$

(124)

By making an additional change of variable ${\rho }^2=v^2-(u^2+{\left(ma\right)}^2)$, and transforming the coordinates in the plane $(u,\rho )$ to polar ones, we arrive at

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }\frac{{\left[{\sigma }^2+{\left(ma\right)}^2\right]}^{ϵ-\frac{1}{2}}{\sigma }^2\mathrm{d}\sigma }{e^{2{({\sigma }^2+{\left(ma\right)}^2)}^{\frac{1}{2}}}+1}.\end{array}$$

(125)

For the case of the massless field, we take $m=0$. In this case, the integral reads,

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_0^{\infty }\frac{{\sigma }^{2ϵ+1}\mathrm{d}\sigma }{e^{2\sigma }+1}.\end{array}$$

(126)

By considering the same integral as in the last subsection, we find

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{1}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{\left(1-2^{-(2ϵ+1)}\right)\Gamma (2ϵ+2)\zeta (2ϵ+2)}{2^{(2ϵ+2)}}.\end{array}$$

(127)

Let us analyze the cases for which $ϵ=2,3$:

• For $ϵ=2$:

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{31{\pi }^4}{64,512a^3}{\left(\frac{l}{a}\right)}^2;\end{array}$$

(128)

• For $ϵ=3$:

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{127{\pi }^6}{245,760a^3}{\left(\frac{l}{a}\right)}^4.\end{array}$$

(129)

Once again, let us evaluate the integral in (125) in its asymptotic limits. To conduct this, we make the following changes of variables ${\xi }^2={\sigma }^2+{\left(ma\right)}^2$ and $\xi =mav$; thus, we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}{\int }_1^{\infty }\frac{v^{2ϵ}{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}+1}.\end{array}$$

(130)

Expressing the denominator in a geometric series as given in (104), we obtain

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\sum _{j=1}^{\infty }{\left(-1\right)}^{j+1}{\int }_1^{\infty }{v}^{2ϵ}{\left(v^2-1\right)}^{ϵ+\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(131)

In addition, using the identity below,

$$\begin{array}{c}\hfill \frac{1}{{\left(2mj\right)}^{2ϵ}}\frac{{\mathrm{d}}^{2ϵ}\left(e^{-2amvj}\right)}{\mathrm{d}{a}^{2ϵ}}=v^{2ϵ}{e}^{-2amvj},\end{array}$$

(132)

the LV Casimir energy can be expressed as

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{1}{{\left(2m\right)}^{2ϵ}}\sum _{j=1}^{\infty }\frac{{\left(-1\right)}^{j+1}}{j^{2ϵ}}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}{\int }_1^{\infty }{\left(v^2-1\right)}^{\frac{1}{2}}{e}^{-2amvj}\mathrm{d}v.\end{array}$$

(133)

Using the integral representation for the modified Bessel function, the above expression becomes

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^{2(ϵ+1)}}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{1}{{\left(2m\right)}^{2ϵ+1}}\sum _{j=1}^{\infty }\frac{{\left(-1\right)}^{j+1}}{j^{2ϵ+1}}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}\left(\frac{K_1\left(2amj\right)}{a}\right).\end{array}$$

(134)

Let us now analyze the asymptotic limits for $ma>>1$ and $ma<<1$:

(i)

For large arguments, $am>>1$, we can use the asymptotic form for the modified Bessel function, Equation (39), and taking the dominant term, $j=1$, we obtain that

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx -\frac{{\left(am\right)}^{2(ϵ+1)}}{8{\left(\pi \right)}^{\frac{3}{2}}{a}^3{m}^{\frac{1}{2}}{\left(2m\right)}^{2ϵ+1}}{\left(\frac{l}{a}\right)}^{2(ϵ-1)}\frac{{\mathrm{d}}^{2ϵ}}{\mathrm{d}{a}^{2ϵ}}\left(\frac{e^{-2am}}{a^{\frac{3}{2}}}\right).\end{array}$$

(135)

We additionally want to consider the two cases $ϵ=2,3$:

• For case $ϵ=2$:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx -\frac{{\left(am\right)}^{\frac{9}{2}}}{16{\left(\pi \right)}^{\frac{3}{2}}{a}^3}{\left(\frac{l}{a}\right)}^2{e}^{-2am};\end{array}$$

  (136)
• For case $ϵ=3$:

$$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}\approx -\frac{{\left(am\right)}^{\frac{13}{2}}}{16{\left(\pi \right)}^{\frac{3}{2}}{a}^3}{\left(\frac{l}{a}\right)}^4{e}^{-2am}.\end{array}$$

(137)

(ii)

For $am<<1$:

• In the case $ϵ=2$, expression (130) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^6}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{v^4{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}+1}.\end{array}$$

  (138)

  Expanding the integrand in powers of v, we obtain an expression that allows us to evaluate the integral:

  $$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & -\frac{{\left(am\right)}^6}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2{\int }_1^{\infty }\frac{\left(v^5-\frac{1}{2}{v}^3-\frac{1}{8}v\right)\mathrm{d}v}{e^{2amv}+1}\hfill \\ & \approx & -\frac{31}{64,512{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^2\left[{\pi }^6-\frac{147}{155}{\pi }^4{\left(am\right)}^2-\frac{42}{31}{\pi }^2{\left(am\right)}^4\right].\hfill \end{array}$$

  (139)

  As for the pressure, we have,

  $$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=-\frac{1}{322,560a^4}{\left(\frac{l}{a}\right)}^2\left[775{\pi }^4-441{\pi }^2{\left(am\right)}^2-210{\left(am\right)}^4\right].\end{array}$$

  (140)
• In the case $ϵ=3$, expression (130) becomes:

  $$\begin{array}{c}\hfill \frac{\widetilde{E_C}}{L^2}=-\frac{{\left(am\right)}^8}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{v^6{\left(v^2-1\right)}^{\frac{1}{2}}\mathrm{d}v}{e^{2amv}+1},\end{array}$$

  (141)

  the series expansion in this case provides

  $$\begin{array}{ccc}\hfill \frac{\widetilde{E_C}}{L^2}& \approx & -\frac{{\left(am\right)}^8}{4{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4{\int }_1^{\infty }\frac{\left(v^7-\frac{1}{2}{v}^5-\frac{1}{8}{v}^3-\frac{1}{16}v\right)\mathrm{d}v}{e^{2amv}+1}\hfill \\ & \approx & -\frac{127}{245,760{\pi }^2{a}^3}{\left(\frac{l}{a}\right)}^4\left[{\pi }^8-\frac{1240}{2667}{\pi }^6{\left(am\right)}^2-\frac{28}{127}{\pi }^4{\left(am\right)}^4\right].\hfill \end{array}$$

  (142)

Consequently, the pressure reads,

$$\begin{array}{c}\hfill {\tilde{P}}_C\left(a\right)=-\frac{1}{5,160,960a^4}{\left(\frac{l}{a}\right)}^4\left[18,669{\pi }^6-6200{\pi }^4{\left(am\right)}^2-1764{\pi }^2{\left(am\right)}^4\right].\end{array}$$

(143)

Once more, as an illustrative example, in Figure 5 we present the behavior of the Casimir energy per unit of area as function of $ma$, considering $\frac{l}{a}=0.01$, for two distinct values of $ϵ$. In the left panel, we consider $ϵ=2$, whereas in the right panel $ϵ=3$.

4. Concluding Remarks

In this work, we investigated the influence of the Lorentz symmetry violation on the Casimir energy associated with a real massive scalar quantum field. We considered the situation in which the field is confined between two parallel plates and assumed that the field obeys boundary conditions of the types of Dirichlet and Neumann, and is mixed on the plates of area $L^2$ separated by the distance $a(a\ll L)$.

The Lorentz symmetry violation is implemented, admitting a direct coupling between a constant space-like four-vector, $u^{\mu }$, in a scenario aether-like CPT-even theoretical model, with higher-order derivative of the field, represented by $l^{2(ϵ-1)}{(u\partial )}^{2ϵ}\varphi \left(x\right)$, as exhibited in the modified Klein–Gordon equation, Equation (3), being $ϵ$ an integer greater or equal to 2, and l is a parameter of order inverse the energy scale where the Lorentz symmetry is broken. As for the constant four-vector, two distinct directions are considered: the vector parallel to the plates and the vector perpendicular to them. We have verified that the combined modifications in the dynamic of the quantum field produce important corrections on the corresponding dispersion relations, and consequently on the Casimir energies. In the calculation of the Casimir energy, we have to develop an integral over the bidimensional space in the $\overrightarrow{k}$ plane associated with the momentum of the field parallel to the plates, and a sum over discrete momentum orthogonal to them, as exhibited in Equation (14), for example. To develop the summation over the discrete momentum, we adopted the Abel–Plana summation formula for integer, Equation (19), and half-integer, Equation (90); quantum numbers, respectively. Because the integrals involved in the obtainment of the Casimir energy do not provide very enlightening results, in our analysis we decided to develop an expansion on the parameter $\frac{l}{a}<<1$ in the integrand, keeping up to the first order term. Adopting this procedure, we were able to provide the corrections on the Casimir energy and pressure due to the LV term, caused by specific boundary conditions obeyed by the quantum field on the plates, the direction of the space-like constant vector, and the order of the derivative. Because the dispersion relation for the Neumann boundary condition is analogous to the one in the Dirichlet condition, we only provide a brief discussion for this case. By our results for the Casimir pressures induced by the LV term, ${\tilde{P}}_c\left(a\right)$, we would like to emphasize that negative values that correspond to attracted forces between the plates, and positive values that correspond to repulsive forces between the plates, depend on the direction of the constant vector, $u^{\mu }$, with respect to the plates, on the value assumed for the parameter $ϵ$ and also depend on the boundary condition imposed on the field at the plates.

In our analysis, the LV Casimir energy is expressed in terms of an integral representation for a massive field. Thus, in order to furnish some quantitative information about this energy, we provided its asymptotic expressions for $am>>1$ and $am<<1$. In the former, the Casimir energy decays exponentially as $e^{-2am}$, whereas in the opposite limit, it presents a term that corresponds to the massless case with additional corrections proportional to some power of the product $am$. In addition, we also presented four graphs for the Casimir energies as a function of $am$ considering $ϵ$ equal to 2 and 3, contemplating all the possible scenarios. Of course, the intensity of the LV Casimir energy depends on the order of the higher-derivative term. It is smaller for a higher value of $ϵ$. Considering the Dirichlet condition, and the vector parallel to the plates, the LV Casimir energy and pressure are positive for $ϵ=2$ and negative for $ϵ=3$; however, for vector perpendicular to the plates, the LV Casimir energy and pressure present the same positive sign for both values of $ϵ$. For a mixed boundary condition, the same behavior related to the sign of the LV Casimir energies and pressure is observable. We would like to say that there are changes in the sign for the Casimir energies when the field obeys the Dirichlet boundary condition and mixed one for each specific situation.

To finish this section, we want to make a few comments about the results obtained in Section 3.1 and Section 3.3. The corrections induced by the Lorentz violation in the Casimir pressure never vanish. Accepting that the Lorentz violation is part of the source in the $1\%$ experimental error estimated in [36], it is possible to infer an upper bound for the parameter l. Considering that the distance a between the parallel plates are of order ${10}^{-8}$ m, the upper limit for l is of order ${10}^{-9}$ m for the case of $ϵ=2$ and of order ${10}^{-8}$ m for the case of $ϵ=3$. In this paper, we have analyzed the Casimir effect associated with the scalar field in a general scenario of LV. In this sense, the results obtained for the LV Casimir energies include two different ingredients: the presence of a background constant vector and spatial higher-order derivative terms. The analysis of the Casimir energy associated with electromagnetic fields in an LV scenario was developed in [37]. In this work, the corrections on the standard Casimir energy due to two different LV approaches are discussed.