This site uses cookies. By continuing to use this site you agree to our use of cookies. To find out more, see our Privacy and Cookies policy.
Brought to you by:
Paper

Casimir effects in systems containing 2D layers such as graphene and 2D electron gases

Published 12 May 2015 © 2015 IOP Publishing Ltd
, , Citation B E Sernelius 2015 J. Phys.: Condens. Matter 27 214017 DOI 10.1088/0953-8984/27/21/214017

0953-8984/27/21/214017

Abstract

We present a variety of methods to derive the Casimir interaction in planar systems containing 2D layers. Examples where this can be of use is graphene, graphene-like layers and 2D electron gases. We present results for two free standing layers and for one layer above a substrate. The results can easily be extended to systems with a larger number of layers.

Export citation and abstract BibTeX RIS

1. Introduction

The dispersion interaction is a correlation effect. The interaction is a result of the correlated motion of the charged carriers in different parts of the system. The dispersion energy between two objects is the inter-object correlation energy. In the field of condensed matter one often uses the Coulomb gauge. In the Coulomb gauge the scalar potential is instantaneous and the vector potential is purely transverse. The scalar potential represents Coulomb interactions between the carriers while the vector potential represents photon interactions. If the two objects are not too far away from each other one may neglect retardation effects. Retardation effects are effects from the finite speed of light in vacuum. A rough estimate of the separation, d, between the objects where the retardation effects start to become important is d ≈ c/ω0, where ω0 is some dominating frequency of the objects (eg the plasma frequency for metallic objects). If retardation effects can be neglected only the Coulomb interactions contribute and the properties of the materials enter in terms of the density–density correlation functions. If retardation has to be taken into account the photon interactions have to be included and, in this additional contribution, the properties of the materials enter in terms of the current–current correlation functions. The interaction can be found using many-body theory [1, 2] with diagrammatic perturbation theory based on Feynman diagrams.

Alternatively the interaction can be viewed as caused by fluctuations: fluctuations in the charge and current densities or fluctuations in the electromagnetic fields. Then one way to find the interaction is in terms of the electromagnetic normal modes of the system [3]. These normal modes are massless bosons and at zero temperature the interaction energy is the sum of the zero-point energy of all these modes.

In a translational invariant and homogeneous system there are three types of normal mode; two with transverse fields and one with longitudinal fields. In the present system there are only two types of normal mode; one with p-polarized fields, i.e. with electric vector in the plane formed by the in-plane momentum k and the normal to the layers, experienced as longitudinal in the two-dimensional (2D) sheets; and one with s-polarized fields, i.e. with electric vector perpendicular to this plane, experienced as transverse in the 2D sheets. The first type of mode is called TM (transverse magnetic) and the second TE (transverse electric). The dispersion curves for each mode type can have one or several branches.

For planar structures the quantum number that characterizes the normal modes is k, the 2D wave vector in the plane of the interfaces or planar structures and the interaction energy per unit area for each mode type can be written as

Equation (1)

where A is the area of a planar interface and the summation over i is a summation over the branches. To get the expression on the right-hand side we have let A be large so that the summation over the discrete k variable can be replaced by an integral over a continuous k variable.

In a simple system with a small number of well-defined modes this summation may be performed directly. In most cases it is more complicated. The complications can e.g. be that the modes form continua or that it is difficult to find the zero-point energies explicitly. An extension of the so-called Argument Principle [35] can then be used to find the results. Let us study a region in the complex frequency plane where two functions are defined: one, φ (z), is analytic in the whole region; and one, f(z), has poles and zeros inside the region. The following relation holds for an integration path around the region:

Equation (2)

where z0 and z are the zeros and poles, respectively, of function f(z). If we choose the function f(z) to be the function in the defining equation for the normal modes of the system, f(ωi) = 0, the function φ(z) as z/2 and let the contour enclose all the zeros and poles of the function f(z), then (2) produces the energy in (1). After deforming the integration path in the complex frequency plane in a proper way one arrives at the expression [3]

Equation (3)

where ${f_{{\bf k}}}\left( {\omega _{{\bf k}}^i} \right) = 0$ is the condition for electromagnetic normal modes. Equation (3) is valid for zero temperature and the interaction energy is the internal energy.

At finite temperature the thermal population of the modes also affects the interaction. The interaction energy is now Helmholtz free energy and can be written as

Equation (4)

where β = 1/kBT. The integral over frequency has been replaced by a summation over discrete frequencies, the so-called Matsubara frequencies [1, 2]. The prime on the summation sign indicates that the n = 0 term should be divided by two. Note that the temperature enters in the discrete frequencies ξn and in the temperature dependence of the dielectric functions.

The normal modes can be found in different ways. One way is to solve Maxwell's equations (in an absence of any external sources) in all regions of the system and use the proper boundary conditions at the interfaces. Another way is to find self-sustained solutions using induced charge- and current-densities in the objects.

