Unifying Theory for Casimir Forces: Bulk and Surface Formulations

The principles of the electromagnetic fluctuation-induced phenomena such as Casimir forces are well understood. However, recent experimental advances require universal and efficient methods to compute these forces. While several approaches have been proposed in the literature, their connection is often not entirely clear, and some of them have been introduced as purely numerical techniques. Here we present a unifying approach for the Casimir force and free energy that builds on both the Maxwell stress tensor and path integral quantization. The result is presented in terms of either bulk or surface operators that describe corresponding current fluctuations. Our surface approach yields a novel formula for the Casimir free energy. The path integral is presented both within a Lagrange and Hamiltonian formulation yielding different surface operators and expressions for the free energy that are equivalent. We compare our approaches to previously developed numerical methods and the scattering approach. The practical application of our methods is exemplified by the derivation of the Lifshitz formula.


Introduction
The interaction induced by quantum and thermal fluctuation of the electromagnetic 15 field is an everyday phenomenon that acts between all neutral objects, both on atomic and 16 macroscopic scales [1][2][3][4][5]. For the Casimir interaction between macroscopic bodies, the 17 last two decades have witnessed unparalleled progress in experimental observations and 18 the development of novel theoretical approaches [6,7]. In most of the recent theoretical 19 approaches, the computation of Casimir forces between multiple objects of different shapes 20 and material composition has been achieved by the use of scattering methods or the 21 so-called TGTG formula [8][9][10][11][12][13][14][15][16]. These approaches have the advantage of relatively low 22 numerical effort; they are rapidly converging and can achieve in principle any desired 23 precision [17][18][19][20][21][22]. Another merit of these methods is the exclusion of UV divergencies 24 by performing the subtraction analytically before any numerical computation. Other 25 efficient approaches that have been developed before the scattering approaches include 26 path integral quantizations where the boundary conditions at the surfaces are implemented 27 by delta functions [23]. These approaches are limited to scalar fields with Dirichlet or 28 Neumann boundary conditions [24], or the electromagnetic field with perfectly conducting 29 boundary conditions [25], with the exception of a similar approach for dielectric boundaries 30 [26]. However, such analytical (and semi-analytical) methods have been restricted to 31 symmetric and simple shapes, like spheres, cylinders or ellipsoids [27][28][29][30][31]. Geometries 32 Version August 17, 2021 submitted to Universe https://www.mdpi.com/journal/universe introduced in the literature as a method for a purely numerical computation of Casimir 48 interactions [32]. There are two different methods to implement the idea of computing 49 Casimir forces from fluctuating currents. One can either integrate the Maxwell stress tensor 50 over a closed surface enclosing the body, directly yielding the Casimir force, or integrate 51 over all electromagnetic gauge field fluctuations in a path integral, yielding the Casimir 52 free energy. We shall consider both approaches here. 53 Compared to scattering theory-based approaches, the surface formulation has the 54 advantage that it does not require the use of eigenfunctions of the vector wave equation 55 that are specific to the shape of the bodies. Hence, our approach is applicable to general 56 geometries and shapes, including interpenetrating structures. In fact, the power of the 57 surface approach has been demonstrated by numerical implementations in Reference [32], 58 where it was used to compute the Casimir force in complicated geometries.

59
In this paper, we present both the Maxwell stress tensor and path integral-based 60 approaches for the Casimir force and free energy in terms of bulk or surface operators. Our the Casimir free energy. 67 We also compare the approaches presented here to methods existing in the literature. 68 For the special case of bodies that can be separated by non-overlapping enclosing surfaces, 69 along which one of the coordinates in which the wave equation is separable is constant, 70 our approach is shown to be equivalent to the scattering approach. Our approaches also 71 show the general equivalence of the use of the Maxwell stress tensor in combination with 72 the fluctuation-dissipation theorem on one side and the path integral representation of the 73 Casimir force on the other side. As the most simple application of our approaches, we 74 re-derive the Lifshitz formula for the Casimir free energy of two dielectric slabs. Other 75 analytical applications of our approach will be presented elsewhere. 76 The geometries and shapes to which our approaches can be applied are shown in