We here derive the dispersion interaction in a variety of ways. Since the non-retarded treatment is much easier to perform we present both the non-retarded and fully retarded derivations in separate sections. In the derivations involving normal modes we derive the mode condition function, fk(ω). The interaction energy is then found by inserting this function into (3) or (4). In the derivations involving diagrammatic perturbation theory the interaction energy is obtained directly.

We do not present any numerical results, only the final expressions and different ways to arrive at these. However, we give the analytical expressions for the zero-temperature dielectric function of graphene and of the 2D electron gas to enable the reader to perform calculations for other geometries. The finite-temperature dielectric functions cannot be found in simple analytic form and due to the limited space we refrain from discussing how to proceed at finite temperature. Numerical results for graphene at zero and finite temperature can be found in the literature [624]. These references are further discussed in section 7.

We begin in section 2 by introducing our notation. The following two sections are devoted to background materials, electromagnetic in section 3 and 2D in section 4. The main results when retardation effects are neglected are derived in section 5 and when retardation effects are taken into account in section 6. We end in section 7 with a summary and some conclusions.

2. Notation

In this work we discuss planar structures. We let all surfaces, interfaces and 2D films be parallel to the xy-plane. We let q be a three-dimensional (3D) wave vector and k its 2D in-plane component. The corresponding notation for the spatial coordinates are R and r, respectively. Thus we have

Equation (5)

We will end up with a formalism depending only on k and r apart from the distances between the various interfaces. Thus the formalism is 2D.

We will make extensive use of the Fourier transform, both the 3D and 2D versions and since they can be defined in different ways we list our definitions here. The 3D Fourier transform and inverse Fourier transform of a function F with only spatial dependence are

Equation (6)

If the function also has a temporal dependence the corresponding relations are

Equation (7)

In 2D the Fourier transform and inverse Fourier transform of a function are

Equation (8)

if the function has only spatial dependence. If the function also has a temporal dependence the relations are

Equation (9)

We will need yet another Fourier transform. We let the function F be a function of the 3D spatial variable R but we only transform with respect to r. Then we have for a function with only spatial dependence

Equation (10)

If the function also has a temporal dependence the relations are

Equation (11)

Note that we have refrained from introducing separate function names for the different Fourier transformed versions of a function. By looking at the argument of the function we can tell if it is a 3D or 2D Fourier transform, a partial transform, or not a transform at all.

3. Electromagnetic background material

This section has two subsections. In the first, section 3.1, we present the proper version of Maxwell's equations and also the boundary conditions at interfaces in the system. In the second, section 3.2, we give the scalar and vector potentials in Coulomb and Lorentz gauges and show how the fields are obtained from these.

3.1. Maxwell's equations

There are different formulations of electromagnetism in the literature. The difference lies in how the conduction carriers are treated. In one formulation these carriers are lumped together with the external charges to form the group of free charges. Then only the bound charges contribute to the screening. We want to be able to treat geometries with metallic regions. Then this formulation is not suitable. In the formulation that we use, the conduction carriers are treated on the same footing as the bound charges. Thus, both bound and conduction charges contribute to the dielectric function. In the two formalisms the E and B fields, the true fields, are of course the same. However, the auxiliary fields, the D and H fields, are different. To indicate this we use an alternative formulation in which we put a tilde above the D and H fields and also above the dielectric functions, polarizabilities and conductivities.

With this formulation Maxwell's equations (MEs) have the form

Equation (12)

and the boundary conditions at an interface between the two media 1 and 2 are

Equation (13)

where ρ and J are volume charge and current densities, respectively, while ρs and K are surface charge and current densities, respectively. We note that the sources to the fields in MEs are the external charge and current densities. Also, in the boundary conditions discontinuities in the normal component of the ${{\bf \tilde D}}$ fields and tangential component of the ${\bf {\tilde H}}$ fields are caused by external surface charge densities and external surface current densities, respectively. The unit vector n is the surface normal pointing into region 2.

3.2. Potentials

Often it is easier to deal with potentials instead of the fields. The potentials are auxiliary functions, or help functions. The potentials are not unique, they depend on which gauge we have chosen. Two very common gauges used in condensed matter theory are the Coulomb gauge (mentioned in the Introduction) and the Lorentz gauge. These are obtained by using different gauge conditions to decouple the two coupled inhomogeneous differential equations relating the potentials to the sources.

In the Coulomb gauge the gauge condition is ∇ · A = 0 and the resulting potentials are

Equation (14)

Here the scalar potential, Φ(R, t), is instantaneous; it describes the Coulomb interaction between the charged carriers; the Coulomb interaction is purely longitudinal. The vector potential, A(R, t), is retarded and transverse; it describes the photon interactions between the charged particles; the photon interactions are purely transverse.

In the Lorentz gauge the gauge condition is ∇ · A + (1/c)∂Φ/∂t = 0 and

Equation (15)

These relations are Lorentz covariant, ie they look the same in all inertial frames. The potentials are both solutions to the inhomogeneous wave equation, where the spatial and temporal coordinates are treated on the same footing. This gauge is the obvious choice in relativistic problems. Here both potentials are retarded.