Stress-Tensor Approach
over any closed surface S r drawn in the vacuum, which surrounds that body (but excludes all other bodies): wheren is the unit normal oriented outside S r , and the angular brackets denote the expectation value taken with respect to quantum and thermal fluctuations. For a system in thermal equilibrium at temperature T, the (equal-time) expectation values of the products of field components (at points x and x in the vacuum region) are provided by the fluctuation-dissipation theorem [34,35]: where ξ n = 2πnk B T/h are the Matsubara imaginary frequencies, and the prime in the summations mean that the n = 0 term is taken with a weight of one half. When the r.h.s of the above equations are plugged into Equation (1), one obtains for T ij (x) a formally divergent expression, since the Green functions G (αβ) ij (x, x ; i ξ n ) are singular in the coincidence limit x = x . This divergence can however be easily disposed of by noticing that the Green's functions admit the decomposition: where G ij (x, x ; i ξ n ) describes the effect of scattering of electromagnetic fields by the bodies. In the coincidence limit, only G ij (x, x ; i ξ n ) attains a finite limit. When the decomposition in Equation (4) is used to evaluate the r.h.s. of Equation (1), one finds that the expectation value of the stress tensor is decomposed in a way analogous to Equation (4): where T (0) ij (x) is the divergent expectation value of the stress tensor in empty space, and Since the divergent contribution T ij (x) is independent of the presence of the bodies, one can just neglect it and then one obtains the following finite expression for the Casimir force on body r due to the presence of the other bodies, The further development of the theory starts from the observation that the dyadic ij (x, x ) (for brevity, from now on we shall not display the dependence of the Green's functions on the Matsubara frequencies ξ n ) can be expressed in two distinct possible forms.
The first representation is general, since it is valid for arbitrary constitutive equations of the magneto-dielectric materials constituting the bodies, which can possibly be nonhomogeneous, anisotropic, and non-local. For all points x and x , it expresses Γ (αβ) ij (x, x ) in the form of an integral of the T-operatorT (for its definition, see Appendix B.1) over the volume V occupied by all bodies: The above formula has a simple intuitive interpretation, if one recalls that according to its definition, the T-operator provides the polarization induced in the volume of the bodies when they are immersed in the electromagnetic field generated by a certain distribution of external sources.
The second representation is less general than Equation (8) because it applies only to magneto-dielectric bodies that are (piecewise) homogeneous and isotropic 1 . It expresses Γ (αβ) ij (x, x ) in the form of an integral of the surface operatorM −1 defined in Equation (A44), over the union Σ = Σ r of the surfaces Σ r of the bodies. For two points x and x both lying in the vacuum region outside the bodies 2 , the surface representation of Γ where F r (y) = 0 is the equation of the surface Σ r . This representation also has a simple intuitive meaning, if one considers that −M −1 (see Appendix B.2 for details) is defined as the operator that provides the fictitious surface polarizations that radiate outside the bodies the same scattered electromagnetic field as the one radiated by the physically induced volumic polarization, in response to an external field. The derivations of Equations (8) and (9) are presented in Appendix B. It is apparent that both representations have the same mathematical structure, consisting of a two-sided convolution of a certain kernel K The only difference between the two representations consists in the expression of K, which in the case of Equation (8) is the three-dimensional kernel T (ρσ) kl (y, y ) supported in the volume V occupied by the bodies: while in Equation (9) K is the two-dimensional kernel − M −1 (ρσ) kl (y, y ) supported on the union Σ of their surfaces: In both cases, the above equation can be concisely written using the operator notation described in Appendix A: In Appendix C, we prove that the structure of the representation of Γ (αβ) ij (x, x ) given in Equation (10), allows to re-express the Casimir force Equation (7) in the following remarkably simple form: A crucial role in the derivation of Equation (14) is played by the fact that the kernel K satisfies a set of reciprocity relations analogous to those satisfied by the Green's functions: where s(E) = 0 and s(H) = 1. It is possible to verify that the reciprocity relations satisfied by G This implies that in Equation (14) the y integral is in fact restricted to V − V r , in accord with one's intuition that the force on body r is due to the interaction with the other bodies. It is possible to present Equation (14) in a more compact and symmetric form, by defining the derivative ∂/∂x r of any kernel A(y, y ), with respect to rigid translations of the r-th body: where ψ r (y) is the characteristic function of V r : ψ r (y) = 1 if y ∈ V r , ψ r (y) = 0 if y / ∈ V r . Using ∂/∂x r , we can rewrite Equation (14) as The expression on the r.h.s. of the above formula can be compactly expressed using the operator notation and the trace operation described in Appendix A: Depending on whether we use for the kernelK the T-operator of Equation (8) or rather the surface operator −M of Equation (9), Equation (19) provides us with two distinct but formally similar representations of the Casimir force, which is expressed either as an integral over the volume V occupied by the bodies or as an integral over their surfaces Σ. One feature of Equation (19) is worth stressing. Since the force is expressed as a trace, Equation (19) can be evaluated in an arbitrary basis, leaving one with complete freedom in the choice of the most convenient basis in a concrete situation. A representation of the Casimir force in the form of a volume integral equivalent to Equation (19) was derived in [11], while the surface-integral representation was obtained in [32]. Equation (19) can be computed numerically for any shapes and dispositions of the bodies, by using discrete meshes covering the bodies. An efficient numerical scheme based on surface-elements methods is described in [32], where it was used to compute the Casimir force in complex geometries, not amenable to analytical techniques.

Casimir Free Energy
In this section, we compute the Casimir free energy F of the system of bodies, starting from the force formula Equation (19). We shall see that the T-operator and the surface-operator approaches lead to two distinct but equivalent representations of the Casimir energy.

T-Operator Approach
Plugging into Equation (19) the expression of the T-operator Equation (A36), we find that the Casimir force can be expressed in the form: where in the last passage, we made use of the fact that the polarization operatorχ defined in Equation (A24) is invariant under a rigid displacement of the body. The r.h.s. of the above equation can be formally expressed as a gradient: where F bare is the bare free energy: Unfortunately, F bare is formally divergent. To obtain the finite Casimir free energy F , one has to subtract from F bare the divergent self-energies F (r) self of the individual bodies: whereχ r =ψ rχψr is the polarizability operator of body r in isolation. In the case of a system composed by two bodies, the renormalized Casimir free energy can be recast in the following TGTG form: whereT is the T-operator of body r in isolation. To prove Equation (24), one notes that for each Matsubara mode the operator identity holds: The above identity in turn allows to prove the following chain of identities: Equating the first line with the last line, we obtain the identity: Upon summing the above identity over all Matsubara modes (with weight one half for the n = 0 term), and then multiplying it by k B T, we find: which is equivalent to Equation (24). The energy formula Equation (24) was derived in [13] using the path-integral method and in [11], using Rytov's fluctuational electrodynamics [36].

Surface Operator Approach
Now we derive the surface-operator representation of the Casimir energy. To do that, we start from the surface-operator representation of the force, which is obtained by replacingK in Equation (19) with minus the inverse of the surface operatorM defined in Equation (A44): This can also be written as: whereΠ is the tangential projection operator defined in Appendix B.2. Now, one notes the identity: which is a direct consequence of Equation (A44) since Plugging Equation (32) into Equation (31), we obtain: The r.h.s. of the above equation can be formally expressed as a gradient: whereF bare is the bare free energy: Similarly to what we found in the T-operator approach, the surface formula of the bare-energyF bare is formally divergent. The finite Casimir free energy is obtained by subtracting fromF bare the limitF (∞) bare of the bare energy when the bodies are taken infinitely apart from each other. From Equation (A44), one sees that in the limit of infinite separations, the operatorM approaches the limitM ∞ whereM Notice that the surface operatorM r is localized onto the surface Σ r of the r-th body. This implies that:M Using Equation (37), we find thatF bare is the formally divergent quantity: The additive character ofF self as representing the (infinite) the self-energy of the bodies in the surface approach. Upon subtracting Equation (40) from Equation (36), we arrive at the following formula for the Casimir energy: An easy computation shows that Substitution of the above formula into Equation (41) results in the following surface formula for the Casimir energy: In the simple case of two bodies, the above formula reduces to: The surface formulas for the Casimir energy given in Equations (44) and (45) were not known before and are presented here for the first time. Comparison of Equation (45) with Equation (24) reveals the striking similarity of the T-operator and surface-approach representations of the Casimir energy. Indeed we see that both formulas can be written in the form: whereK r is the kernel, which gives the scattering Green's function of body r in isolation:

Equivalence of the Surface-Formula with the Scattering Formula for the Casimir Energy
In the previous sections, we have shown that, both in the T-operator and in the surface approaches, the Casimir energy F of two bodies can be expressed by the general Equation (46). This formula is valid for any shape and relative dispositions of the two bodies, and in particular for two interleaved bodies (see Figure 1a). Now we show that when the two bodies can be enclosed within two non-overlapping spheres (see Figure 1b), Equation (46) is the same as the well-known scattering formula [8][9][10][11][12]: where T (r) is the scattering matrix of body r (see Equation (A76) for the definition of T (r) ), U (rs) are the translation matrices defined in Equation (A82) and tr denotes a trace over multipole indices.
To prove equivalence of Equation (46) with Equation (48), one starts from the observation that the trace operation in Equation (46) involves evaluating the Green functions ↔ G (αβ;0) (y, y ) at points y and y , one of which (call it y 1 ) belongs to body 1, while the other (call it y 2 ) belongs to body 2. For two bodies that can be separated by non-overlapping spheres, it is warranted that |y 1 − X 1 | < |y 2 − X 1 | and |y 2 − X 2 | < |y 1 − X 2 |, where X 1 and X 2 are the positions of the centers of the spheres S (1) and S (2) , respectively, and d = |X 2 − X 1 | is their distance. This condition satisfied by y 1 and y 2 ensures that it is legitimate to express ↔ G (αβ;0) (y r , y s ) (with r = s = 1, 2) by the partial-wave expansion (see Appendix D): where Φ (reg/out) plm (y r − X r ) are a basis of regular and outgoing spherical waves with origin at X r . When the above expansion is substituted into Equation (46) and the trace is evaluated, one finds that F can be recast in the form: where N is the matrix of elements: pl,p l (d) U (12) p l ;p l (d) Recalling the formula Equation (A79) for the scattering matrices T (r) of the two bodies, we see that N is the matrix: Upon substituting the above expression into the r.h.s. of Equation (50), and using cyclicity of the trace, we see that Equation (50) indeed coincides with the scattering formula Equation (48).

Path Integral Approach
As in the previous sections, we consider again N dielectric bodies occupying the volumes V r , r = 1, . . . , N, bounded by surfaces Σ r . Their electromagnetic properties are described by the dielectric functions (r) and magnetic permeability µ (r) . The bodies are embedded in a homogeneous medium occupying the outside volume of the bodies, V 0 , with dielectric function (0) and magnetic permeability µ (0) .
In the Euclidean path integral quantization of the electromagnetic field, the Casimir free energy at finite temperature T can be obtained as where the sum runs over the Matsubara momenta κ n = 2πnk B T/hc, with a weight of 1/2 for n = 0. The partition function Z is given by a path integral that we shall derive now. The partition function Z ∞ describes the configuration of infinitely separated bodies and subtracts the self-energies of the bodies from the bare free energy. In the following two sections, we shall derive both a Lagrangian and a Hamiltonian path integral expression of the partition function. In both cases, we employ a fluctuating surface current approach. A path integral approach that is based on bulk currents can be found, e.g., in Reference [10].