In all gauges the relations between the true fields and the potentials are the same, viz.

Equation (16)

When retardation effects are neglected (c → ) both gauges coincide and the only potential that remains is the instantaneous scalar potential.

4. Background material for 2D systems

4.1. Fourier transforms of a particular function

The different versions of the Fourier transform of the function F(R) = 1/R appear many times in the derivations so we compile the results in this subsection. This function is the Coulomb potential between two electrons, v(R) = e2/R, stripped of the charges.

The 3D Fourier transform is

Equation (17)

the 2D Fourier transform performed over a plane the distance z from the xy-plane is

Equation (18)

the 2D Fourier transform performed over the xy-plane is

Equation (19)

These relations appear in expressions involving the scalar potential in Coulomb gauge and when retardation effects are negligible. When retardation is important the function has a temporal dependence, F(R, t) = δ(t − R/c)/R. Then the following Fourier transforms are important: the 3D Fourier transform,

Equation (20)

the 2D Fourier transform performed over a plane the distance z from the xy-plane,

Equation (21)

the 2D Fourier transform performed over the xy-plane,

Equation (22)

where in the last two equations

Equation (23)

We will need one further relation, viz.

Equation (24)

4.2. Method of images

When dealing with a 2D layer above a substrate it is useful to be able to find the image charge- and current-densities produced in the substrate. We start from a given time dependent charge- and current-density in the 2D layer, calculate the potentials and, from these, the E- and B-fields. We make an ansatz for the image densities and determine these by using the standard boundary conditions. Here, we use the Lorentz gauge in which the scalar and vector potentials are given in (15). Both are convolution integrals of the type

Equation (25)

which has the nice property that F (q, ω) = G (q, ω)H (q, ω). In the case of the scalar potential these functions are G(R, t) = ρ(R, t) and H(R, t) = δ(t − R/c)/R, respectively. For the vector potential G(R, t) = J(R, t) and H(R, t) = δ(t - R/c)/R, respectively. Thus using (20) we find the Fourier transforms of the potentials are

Equation (26)

The electric field in terms of the potentials was given in (16). This leads to the following relation between the Fourier transforms

Equation (27)

Thus,

Equation (28)

Now, in our system the charge- and current-densities are surface densities,

Equation (29)

and

Equation (30)

which means that

Equation (31)

Thus, we have for the electric field from surface charge- and current-densities, confined to the xy-plane

Equation (32)

in vacuum. If the source densities are embedded in a medium with dielectric function ${{\tilde \varepsilon}_m}\left( \omega \right)$ the electric field is

Equation (33)

The magnetic induction is

Equation (34)

in vacuum and

Equation (35)

in a medium.

We may use the equation of continuity to eliminate ρs(k, ω) in favor of K(k, ω). The equation of continuity reads

Equation (36)

The relation for the Fourier transforms becomes

Equation (37)

In the present system with surface densities we have from (31)

Equation (38)

and since K is in the plane

Equation (39)

The scalar product picks out the longitudinal part of the surface current density. Thus we have

Equation (40)

We use this relation to eliminate the charge density from the field relations. From now on we assume that the 2D layer is embedded in a medium. Let us now choose the x-axis to point along k. A general vector quantity has a component normal to the planar interfaces. We attach a subscript n to this component. The in-plane part of the vector has a longitudinal and a transverse part. The longitudinal is parallel and the transverse perpendicular to k. We attach the subscripts ∥ and ⊥, respectively, to these components. The current density has no normal component. To summarize we have

Equation (41)

We then have in the source plane

Equation (42)

where ${\gamma _m}\left( {k,\omega} \right) = \sqrt {1 - {{\tilde \varepsilon}_m}\left( \omega \right){{\left( {\omega /ck} \right)}^2}}$ . In a plane parallel to the source plane (xy-plane) we have

Equation (43)

where we have used (40) in the third relation.

Now, we have all we need to determine the image densities. In the spirit of image theory we assume that the fields outside the substrate can be reproduced by the fields from the actual source densities in the 2D layer (the distance d from the interface) plus the fields from the mirror densities (the distance d from the interface on the opposite side), both sources embedded in the medium with dielectric function ${{{\tilde \varepsilon}_m}\left( \omega \right)}$ . The field inside the substrate we assume can be reproduced by another mirror density at the position of the 2D layer embedded in the medium with dielectric function ${{{\tilde \varepsilon}_s}\left( \omega \right)}$ . We mark the mirror densities of the first kind by a prime and those of the second kind by a double prime.

We start with the longitudinal current densities. We first use the boundary condition that the in-plane electric field is continuous across the medium-substrate interface. This gives

Equation (44)

Then we use the condition that the normal component of the displacement field ${{\bf \tilde D}} = \tilde \varepsilon {{\bf E}}$ is continuous across the interface. This gives

Equation (45)

Combining these two equations results in the mirror densities

Equation (46)