Lagrange Formulation
The action of the electromagnetic field coupled to bound sources P ind in the absence of free sources is in general given by In the following, we express the action in terms of the gauge field A choosing the transverse or temporal gauge with A 0 = 0. The functional integral will then run over A only. The electric field is given by E = ikA → −κA and the magnetic field by B = ∇ × A. Then the action in terms of the induced sources at fixed frequency κ is given bŷ for fluctuations A of the gauge field, and induced bulk currents P r inside the objects. The inverse of the kernel of the quadratic part of this action is given by the Green tensor G(x, x ), which is defined by For spatially constant and µ with body r, this yields the free Green's tensor which is symmetric, reflecting reciprocity. From the relation between the gauge field A and the electric field E follows the relation −κ 2 ↔ G , which allows to compare the results below to those of the stress-tensor-based derivation.
Next, we define the classical solutions A r of the vector wave equation in each region We use this definition together with the fact that A has no sources inside V r , i.e., obeys above wave equation with vanishing right-hand side, to rewrite the source terms of Now we have to consider the electric field E = −κA only on the surfaces Σ r . However, the values of the electric field E and its curl ∇ × E are those when the surface is approached from the inside, denoted by E − and (∇ × E) − . It is important to realize that in the above surface integral, A and ∇ × A multiply vectors that are tangential to the surface, and hence only the tangential components of A and ∇ × A contribute to the integral. Hence, we can use the continuity conditions of the tangential components of E and H, to write the source terms as where the first form applies to A inside the objects and the second form to A outside the objects. There is another advantage of having expressed the latter integrals in terms of the values of A and ∇ × A when the surfaces are approached from either the outside or the inside of the objects. In the region V 0 , the field A ≡ A 0 is fully determined by its values on the surfaces Σ r and the dielectric function 0 and permeability µ 0 , which are constant across V 0 . When integrating out A 0 , in fact, one computes the two-point correlation function of A + and (∇ × A) + on the surfaces Σ r , and hence the behavior of A 0 inside the regions V r with r > 0 is irrelevant. Following the same arguments for A ≡ A α inside the objects, the behavior of A r outside of region V r is irrelevant for computing the correlations of A − and (∇ × A) − on the surfaces Σ r . Hence, we can replace in the actionŜ[A, {A r }] the spatially dependent x by 0 when the coupling of A 0 to the surface fields A r is represented by the second line of Equation (60), and similarly replace x by r when the coupling of A r to the surface fields A r is represented by the first line of Equation (60). That this is justified can also be understood as follows. The field A 0 in region V 0 can be expanded in a basis of functions that obey the wave equation with 0 . The same can be done for A r in the interior of each object, i.e., A r can be expanded in a basis of functions that obey the wave equation with r in V r . For each given set of expansion coefficients in V 0 there are corresponding coefficients within each region V r that are determined by the continuity conditions at the surfaces Σ r . The functional integral over A then corresponds to integrating over consistent sets of expansion coefficients that are related by the continuity conditions. The two-point correlations of A + and (∇ × A) + on the surfaces Σ r are then fully determined by the integral over the expansion coefficients of A 0 in V 0 only, and the interior expansion coefficients play no role. Equivalently, the two-point correlations of A − and (∇ × A) − on the surfaces Σ r are then fully determined by the integral over the expansion coefficients of A r in V r only, and now the exterior expansion coefficients are irrelevant. Hence, in the functional integral, the integration of A can be replaced by N + 1 integrations over the fields A r , r = 0, . . . , N, where each A r is allowed to extend over unbounded space with the action for a free field in a homogeneous space with r , µ r . However, it is important that the correct of the two possible forms of the surface integral in Equation (60) is used. The multiple counting of degrees of freedom that results from N + 1 functional integrations poses no problem since the (formally infinite) factor in the partition function cancels when the Casimir energy is computed from Equation (52).
With this representation, we can write the partition function as a functional integral over A, separately in each region V r , and the surface fields A r on body r, leading to the partition function with the action Now, the fluctuations A can be integrated out easily, noting that the two point correlation function A r (x)A r (x ) = 0 for all r, r = 0, . . . , N with r = r , and for equal-region correlations where ∇ always acts on the argument x of . We obtain for the partition function with the kernels where ↔ G (AA;r) (x, x ) is the free Green function of Equation (56), and the arrow over the placeholder · indicates to which side of the kernel M acts. This notation implies that the derivatives are taken before the kernel is evaluated with x and x on the surfaces Σ r .

Hamiltonian Formulation
The representation of the partition function in the previous subsection sums over all configurations of the surface fields A r , and the action depends both on A r and the tangential part of its curl, which is functionally dependent on A r . Hence, the situation is similar to classical mechanics where the Lagrangian depends on the trajectory q(t) and its velocityq(t). The Lagrangian path integral runs then over all of path q(t) withq(t) determined by the path automatically. To obtain a representation in terms of a space of functions that are defined strictly on the surfaces Σ r only, it would be useful to be able to integrate over A r and its derivatives independently. In classical mechanics, this is achieved by Lagrange multipliers that lead to a Legendre transformation of the action to its Hamiltonian form. Here the situation is similar. To see this, it is important to realize that the bilinear form described by L r is degenerate on the space of functions over which the functional integral runs, i.e., Σ r d 3 x Σ r d 3 x A(x)L r (x, x )A(x ) = 0 for all A(x) that are regular solutions of the vector wave equation ∇ × ∇ × A + r µ r κ 2 A = 0 inside region V r . With a basis {A (reg,r) ν (x)} for this functional space, the elements of L r can be expressed as where we used the relations of Equation (58), and defined and made use of the fact that A ν (x) is also a solution of the vector wave equation inside V r . This implies that the kernel L r can be ignored in the above functional integral over regular waves A r inside the objects. However, the appearance of the kernel L r is important in what follows. Let us consider the part of the actionŜ[{A r }, {A r }] which, after functional integration over A r , generates the kernel L r . It is given by The exponential of this action can be written as a functional integral over two new vector fields K r and K r that are defined on the surfaces Σ r and are tangential to the surfaces, where Z r is some normalization coefficient, and we have used DK r DK r to indicate that the functional integral extends only over vector fields that are tangential to the surface Σ r . This representation shows that A r acts as a Lagrange multiplier. Integration over this field removes the imposed constraints between the dependent tangential fields n r × A r , n r × (∇ × A r ) by replacing them with the independent tangential fields K r and K r , respectively. Substituting Equation (72) for each object into the expression for the partition in Equation (61), we obtain with K r = (K r , K r ) the partition function with the kernel from Equation (72), and with the effective action where we have integrated out A r for r = 1, . . . , N, constraining the functional integral over A r to be replaced by the substitutions n r × A r → K r and n r × (∇ × A r ) → K r . Integrating out A 0 , finally yields with the additional kernel It should be noted again that the functional integral in Equation (76) runs over tangential vector fields K r , K r defined on the surfaces Σ r only. The kernelsL andM can be combined into the joint kernelN rr =L r δ rr +M rr .
SinceN acts in the path integral only on tangential vectors, the projections ofN on the tangent space of the surfaces Σ r have to be taken. Let t r,1 (x), t r,2 (x) be two tangent vector fields that span the tangent space of Σ r at x. The 4 × 4 matrix kernels then becomẽ which expresses all kernels in terms of tangential and normal derivatives of the scalar Green function g r (|x − x |) = e − √ r µ r κ|x−x | /|x − x |. These expressions simplify when an orthonormal basis t r,1 (x), t r,2 (x), n r (x) = t r,1 (x) × t r,2 (x) is used. The Casimir free energy is then given by where the determinant runs over all indices, i.e., x, x located on the surfaces Σ r , and r, r = 1, . . . , N. The kernelN ∞ is obtained from the kernelN by taking the distance between all bodies to infinity, i.e, by settingM rr = 0 for all r = r . In the following we shall again denote the form of the partition function in Equation (67) as Lagrange representation, and the one of Equation (76) as a Hamiltonian representation. By a simple computation, one can verify that the Hamiltonian representation of the Casimir free energy in Equation (79) is indeed equivalent to the surface formula Equation (44).

Application: Derivation of the Lifshitz Theory
As a simple example to demonstrate the practical application of the surface formulations, we consider two dielectric half-spaces, one covering the region z ≤ z 1 = 0, with the surface Σ 1 and dielectric function 1 and magnetic permeability µ 1 , and the other covering the region z ≥ z 2 = H, with the surface Σ 2 and dielectric function 2 magnetic permeability µ 2 . We shall consider both the Lagrange and Hamiltonian representation in the following.

Lagrange Representation
We compute the matrix elements of the kernels L and M of Equation (68) in the basis of transverse vector plane waves, given by with p r = r µ r κ 2 + k 2 and the sign of z is fixed so that the waves are regular inside the half-spaces. Note that we include here a z dependence to be able to compute the curl on the surfaces. For the Green tensor, we use the representation x q x q y q x q z q y q x r µ r κ 2 + q 2 y q y q z q z q x q z q y r µ r κ 2 + q 2 which yields after the curl operations We also need the following expressions for the operators that appear in the kernels, acting on the basis functions, which are tangential to the surfaces. On surface Σ 1 we havê and similarly on surface Σ 2 , It is straightforward to show that the matrix elements of L r in the above basis all vanish, as the basis functions are regular solutions of the vector wave equation. This observation is in agreement with the above finding that the kernel L r is degenerate on the space of those solutions. We proceed with the computation of the elements of kernel M. We find for the case r = r , and for the case r = r we get where the sign in e iq z (−1) r H determines upon integration over q z the sign of the terms ∼ iq z . The vanishing of the off-diagonal elements reflects the fact that the two polarizations described by the basis functions M and N do not couple for planar surfaces. The total kernel M can be written as The Casimir free energy is given by in terms of the determinant of the matrix which has four dimensions due to two sets of basis functions (polarisations) M and N per surface. This yields the final result This result is in agreement with the Lifshitz formula [37].