Now, we continue with the transverse current densities. We first use the boundary condition that the in-plane TE field is continuous across the interface. This gives

Equation (47)

Then we use the condition that the in-plane longitudinal component of the B-field is continuous across the interface. This gives

Equation (48)

These two equations result in the mirror densities

Equation (49)

In later sections we will need the single primed mirror densities in (46) and (49).

4.3. Dielectric screening in a 2D system

The relation between the dielectric function, $\tilde {\varepsilon}$ and dynamical conductivity, $\tilde{\sigma}$ , is different in 2D and 3D. The relations are

Equation (50)

where $\tilde {\alpha}$ is the polarizability. The relations between the induced current densities and the electric field are

Equation (51)

Note that both bound and conduction carriers contribute to the polarizabilities and conductivities. It is sometimes convenient to introduce another correlation function, χ, the polarization bubble (in the language of Feynman diagrams), or lowest order contribution to the density–density correlation function,

Equation (52)

We now give the analytical expression for the 2D dielectric function of graphene and of a 2D electron gas at zero temperature. Let us first begin with an undoped graphene sheet. In a general point, z, in the complex frequency plane, away from the real axis the density–density correlation function is [25]

Equation (53)

where v is the carrier velocity which is a constant in graphene (E = ±vk) and g represents the degeneracy parameter with the value of 4 (a factor of 2 for spin and a factor of 2 for the cone degeneracy.) In our numerical calculations [610] we use the value [26] 8.737 23 × 105 m s−1 for v. In the present work we do not give any numerical results.

When the graphene sheet is doped the dielectric function becomes much more complicated. However it has been derived by several groups [7, 24, 26, 27]. In the two next equations we use dimension-less variables: x = k/2kF; y = ω/2EF; $\tilde z = \hbar z/2{E_{\rm F}}$ . The density–density correlation function in a general point in the complex frequency plane, z, away from the real axis is [7]

Equation (54)

where ${D_0} = \sqrt {gn/\pi {\hbar ^2}{v^2}}$ is the density of states at the fermi level and n is the doping concentration. The same result holds for excess of electrons and excess of holes. In the calculations one needs the function at the imaginary frequency axis. We have derived a very useful analytical expression valid along the imaginary axis; an expression in terms of real valued functions of real valued variables [7, 8]:

Equation (55)

where the arcus tangens function is taken from the branch where 0 ⩽ atan < π. The density–density correlation function on the imaginary frequency axis has been derived before in a compact and inexplicit form (see [28] and references therein). Here we have chosen to express it in an explicit form in terms of real valued functions of real valued variables.

The 2D polarizability of an electron gas in RPA (random phase approximation) is given by [29]

Equation (56)

where

Equation (57)

Now we have presented all the background material needed for the actual derivations.

5. Non-retarded treatment, van der Waals interactions

In the non-retarded treatment many terms in MEs drop out and the derivations become much simpler. There are no photons, only Coulomb interactions. The dispersion interactions treated here are called van der Waals interactions or sometimes non-retarded Casimir interactions. We begin in section 5.1 with two parallel 2D sheets and use different methods to derive the results. In section 5.2 we treat a 2D sheet above a substrate.

5.1. Two 2D films in non-retarded treatment

5.1.1. Many-body approach in non-retarded treatment.

For a system of parallel 2D sheets the derivation in the language of many-body theory becomes very simple, especially when retardation effects can be neglected. In the present two-sheet system the interaction energy is only the inter-sheet correlation energy [30]. Using diagrammatic perturbation theory the Feynman diagrams representing the correlation energy are given in figure 1. To get the inter-sheet contribution we can either subtract the intra-sheet part or subtract the total result when the separation between the sheets goes to infinity (at that limit only the intra-sheet contribution remains). To save space we assume that the two sheets are identical. It is straightforward to extend the treatment to different sheets. Each ellipse represents a polarization bubble, χ(k, ω) and the number indicates in which sheet the process occurs. Each of the four diagrams represents an infinite series of diagrams. Since the sheets are identical the diagrams in the second row give contributions identical to the diagrams in the first row. We can use the first row diagrams and multiply the result with a factor of two. The interaction energy per unit area can be written as [10]

Equation (58)

where λ is the coupling constant and the factor of two has been inserted.

Figure 1.

Figure 1. Feynman diagrams (from [30]) for the correlation energy in the two 2D sheet system. The ellipses represent polarization bubbles and the dashed lines the interactions indicated in the figure. The numbers 1 and 2 refer to which sheet the carrier belongs to. See [30] for details.

Standard image High-resolution image

Now v11 = v2D(k, ω) and v12 = exp(−kd)v2D(k, ω), respectively. Note that with our assumptions v22 = v11 and v21 = v12. The distance between the sheets is d and v2D is the 2D Coulomb interaction, v2D(r) = e2/r. The exponential factor in the second potential is the result from taking the 2D Fourier transform of the Coulomb potential in a plane the distance d from the center of the potential, see (18). The interaction lines with double bars represent a series of terms, with zero, one, two ... number of polarization bubbles. This can be expressed as

Equation (59)

where we have closed the two infinite series. This system of equations can be solved and the result is

Equation (60)

where we have used the relation $\tilde \alpha \left( {{{\bf k}},\omega} \right) =- {v^{{{\rm 2D}}}}\left( k \right)\chi \left( {{{\bf k}},\omega} \right)$ in accordance with (52). Interpreting the Feynman diagrams in figure 1 one finds that the square brackets in (58) are $\left[ {{{\tilde v}_{11}}{v_{11}}{\chi ^2} + {{\tilde v}_{12}}{v_{12}}{\chi ^2} - \left( {{v^{2D}}/\varepsilon} \right){v^{2D}}{\chi ^2}} \right]$ , where the last term comes from the subtraction of the intra-layer correlation energy. Each factor of e2 appearing implicitly in the expression should be multiplied by the coupling constant. Performing the integration over coupling constant gives

Equation (61)

Here for simplicity we assumed that the two 2D layers were identical. If they are not the result is

Equation (62)

This is the result from diagrammatic perturbation theory within the RPA. In the next section we derive the same thing using normal modes.

5.1.2. Coupled induced sources in non-retarded treatment.

We will now go through one possible way to find the normal modes and their zero-point energies. Let us assume that we have an induced charge distribution, ρ1(k, ω), in sheet 1. This gives rise to the potential v(k, ω) = v2D(k)ρ1(k, ω) and exp(−kd)v2D(k)ρ1(k, ω) in sheets 1 and 2, respectively. The resulting potential in sheet 2 after screening by the carriers is $\exp \left( {- kd} \right){v^{2D}}\left( k \right){\rho _1}\left( {{{\bf k}},\omega} \right)/\left[ {1 + \tilde {\alpha}_2 \left( {{{\bf k}},\omega} \right)} \right]$ , which gives rise to an induced charge distribution in sheet 2,

Equation (63)

In complete analogy, this charge distribution in sheet 2 gives rise to a charge distribution in sheet 1,

Equation (64)

To find the condition for self-sustained fields, or normal modes, we let this induced charge density in sheet 1 be the charge density we started from. This leads to

Equation (65)

The left hand side of this equation (fk(ω)) is exactly the argument of the logarithm in (62). Inserting this function into (3) we obtain exactly the same energy as in (62). That result was derived using many-body theory and in that derivation there is no discussion of zero-point energies at all.

5.1.3. Maxwell's equations and boundary conditions in non-retarded treatment.

Here we solve MEs in the three regions separated by the two 2D layers. We treat the induced charge and current densities in the 2D layers as external to our system. Since we have vacuum in between the layers our system is vacuum and the D and H fields are equal to the E and B fields, respectively. It is enough to study the E-field in the non-retarded treatment. We follow the procedure in [3] and make the following ansatz

Equation (66)

The two MEs involving the E-field, ∇ · E = 0 and ∇ × E = 0, both give

Equation (67)

The boundary conditions are that the x-components of the fields are continuous across the two interfaces and there is a jump in the z-component equal to 4πρs. Thus we have

Equation (68)

or on matrix form

Equation (69)

and

Equation (70)

Now, we eliminate the z-components in favor of the x-components using (67) and use ${\rho _s} = kK/\omega = k\tilde {\sigma} {{{E}}_x}/\omega$ and assume for simplicity that the two layers are identical. We find

Equation (71)

or on matrix form

Equation (72)

Combining (69) and (72) we find that the condition for having normal modes is that the determinant of the matrix

Equation (73)

vanishes. The mode condition function becomes

Equation (74)

We have taken as a reference system the system when the two sheets are at infinite separation. This means that we have divided the mode condition function with the corresponding result when d is infinite. If the two layers are different we instead find

Equation (75)

which agrees with the results of the previous sections.

5.1.4. Coupled potentials in non-retarded treatment.

Another way to find the normal modes is to find the frequency where the potential between carriers in the same layer or in different layers diverges. The condition for modes is that the denominators in (60) are equal to zero. This leads to the same mode condition function as in (65) and (75) when we use as a reference system the system when the 2D layers are at infinite separation.

5.2. One 2D film above a substrate in non-retarded treatment

5.2.1. Coupled induced sources in non-retarded treatment.

If we have a 2D layer (like a graphene sheet) with the distance d above a substrate the procedure is very similar to that in section 5.1.2. We start with an induced mirror charge density, ρ1(k, ω), in the substrate. The induced charge density in the 2D layer is given by the expression in (63) except that now the distance between the mirror charge and the 2D layer is 2d instead of d. The mirror charges can be obtained from (40) and (46) and agree with the traditional results for the mirror charges at a planar vacuum dielectric interface [31]. Equation (64) is then replaced by

Equation (76)

and the condition for normal modes becomes

Equation (77)

and hence the mode condition function is

Equation (78)

5.2.2. Maxwell's equations and boundary conditions in non-retarded treatment.