Hamiltonian Representation
Now we derive the Lifshitz expression for the free energy of two dielectric half-spaces in the Hamiltonian representation. Since the kernelsL andM of Equations (74) and (77) act on vector fields that are tangential to the surfaces, we need to compute the matrix elements in a basis of tangential vectors t r,1 (x) and t r,2 (x) that span the tangent space of surface Σ r at position x. For a planar surface, one can simply set t r,1 (x) =x 1 and t r,2 (x) =x 2 . For a given pair of positions x, x on the surface and fixed surface indices r, r we obtain the following 4 × 4 dimensional matriceŝ where we set x = (x , 0) and x = (x , H) for surfaces 1 and 2, respectively. We determined the sign of the terms ∼ iq z from the q z -integration by the observation that z, z have to be taken to the surface with z − z staying inside the object. The upper (lower) sign of p r refers to r = 1 (r = 2). Analogously, for kernel M we get for the case r = r and for r = r , where again the upper (lower) sign everywhere refers to r = 1 (r = 2). For the kernelM we determined the sign of the terms ∼ iq z from the q z -integration by the observation that z, z have to be taken to the surface with z − z staying outside the object. When combining the kernelsL andM into the joint kernel N rr =L r δ rr +M rr , it is diagonal in q -space with the blocks N(q ) on the diagonal given by the 8 × 8 matrix shown in Figure 2.

Figure 2.
Matrix N(q ) forming the diagonal blocks of the matrix N.
The Casimir free energy is given by in terms of the determinant of the above matrix where N ∞ is the matrix with H → ∞, i.e., the matrix N with the off-diagonal 4 × 4 blocks ∼ e −p 0 H vanishing. This yields the final result which is again identical to the Lifshitz formula. Note that in the Hamiltonian approach there is no need to express the kernels in a basis for the space of functions that are regular solutions of the wave equation inside the objects. However, the number of fields per surface is now doubled compared to the Lagrangian approach, resembling the situation in quantum mechanics where the Hamiltonian path integrals run over the canonical coordinates q and p independently.

Conclusions
To date, analytical and purely numerical approaches to compute Casimir interactions have been developed independently, and it remained an open question if and how they are related. Analytical methods build on ideas from scattering theory and hence require an expansion of the Green function and bulk or surface operators in terms of special functions that are solutions of the wave equation. Hence, the very existence of such functions and the convergence of the expansion limit these approaches to sufficiently symmetric problems. Purely numerical approaches, such as that developed in [32], can be applied to basically arbitrary geometries but the numerical effort can be extremely high. Hence, it appeared useful to us to study the relation between these approaches in order to develop methods that can serve as semi-numerical approaches that combine the versatility of the purely numerical approaches with the smaller numerical effort of analytical methods. Hence, we have presented in this work a new compact derivation of formulas for the Casimir force, which is based both on bulk and surface operators that also enable analytical evaluations. This we have demonstrated for the simplest case of two dielectric slabs. Further semianalytical implementations of our approaches are underway.
Our Hamiltonian path integral representation is equivalent to the one derived by Johnson et al. as a purely numerical approach using Lagrange multipliers to enforce the boundary conditions in the path integral. Interestingly, the here-presented derivation of this representation from a Lagrangian path integral demonstrates the relation of this approach to the scattering approach when the T-matrix is defined, as originally by Waterman, by surface integrals of regular solutions of the wave equation over the bodies' surfaces [38]. This shows the close connection of these approaches, motivating further research in the direction of new semi-analytical methods to compute Casimir forces.