In section 5.1.3 we had a geometry with two interfaces, one at z = −d and one at z = 0. There there was a 2D layer at both interfaces. We use the same geometry here but there is a 2D layer only at z = −d. At z = 0 we place the vacuum-substrate interface. The derivation and equations are equal up until (70). When the media are no longer vacuum there is a jump in ${{{{\tilde D}}}_z}$ at the 2D layer instead of in Ez. Note that here there is only a 2D layer at the left interface so the normal component of the $\tilde{D}$ -field is continuous at the substrate surface. The new equation is

Equation (79)

Now we eliminate the z-components in favor of the x-components and use ${\rho _s} = kK/\omega = k\tilde {\sigma} {{{E}}_x}/\omega$ so that

Equation (80)

or on matrix form

Equation (81)

Combining (69) and (81) we find that the condition for having normal modes is that the determinant of the matrix

Equation (82)

is zero. Thus, the mode condition function becomes

Equation (83)

which agrees with the results of the previous sections.

6. Retarded treatment, Casimir interactions

Here we take retardation effects fully into account. The formalism becomes much more involved. The dispersion interactions treated here are called Casimir interactions or sometimes retarded van der Waals interactions. We begin in section 6.1 with two parallel 2D sheets and use different methods to derive the results. In section 6.2 we treat a 2D sheet above a substrate.

6.1. Two 2D films in fully retarded treatment

6.1.1. Many-body approach in fully retarded treatment.

We do not perform a detailed derivation here, we just describe how the derivation is performed. We follow the derivation in [30] and refer to that work for more details. We have chosen to work in Coulomb gauge. In this gauge one part of the interaction between the electrons is in the form of the instantaneous, longitudinal Coulomb interaction (scalar potential) and the other part is via transverse photons (vector potential). The Hamiltonian for the system is

Equation (84)

where the first term is the kinetic energy that contains the interactions via the vector potential A. The second term represents the scalar potential interaction and the last term is the free photon Hamiltonian. The non-retarded results are obtained by letting the speed of light tend to infinity. From the Hamiltonian we see that this corresponds to neglecting the vector-potential interaction completely. Now, in the specific system we consider here, the electrons are only free to move in a plane. In a plane only one longitudinal and one TE field can exist. It turns out that the longitudinal electric field and the p-polarized photons combine into one field that is longitudinal in the plane. The induced charge and current densities from such a field will also produce a field of that kind. The s-polarized photons produce an electric field that is transverse in the plane. The current induced by such a field produces a field of the same type.

The two types of interaction involving longitudinal and TE fields in the planes defined by the two sheets will not mix and the energy is given by a sum of two sets of diagrams of the type given in figure 1, one set for the longitudinal interaction and one containing only transverse interactions.

Let us now discuss the longitudinal interaction. The coupling to the longitudinal field and p-polarized photons is given by the scalar potential and the p · A terms in the Hamiltonian. The resulting interaction energy per unit area is

Equation (85)

Compared with the non-retarded result for the correlation energy we gave earlier in (61), we find that the polarizability has attained a factor γ(0)(k, iω) in front, which is in agreement with the polarizability obtained by Stern [32]. Furthermore, the factor appearing in the interlayer interactions is modified with the same factor in the exponent.

The response to the s-polarized fields is a transverse current. It is described by the transverse conductivity, which is dominated by the contribution originating from the A2-term in the Hamiltonian. Also the p·A terms in the Hamiltonian give a contribution to the transverse conductivity, in the form of the current–current correlation function. The energy per unit area from the s-photon interaction is

Equation (86)

6.1.2. Coupled induced sources in fully retarded treatment.

We could here proceed in a way analogous to how we derived the modes in section 5.1.2 by starting with a current density (longitudinal or transverse) in layer 1. This current density gives rise to an electric field in layer 2; this field induces a current density in layer 2; this current density in layer 2 gives rise to an electric field in layer 1 which causes an induced current density in layer 1. If we let this current density be the one we started with we have closed the loop and produced self-sustained fields, or normal modes. When the current densities and electric fields are longitudinal, the modes are of TM type and the TE modes are the results when the current densities and fields are transverse.

To vary the derivations we instead proceed in a slightly different way. We have two layers, 1 and 2. The electric field in layer 1 has contributions from both layers and the induced current density in layer 1 is the conductivity of layer 1 times the electric field in layer 1. The same is true for layer 2. We begin with the TM modes. The longitudinal field in layer 1 is from (42) and (43)

Equation (87)

and the induced current density is

Equation (88)

In analogy we get the result for the induced current density in layer 2 by interchanging the superscripts. Thus we have

Equation (89)

and on matrix form

Equation (90)

where we have omitted all function arguments (k, ω). This system of equations has a trivial solution, in which both current densities are zero. For a non-trivial solution the determinant of the matrix is equal to zero. Thus we get the condition for modes by letting the determinant of the matrix be equal to zero. The mode condition function becomes

Equation (91)

where we have let the reference system be when the layers are at infinite separation. For the TE modes we need the transverse fields. They are from (42) and (43)

Equation (92)

and the induced current densities

Equation (93)

or on matrix form

Equation (94)

The mode condition function for TE modes is

Equation (95)

where we have chosen as a reference system the geometry when the layers are infinitely apart.

6.1.3. Maxwell's equations and boundary conditions in fully retarded treatment.

The normal modes can be found from self-sustained charge and current densities, from self-sustained potentials or from self-sustained fields. Here, we will find self-sustained fields. To do that one solves the MEs in all regions of the geometry and makes use of the standard boundary conditions at all interfaces between the regions. For the present task we need a geometry consisting of three regions and two interfaces, 1|2|3. This geometry gives the mode condition function

Equation (96)

where rij is the amplitude reflection coefficient for a wave impinging on the interface between medium i and j from the i side, d is the thickness of region 2 and ${\gamma _i} = \sqrt {1 - {{\tilde \varepsilon}_i}\left( \omega \right){{\left( {\omega /ck} \right)}^2}}$ . The function ${{{\tilde \varepsilon}_i}\left( \omega \right)}$ is the dielectric function of medium i. At an interface where there is no 2D sheet the TM and TE amplitude reflection coefficients are [3]

Equation (97)

and

Equation (98)

respectively. Note that rji = −rij holds for both mode types.

The amplitude reflection coefficient gets modified if there is a 2D layer at the interface. We treat the 2D layer at the interface as external to our system. The fields will then induce external surface charge and current densities at the interface. Two of the boundary conditions are enough to obtain the modified Fresnel coefficients. The other two are redundant. We chose the following two,

Equation (99)

where $\tilde {\sigma}$ is the conductivity of the 2D sheet. The modified amplitude reflection coefficient for a TM mode is [8]

Equation (100)

where the polarizability of the 2D sheet is obtained from the dynamical conductivity according to (50). For TM modes the tangential component of the electric field, which will induce the external current, is parallel to k, so the longitudinal 2D dielectric function of the layer enters. The bound charges in the 2D sheet also contribute to the dynamical conductivity and the polarizability.

The modified amplitude reflection coefficient for a TE mode is [8]

Equation (101)

where the polarizability of the 2D sheet is obtained from the dynamical conductivity according to (50). For a TE wave the electric field is perpendicular to k, so the transverse 2D dielectric function of the layer enters. The bound charges in the 2D sheet also contribute to the dynamical conductivity and the polarizability.

Now, when we have geometrical structures with 2D layers at some of the interfaces we just substitute the new reflection coefficients at the proper interfaces. Note that r21 = −r12 no longer applies so one has to be careful when making the substitutions. To arrive at the starting expression one might have used this relation. For the three regions system 1|2|3 we find the most general mode condition functions as

Equation (102)

where we assumed there are 2D sheets at both interfaces. These sheets may be different so we have put the subscripts L and R on the polarizability of the left and right sheet, respectively.

We obtain two freestanding 2D sheets from our geometry by letting all three media be vacuum. Then the mode condition functions in (102) reduce to

Equation (103)

6.2. One 2D film above a substrate in fully retarded treatment

6.2.1. Coupled induced sources in fully retarded treatment.

When a 2D layer is placed next to a substrate there will be induced current densities at the interface. The effect of these densities can be obtained in a much simpler way using mirror images. What these mirror images look like was derived in the section 4.2 for a more complex system where the 2D layer is embedded in a medium. For our system the mirror-image current-densities are

Equation (104)

when K(k, ω) and K(k, ω) are the actual longitudinal and transverse, respectively, current densities in the 2D layer. These mirror images are at a distance d from the interface and inside the substrate. We have suppressed the argument (k, ω) in all γ-functions. The coupled longitudinal current densities give rise to TM modes and the transverse to TE modes.

We start with the TM modes. The resulting longitudinal electric field from a longitudinal current density was derived in (42) and (43). The field in the 2D layer is

Equation (105)

We close the loop by noting that the current density in the 2D layer is the conductivity times the electric field. Thus,

Equation (106)

where we have used (104) in the last step to eliminate the mirror density. We may now use the first and last parts of this equation to find

Equation (107)

This has the trivial solution that the current density is zero. There are non-trivial solutions, the normal modes, if what is inside the curly brackets is zero. This leads to the mode condition function for TM modes,

Equation (108)

where we have chosen as a reference system the system when the 2D layer is at infinite distance from the substrate.

Now, we proceed with the TE modes. The resulting TE field from a transverse current density was derived in (42) and (43). The field in the 2D layer is

Equation (109)

We close the loop by noting that the current density in the 2D layer is the conductivity times the electric field. Thus,

Equation (110)

where we have used (104) in the last step to eliminate the mirror density. We may now use the first and last parts of this equation to find

Equation (111)

This has the trivial solution that the current density is zero. There are non-trivial solutions, the normal modes, if what is inside the curly brackets is zero. This leads to the mode condition function for TE modes,

Equation (112)