Appendix A. Dyadic Green's Functions and Notations
It is well known [34] that the expectation values of quantum fields in systems in thermal equilibrium can be expressed, via the Matsubara formalism, in terms of the analytic continuation to the positive imaginary axis of the appropriate Green's functions. Along the imaginary axis, response functions have a simpler behavior than along the real frequency axis. For example, it can be shown that the electric and the magnetic permittivities of a magneto-dielectric medium are real-valued and positive definite for imaginary frequencies [39]. In the same manner, retarded Green's functions are real-valued for positive imaginary frequencies. In this Appendix, we briefly review their definitions and main properties.
To be specific, consider a system consisting of a collection of possibly spatially dispersive magneto-dielectric bodies in vacuum, occupying the (non-overlapping) regions of space V 1 , . . . , V N . We let V 0 be the vacuum region that separates the bodies. The electromagnetic response of the system is described by the following real-valued and positive definite frequency-dependent electric and magnetic permittivities The principle of microscopic reversibility [40] requires that the permittivities satisfy the symmetry relations: The imaginary-frequency Maxwell's Equations are: where κ = ξ/c, P ext and M ext are the external polarization and magnetization, and At the boundary separating media r and s, the tangential components of E and H are continuous: wheren is the unit normal to the boundary. The dyadic Green's functions G (αβ) ij (x − x ; i ξ n ) are defined such that: It is convenient to introduce a compact notation for the fields and the Green's functions. We collect the fields and the sources into six-rows column-vectors: The defining Equations of the Green's functions Equation (A6) can then be interpreted as defining the action of the linear operatorĜ onto the vector K ext : (A19) With this notation, the reciprocity relations Equation (A7) are expressed by the operator Equation: whereŜ is the 6 × 6 diagonal matrix of elements S where the operatorψ r isψ The constitutive Equations (A4) are expressed as: For later use, it is convenient to define the polarization operatorχ: Equation (A1) implies that the polarizability χ of a collection of bodies is: whereχ (r) is the polarizability of body r. Clearlyχ (r) is supported in the volume V r . This implies:χ (r)χ(s) = 0 , if r = s .
Finally, we define the trace operation Tr of any operatorÂ: ii (y) .
whereΠ r is the surface operator, which acts on the field Φ(x) by projecting it onto the tangential plane to Σ r at x: ↔ Πr(x, x ) = δ (3) (x − x ) δ(F r (x)) ↔ 1 −n(x) ⊗n(x) , r = 1, . . . , N (A42) Equation (A41) can be recast in the form: whereM is the surface operator whereΠ = ∑ N r=1Π r . The operatorM is invertible [41]. Equation (A43) then determines the unique surface current K surf that solves the boundary-value problem: Note that according to its definition,M −1 acts on tangential fields defined on the union Σ of the surfaces Σ r of the bodies and returns as a result a tangential polarization field K surf defined on Σ. Suppose now that the external sources K ext are localized in the vacuum region V 0 : We thus see that in the vacuum outside the bodiesΓ has precisely the form of Equation (13) In this Appendix, we demonstrate the force formula Equation (14). The simple proof presented below holds for smooth kernels K (ρσ) kl (y, y ). Unfortunately, the kernels K (ρσ) kl (y, y ) involved in both the volume and the surface representations of the scattering Green's functions lack the necessary smoothness. Indeed the T-operator in Equation (8) is in general discontinuous on the surfaces of the bodies. As to the surface operatorM in Equation (9), it is a singular kernel supported on the surfaces of the bodies. Fortunately, it is possible to remedy this difficulty by representing the kernel K (ρσ) kl (y, y ) as the limit of a one-parameter family of smooth real kernelsK m|i (x) are well defined in all space, and by construction they satisfy the following euclidean-time Maxwell Equations: where we set Φ m ≡ (E m , H m ), K m ≡ (P m , M m ) and κ = ξ n /c. When the expansion Equation (A57) is substituted into Equation (6), we find for the dyad ↔ Θ the following expression: It is worth noticing that the minus sign that multiplies ↔ Θ (HH) m in Equation (A61) is a direct consequence of the factor (−1) s(β) in the r.h.s. of Equation (A57). At this point, we use the divergence theorem to convert the surface integral in Equation (7), giving the force F (r) , into a volume integral. That gives: where in the last passage, we made use of Equation (A56), and restored the dependence of K on the parameter λ. By taking the limit λ → 0 of the r.h.s., we see that Equation (A67) reproduces the formula for force in Equation (14).
containing the body in its interior and centered at the point X r ,