where we have chosen as reference system the system when the 2D layer is at infinite distance from the substrate.

6.2.2. Maxwell's equations and boundary conditions in fully retarded treatment.

We obtain a freestanding 2D sheet above a substrate from our geometry in section 6.1.3 by letting media 1 and 2 be vacuum and medium 3 be the substrate. Then the mode condition functions in (102) reduce to

Equation (113)

If retardation effects are neglected there are no TE modes and the mode condition function for TM modes is reduced further to

Equation (114)

which is fully in line with what we obtained starting from the non-retarded theory.

7. Summary and conclusions

We have derived the expressions for the dispersion interaction in systems containing 2D sheets in several different complementary ways, all with a condensed matter physics perspective. The geometries considered were two free standing 2D sheets and one 2D sheet above a substrate. The derivations were either performed using many-body theory with diagrammatic perturbation expansions based on Feynman diagrams or were based on electromagnetic normal modes. The normal modes were generated in several different ways. The zero-point energy is a very important concept in the normal mode derivation of the dispersion interaction. We note that the same results are obtained using many-body theory without invoking zero-point energy. Both approaches, although they appear to be quite different, derive the same thing and give identical results. We refer the reader to chapter 3 of [3] where it is shown that the exchange and correlation energy of a metal is just the change of the zero-point energy of the longitudinal normal modes of the system when the electrons are brought together at the formation of the solid.

Our derivations here were all made with a condensed matter physics perspective and were based on the density–density or current–current correlation functions. Numerical works with this perspective can be found in the literature in [616, 30, 33, 34]. In [6] it was shown how a strictly 2D film could be represented by a film of finite thickness with the proper choice of 3D dielectric function. In [33] it was noted that the interaction between two thin metal films approached the interaction between two 2D metal films when the film thickness goes towards zero. In [7] the non-retarded results for two freestanding graphene sheets and for a graphene sheet above a gold substrate at zero temperature were determined. Undoped and doped graphene sheets were treated. In [8] the results from a fully retarded calculation were presented for zero and room temperature, for undoped and doped graphene sheets. In [9] it was demonstrated that when a bilayer of graphene on a silica substrate was exposed to hydrogen, the gas atoms had a tendency to enter in between the graphene layers and separate them. In [10] the non-retarded interaction between doped and undoped graphene sheets and between metal films were studied at zero temperature. It was noted that the effect of spatial dispersion found in the metal film geometry [30, 33] was absent for undoped graphene but present in the contribution from the dopant carriers in the doped graphene case. Furthermore, the fractional power law discovered [30, 35] in the metal film geometry was also obtained for doped graphene sheets in the separation region where the doping contribution dominated. In [11] one studied van der Waals interactions between graphene layers at finite temperature and found that retardation effects were negligible and that the thermal effects set in much earlier than expected, already at distances of tens of nanometers at room temperature. In [12] the focus was on the change in force between a silica membrane and a metal plate when the membrane was covered by a graphene sheet. In [13] the force between undoped graphene sheets at zero temperature was studied. They used a model where the conductivity is a constant, σ(k, ω) = σ0 = e2/4. This is the limiting form for vk < < ω,

Equation (115)

One has to be careful with this approximation, since it is not valid in the finite temperature n = 0 term, which dominates at large separations: in [14] and [15] how a graphene sheet or a 2D electron gas could be used to modify the reflection of atoms was studied. In [16] one presented an alternative general way to derive the electromagnetic normal modes and Casimir effects in layered structures, also including 2D films. Special emphasis was put on planar, spherical and cylindrical structures. In [34] one derived the thermal Casimir effect for two parallel 2D metal films.

Other derivations in the literature were made with a quantum field theory perspective and were based on the polarization tensor defined in a (2 + 1)-dimensional space time [1723]. In principle this approach produces the same results as the other. It was shown in [24] that the two approaches are indeed equivalent. Deviations may occur due to different approximations used for the dielectric functions of the 2D systems. In [17, 18] they let graphene have a small but finite band gap. Thus, the results are different. In [1923] a zero band gap was assumed and the results were as expected. In [23] they calculated the gradient of the Casimir force between an Au-coated sphere and a graphene sheet deposited on a SiO2 film covering a Si plate, which is the configuration of a recent experiment [36] performed by means of a dynamic atomic force microscope. The theoretical results were found to be in very good agreement with the experimental data.

One characteristic feature for a system of two parallel graphene sheets is that the effect of finite temperature appears at a much shorter distance than in other planar systems made by conventional materials [11, 12, 1822].

In the reference list we have up till now only included references related to the specific geometries treated in this work. If the reader is interested in the Casimir Effect in general we suggest the latest review book [37].

To summarize, the Casimir effects can be approached from several fields of physics and the derivations can appear to be very different but the results are, at the end, the same. We have here shown how the Casimir effects in planar systems containing 2D layers can be derived in a variety of ways, all with the condensed matter physics perspective.

Please wait… references are loading.
10.1088/0953-8984/27/21/214017