Quantum electrodynamics of inhomogeneous anisotropic media

In this work we calculate the closed time path (CTP) generating functional for the electromagnetic (EM) field interacting with inhomogeneous anisotropic matter. For this purpose, we first find a general expression for the electromagnetic field's influence action from the interaction of the field with a composite environment consisting in the quantum polarization degrees of freedom in each point of space, at arbitrary temperatures, connected to thermal baths. Then, we evaluate the generating functional for the gauge field, in the temporal gauge, by implementing the Faddeev-Popov procedure. Finally, through the point-splitting technique, we calculate closed expressions for the energy, the Poynting vector and the Maxwell tensor in terms of the Hadamard propagator. We show that all the quantities have contributions from the field's initial conditions and also from the matter degrees of freedom. Throughout the whole work we discuss and give insights about how the gauge invariance must be treated in the formalism when the EM field is interacting with inhomogeneous anisotropic matter. We study the electrodynamics in the temporal gauge, obtaining the EM field's equation and a residual condition. Finally analyze the case of the EM field in bulk material and also discuss several general implications of our results in relation with the Casimir physics in a nonequilibrium scenario.


Introduction
The main subject of this paper is to develop a CTPintegral formulation of nonequilibrium quantum electrodynamics in an inhomogeneous, anisotropic real medium. The CTP-method has been used in quantum field theory as a tool to nonequilibrium descriptions of dynamical problems, where dissipative effects arise at the macroscopic level after coarse graining the detailed information in one or more subsystems, by tracing out those degrees of freedom. In fact, this approach presents a combination of both quantum field theory and nonequilibrium statistical mechanics. The former is needed by the quantum electromagnetic (EM) field and the latter for treating processes involving quantum dissipation and noises. As far as we know, this complete formulation has not been applied previously to the quantum theory of the electromagnetic field in interaction with matter capable to dissipate and absorb energy to and from an external thermal environment, at a microscopic level. The general problem is of both basic and practical importance. The polarization properties of real media change the properties of the EM field dramatically in comparison with the situation in vacuum, especially in the case of a dispersive and absorbing dielectric. Therefore, Casimir effect appears as one macroscopic manifestation of these effects, for which it is important to have a complete theoretical formulation to contrast first principle based models with experiments.
The very well known Lifshitz formula [1] describes the forces between dielectrics in terms of their macroscopic EM properties, in thermal equilibrium. Its original derivation is not based on a first principle quantum framework, but on a macroscopic approach, starting from stochastic Maxwell equations and using thermodynamical properties for the stochastic fields. The connection between this approach and one based on a fully quantized model including lossy dielectrics is not completely clear [2,3,4].
When dealing with a composite system, in which there are noise, fluctuations, and also dissipative effects between different parts of the full system (mirrors, vacuum field and environment), the theory of open quantum systems [5] is the most appropriate framework to clarify the role of these effects in Casimir physics. Indeed, in this description, dissipation and noise appear in the effective theory of the relevant degrees of freedom (the EM field) after integration of the matter and other environmental degrees of freedom.
The quantization at the steady situation can be performed starting from the macroscopic Maxwell equations, and including noise terms to account for absorption [6] (one can also to couple the EM field to an external reservoir [3], following the standard route to include dissipation). Regarding microscopic models, the fully canonical quantization of the EM field in dispersive and lossy dielectrics has been performed by Huttner and Barnett (HB) [7]. In their model, the EM field is coupled to matter (the polarization field), and the matter is coupled to a reservoir that is included into the model to describe the losses. In the context of the theory of quantum open systems, one can think the HB model as a composite system in which the relevant degrees of freedom belong to two subsystems (the EM field and the matter), and the matter degrees of freedom are in turn coupled to an environment (the thermal reservoir). The indirect coupling between the electromagnetic field and the thermal reservoir is responsible for the losses.
In Ref. [8], we have followed an equilibrium canonical quantization program similar to that of Ref. [9], generalizing it by considering a general and well defined open quantum system. Recently, in Ref. [10], we have considered two simplified models analogous to the one of HB, both assuming that the dielectric atoms in the slabs are quantum Brownian particles, and that they were subjected to fluctuations (noise) and dissipation, due to the coupling to an external thermal environment, also generalizing the constant dissipation model of Ref. [8]. Indeed, after integration of the environmental degrees of freedom, it was possible to obtain the dissipation and noise kernels that modify the unitary equation of motion of the dielectric atoms. Using the Schwinger-Keldysh formalism (or closed time path (CTP) -in-in -formalism) [11,12], we have studied the time evolution of the expectation value of the energymomentum tensor of a scalar field in the presence of real materials. The work included a fully nonequilibrium scenario done for two different couplings between the scalar field and the polarization degrees of freedom of matter. There it was shown that two contributions always take place in the transient evolution of the energy-momentum tensor: one associated with the material, and other one only related to the field. Therefore, we have shown that the material always contribute unless is non-dissipative. Conversely, the proper field contribution vanishes unless the material is non-dissipative or, moreover, at least for the 1 + 1 case, if there are regions without material. We finally concluded that any steady quantization scheme in 1 + 1 dimensions must consider both contributions and, on the other hand, we argued why these results are physically expected from a dynamical point of view, and also could be valid for higher dimensions based on the expected continuity between the non-dissipative and real material cases [10].
In this work, we will extend Ref. [10] in order to calculate the CTP generating functional for the EM field (abelian gauge field) interacting with inhomogeneous and anisotropic matter through the open system framework. For this purpose, we will firstly calculate a general expression for the EM field's influence action from the interaction of the field with a composite environment consisting in quantum polarization degrees of freedom in each point of space (having arbitrary temperatures) and connected to thermal baths (with arbitrary temperatures too). Then, we will evaluate the CTP EM-field generating functional in the temporal gauge by implementing the Faddeev-Popov procedure. Moreover, we will calculate closed expressions for the EM-energy, the Poynting vector and the Maxwell tensor in terms of the Hadamard propagator, showing that all of these quantities present contributions from the field's initial conditions and also from the matter degrees of freedom in the mirrors. Then, we analyzed the dynamics of the EM field in the temporal gauge and study the case of an infinite homogeneous and isotropic material, connecting our results with previous ones. A detailed analysis is performed in relation to how the gauge invariance must be treated in the CTP formalism when the EM field is interacting with inhomogeneous anisotropic matter. We will also discuss several general implications of our results in relation to the nonequilibrium calculation of energies and forces in Casimir physics. In fact, a detailed study of the Casimir-Lifshitz problem in a fully nonequilibrium situation will be presented elsewhere using the results of this paper [13].
This paper is organized as follows: In the next Section we will discuss the CTP integration of the polarization degrees of freedom in interaction with both, the EM-field and the thermal environment. In Section III, we evaluate the CTP-generating functional for the gauge field in the composite system. Section IV contains the formal calculation of the energy, poynting vector and mean value of the Maxwell tensor. In Section V, we give insights about the electrodynamics in the temporal gauge and analyzes a concrete example. Finally, in Section VI we will present the final remarks.

CTP Integration for the Electromagnetic Field -Matter Interaction
The main goal will be to calculate the generating functional. Therefore, the first step in this direction is to calculate the influence action over the field. We will describe the ordinary polarizable matter by (non-relativistic) quantum degrees of freedom associated to the polarization vector P of each volume element of the polarizable body, each one subjected to an independent bath generating an influence action as in the quantum Brownian motion (QBM) theory [5], i. e., a quantum harmonic oscillator interacting linearly with a thermal bath consisting in a set of quantum harmonic oscillators. On the other hand, the EM field is described by a massless spin-1 gauge field A µ . The interaction term can be taken as a coupling between the field and the current generated in the polarizable matter or, equivalently, as the dipolar interaction between the polarization dipoles and the field. In other words, the interaction can be considered in two ways, depending which degree of freedom is differentiated. However, since we have nonrelativistic matter degrees of freedom coupled to a gauge field which is, from the beginning, a relativistic system, differentiation is not only a time derivative, as it happens for the scalar field in Ref. [10].
In other words, the coupling terms can be proportionals to d 4 x P j (x) E j (x), coupling the polarization vector to the field's canonical momentum (i. e., the electric field) in a dipole interaction way; or, analogously, , coupling the current four-vector of the (non-relativistic) polarizable matter to the (relativistic) EM field. From now on, our sum notation will be Einstein's for greek sub-and superscripts (a subscript sum with a superscript), but latin ones may sum without being opposite scripts (the covariant or contravariant nature of spatial coordinates will be adjusted by introducing minus signs when passing from a subscript to a superscript or viceversa).
Since the electric field is E j = −∂ j A 0 − ∂ 0 A j , it is clear that the interaction term is gauge invariant, while from the second way it is clear that the interaction term is a Lorentz scalar. In order to keep gauge invariance and considering also that the current J µ must be a conserved current (its four-divergence must vanish), the current fourvector have to be written as J µ = ∇ · P, −Ṗ .
After the integration over P, we will eventually integrate over the EM field. Then, we will consider the second form for the interaction term and leave a longer discussion for the next Section. However, as we will integrate over P firstly, we integrate by parts in the respective space coordinate for each derivative of the term involving the divergence. Therefore, considering that the field paths vanish at infinity in every direction, the interaction term reads The model for the total system can be described by the initial total action: are the free actions for the EM field, the polarization vector and the degrees of freedom of the thermal baths which affects the polarization vectors in each spatial point, respectively. The spatial labels denote the fact that the properties change with the position while the degree of freedom presents no spatial differentiation in its dynamics, i. e., the degrees of freedom are a spatial continuous set of 0 + 1 fields. The last two actions are the interaction actions between the EM field and the polarization vectors (as we have discussed above) and between the polarization vectors and the thermal baths in each point (which are linear couplings).
Therefore, the first step would be CTP integration over the thermal baths ({q n,x }). However, from the well-known QBM theory we know that this is already done since, as we mentioned before, it gives that the polarization vectors under the influence of the thermal baths behave effectively as Brownian particles (see Refs. [5,14]). Then, for each j−component in each point of space, the polarization vector will have it unitary evolution modified by the QBM influence action: where D x and N x are the QBM's dissipation and noise kernels respectively, while ∆P = P ′ − P and ΣP = (P + P ′ )/2 are called difference and sum variables each one. We are interested in calculating the CTP integral associated to the calculation of the influence action which acts over the EM field due to the interaction with matter. Considering a coupling constant λ 0 between the field and the polarization degrees of freedom, the integral is given by: where A·B ≡ t f ti dt A(t) B(t) and, for simplicity, products between vectors and matrices will be omitted in the vectorial form, and then for example and so on. The matter distribution function g, which takes binary values (1 or 0) whether or not there is matter in each point of space x, is introduced to denote the fact that this calculation takes place in every point of space that contains polarizable material. Therefore, the influence action that acts over the EM field will be defined in space by the matter distribution function g, which will define in which points the influence action does not vanish.
As we are dealing with a three dimensional problem, the polarization vectors P described as three-dimensional harmonic oscillators can be decomposed in cartesian components, each one suffering the action of different baths. Therefore, we trivially have that S 0 [P] = 3 j=1 S 0 [P j ] and S IF [P, Considering a separable initial state for the polarization vector, the initial density matrix is the product of density matrices for each component of the polarization vector, ρ P (P i , Then, finally noting that the interaction term also separates in each component, the CTP integral can be written as the product of three integrals for each component or direction. One step further can be done by taking advantage from the fact that the material is conceived as a continuous of independent degrees of freedom (there is no interactions between them), giving the chance for skeletonizing the spatial grid in volume elements ∆x and writing the integrals for each component as a product over the spatial points ( x ) of a generic CTP-integral with a spatial label x.
We now proceed to write the actions in terms of the difference and sum variables. For this purpose, we consider that in each point of space for the interaction terms we have and for the free actions we can write where we have allowed each degree of freedom, at each point x, to have its own properties (mass and natural frequency). The functional integrations over ∆P j are Gaussians in each point of space, and can be performed straightforwardly considering the noise kernels in each direction and point, and defining the linear sources as R j (t) = ∆x . At this point, the remaining functional integration is over ΣP j in each point of space. However, as it was done in Refs. [10] and [15], we can write every path ΣP j in terms of an homogenous solution P j 0 (t) satisfying the initial conditions and a particular solution P j ξ (t) describing the deviation of the paths from the homogeneous ones, i. e., we write ΣP j (t) = P j 0 (t) + P j ξ (t). Considering that, from the initial actions, the canonical momentum associated to P j is given by MṖ j + λ 0 A j , the solutions can be written as with Π j i = M ΣṖ j i + λ 0 ΣA j i and being G Ret,x (t − t ′ ) the retarded Green function for the linear integro-differential operator, associated to the kernel L of the degree of freedom at x. Therefore, we can replace the functional integration over ΣP j with integration limits, by the functional integration over ξ without them, and an ordinary integration over all the values of the initial canonical momentum Π j i . Then, after realizing the replacements it turns out that the functional integration over ξ is straightforward and, by omitting for simplicity the spatial labels, we obtain for each integral in each point of space: where we have introduced a delta function in order to take into account the restriction on the final points, and the Wigner function for the j-component of the polarization vector in the point x defined as in Refs. [10] and [15] by (10) Now, on one hand, we note that the homogeneous solution can be written as P j is the homogeneous solution for a bi-lineal coupling (where the canonical momentum is only related to the time derivatives of the degree of freedom, see for example Ref. [10]). On the other hand, we can take the exponent of the last factor and integrate by parts on the second time integration, obtaining where we have used that G Ret,x (t − t f ) = 0 because it is a retarded Green function involving a Heaviside function. By replacing into the CTP integral, we find that one factor constituted by the second term of the last expression of P j 0 in Eq.(8) cancels with a factor constituted by the first term of the right hand side of last equation.
As a last step to cancel out the factors involving initial conditions, we must consider the first term in the exponent of the first factor in Eq. (9). Again, by integrating by parts and using the CTP condition for the EM field, we have Then, the first term in this expression cancels out with another factor containing initial conditions. Finally, it is straightforward to integrate over ∆P j f by using the delta function. After this step, it is easy to evaluate the integral over ΣP j f by using the other delta function. Considering that the result is for the j-component of the polarization vector at the point x, then, we obtain the product over all the positions where there is material. Grouping each factor and taking the continuum limit for the spatial grid, the exponents result in integrals limited by the matter distribution g(x), thus i,x , t i and we let each material's property (charge, noise kernel and Green function) to depends on position by introducing a subscript x as a label, allowing the material to be inhomogeneous. Now, we have to consider that this last expression is for the j-component of the polarization vector. The final expression for the field's influence action on Eq. (3) results from the product of this expression for each component. At this point, we can give the material another freedom, which is to be anisotropic (birefringent material). Each direction j of the polarization vector can be matched with each of the three principal axes of Fresnel's (or refractive index) ellipsoid in each point of the material (we give more insights of this in Sec.5.1). Therefore, each component can have different properties (except for the charge), and by introducing delta functions in space in the last two factors, we can write, more compactly Ret, which are the bi-linear dissipation matrix-kernel associated to the EM field -matter interaction model (see Ref. [10]) and the bilinear noise matrix-kernel which is only related to the contribution of the baths. The superscripts [j] denote the dependence on direction (anisotropic material).
To give a closed expression for the influence action on the field, we have firstly to move the temporal and spatial derivatives on the field's components to the bi-linear matrix-kernels to define correctly the 'current' matrixkernels which act over the EM field. This is in fact straightforward due to the causal behavior of the Green function, the CTP condition on the EM field's components and the convergence of each path of the EM field whenever any spatial coordinate goes to infinity. These allow us to integrate by parts both in time and spatial coordinates resulting in a transfer of the derivative on the field to the matrix-kernels. Then, we can write a covariant form for the last two factors with the (covariant) EM dissipation kernels D µν and the (also covariant) EM noise kernel N B µν associated to the contribution of the baths are given by where the operator Γ µν tt ′ with the prime denoting derivation on the respective coordinate of the point x ′ and the covariant delta is introduced with Einstein's notation, unlike all delta in matrix notation employed at the moment for spatial sub and superscripts.
It is clear that the first factor on the right hand of Eq.(18) is entirely related to the initial state of the polarization degrees of freedom. However, in order to obtain an expression for the field's influence functional, we have to calculate the factor for the chosen initial state. This can be easily done for the case that the initial state is a thermal one for each direction of the polarization in each volume element. Therefore, considering temperatures β P j x for each direction in each point, the integrals over ΣP j i and Π j i in each point of space are Gaussian, obtaining by discarding the normalization factor: where the EM noise kernel associated to the polarization degrees of freedom is given by with: Ret, Ret, Ret, Therefore, the EM field's influence action reads as First, it is worth noting that it has the form of Eq.(2). This is in principle not surprising, however the present influence action for the EM field is the result of two CTP integrations since the system is a composite of three parts, and the fact that it has this closed form is a merit of the choice of a thermal initial state for the polarization degrees of freedom. For a non-thermal initial state, this form is not achieved. Moreover, while the kernels in the QBM theory depend on the difference of the arguments, the present EM kernels do not.
It is remarkable that this last expression is analogous to the one found for the scalar field case for the polarization degrees of freedom contribution in Ref. [10]. However, since in this case we are dealing with an abelian gauge field, all the kernels contain the differential operator Γ jk µν , which basically ensures that the influence action is gauge invariant. In fact, it is easy to see that This is the property needed to ensure gauge invariance, so it is the physical requirement for every EM field's influence action of the form of Eq.(23) in order to be gauge invariant. In other words, every quadratic influence action must contain EM dissipation and noise kernels with null four-divergence in both subscripts.
Nevertheless, this is indeed an expected mathematical requirement which cames from the physical fact that the kernels are basically correlations functions of the current four-vector J µ , which is a conserved current thanks to the gauge invariance of the EM theory. If we break gauge invariance, the current is not neccesarily conserved as it happens for the Proca field theory (see Ref. [19]), so all the obtained properties are clearly expected as requirements to satisfy by a physically consistent EM theory.

CTP Generating Functional for a Gauge Field
In this Section we will calculate the CTP generating functional for the gauge field. Clearly, the result obtained in the last Section for the EM-field influence action will be the starting point. We are dealing with a spin-1 abelian gauge field A µ = (A 0 , A) (being A 0 and A the electric and vector potentials respectively) under the influence of matter degrees of freedom, modeled through the EM dissipation and noise kernels in the influence action obtained in the last Section. Although we will use the main result in Eq.(23), we resume briefly some points of the previous discussion about the interaction terms.
It is worth noting that, if we consider a typical spin-0 scalar field instead of the EM field, the calculation of the generating functional can be easily done as it was shown in Ref. [10] for both bilineal and current-type couplings, without critical changes and obtaining the same formal result. In contrast, for the case of the EM field, a few subtle points must be taken into account.
Considering the gauge symmetry of the EM theory, the interaction term with matter must be gauge invariant. Thus, bilineal coupling models (as the one considered in Ref. [10]) are forbidden from the very beginning because the interaction term in the initial actions is not gauge invariant. In fact, the influence action for that case is not gauge invariant since the kernels do not verify the correct properties that ensures gauge invariance.
However, the current-type coupling presents a subtle conceptual variation. This type of coupling cannot be on the field's time derivative as it happens for the scalar field in Ref. [10] because it also breaks gauge symmetry. Then, the interaction term must be in terms of the electric and magnetic fields in order to keep the gauge invariance of the whole theory. The general rule that is behind all of these choices is that the interaction for a current-type coupling must be on the canonical conjugate momentum. In the classical free EM field theory, the canonical momenta are defined as From a formal point of view, it is well known that Π 0 ≡ 0 and, consequently, the canonical momentum for the temporal component of A µ is not well-defined, implying that a quantization procedure will not be so straightforward. However, we have no need to focus on this point at this time, and we will see how to deal with that problem in our formalism.
By considering the linearity of the interaction terms, we can write directly the CTP generating functional for the gauge field as a CTP-Feynman path integral: We can proceed by performing the Faddeev-Popov procedure to extract the redundant sums over paths on the same gauge class. Therefore, this includes a gauge fixing term (depending on the gauge choice defined by the gauge condition F [A µ ] = 0 which will be taken linear) and the determinant which gives the ghost terms in the Lagrangian, which, in our case, can be discarded since ghosts do not couple the EM field for every linear gauge condition due to the abelian nature of the field.
We re-write the delta-functionals as lagrangian gauge fixing terms in the exponential containing the free field actions. For this purpose, we have to continue the typical Faddeev-Popov procedure one step further. Then, we change the gauge condition on both delta functionals to F [A µ ] = C(x) and F [A ′ µ ] = C ′ (x) respectively, where C and C ′ are arbitrary functions of the coordinates (the determinant giving the ghosts do not change anyway). As the generating functional is independent of C and C ′ , we can multiply by a CTP-type distribution e − i 2α (C * C−C ′ * C ′ ) (having an unique gauge fixing parameter α for both fields) and then integrate over the arbitrary functions. As we need to keep Eq. (24) valid, we must consider the Landau gauge for a given choice of the gauge condition, where α → 0, and we effectively re-obtain the delta-functionals from the gauge fixing exponentials [17]. Thus, we get being the new action for the EM, field containing the typical gauge fixing term, which breaks gauge invariance for every α, but that allow us to treat the component of the field as independent variables. It is clear that it turns out to be a crucial point that the gauge fixing parameter is unique, which is indeed very natural because both fields in the CTP formalism are needed to describe an unique quantum field, so the effective action must depends only on one parameter. On the other hand, it can be said that α must be unique since a CTP integral can always be written in terms of an unique CTP field, by parametrizing both branches of the CTP contour. Therefore, applying Faddeev-Popov procedure at this point, it is clear that α must be the same for both branches.
It is worth noting that we have properly included subscripts α as a reminder that when we take the Landau gauge (taking the limit α → 0), the initial and final field configurations must be written satisfying the gauge condition, since in that limit, the result must be the same as the one obtained by evaluating the deltas on the gauge condition at the very beginning. Indeed, we have to take care of this only for the initial and final points, and we do not have to introduce subscripts in the functional integrations since, as we will show below, the correct result will be naturally ensured by taking the limit on the Green function associated to the effective CTP action including the gauge fixing terms. The limit α → 0 or Landau gauge imposes the vanishing gauge condition (it sets the arbitrary function C), i. e., taking the limit will imply that, for each point, the field satisfies the vanishing gauge condition. Now, we change variables to ∆A µ = A ′ µ − A µ and ΣA µ = (A µ + A ′ µ )/2, what allow us to treat each component as an independent variable thanks to gauge fixing terms which break gauge invariance. To continue the calculation, we first have to integrate by parts the free field actions, including the gauge fixing term. Therefore, we have to choose a gauge condition.
At this point, as we are considering linear gauge conditions, if the gauge condition actually contains or not a derivative on the field will generate crucial differences. If we take the temporal or axial gauge, the gauge condition can be written in general for both cases as F [A µ ] = t µ A µ , allowing t µ to be a temporal or spatial four-vector for each case. For now, we will keep generality on this four-vector without choosing an specific gauge condition. Thus, in this case we do not have to integrate by parts because we directly have: For the free field actions, we get . This way, by defining the Wigner functional for the EM field as a natural (gauge dependent) extension of Eq.(10) (see Ref. [18]), we obtain: The next step is to write the paths of ΣA µ in terms of it initial conditions and the retarded Green function associated to the operator L µν . Given the EM field equation of motion associated to the CTP effective action by: then, the retarded Green function for t > t ′ is defined by: It is worth noting that we have remarked the fact that the retarded Green function depends on the gauge parameter α and on the time difference.
Thus, the solutions for each component of ΣA µ , splitting in homogeneous A µ 0,α and inhomogeneous A µ ξ,α solutions, can be written as: where we have clearly denoted the fact that both initial conditions on the field and its time derivative, depend on the gauge parameter, while the four-function ξ ν does not, since it represents a deviation of the solutions from the classical ones and, therefore, each component can be treated independently without α−dependence. However, the present step is quiet subtle, because it consists in the replacement of the functional integration over possible paths by two ordinary integrations over the initial field and momentum configurations (which involve the homogeneous solution) and a functional integration over ξ µ , which represents the shift in the path from the classical trajectory (and it includes the inhomogeneous solution).
In the present case, the replacement is not so easy since the canonical momenta are not proportional to the field components' time derivatives, as it always happens until this Section. Moreover, the choice on the gauge condition is a crucial point on defining these replacements. As we have noted in the beginning of this Section, it is well known that the canonical momentum for the temporal component A 0 is not well defined since Π 0 ≡ 0. This is intimately related to the problematic position of the EM theory in order to be quantized, since the canonical quantization procedure cannot be developed in a straightforward way [19]. The same happens for a (non-CTP) path integral quantization, but the Faddeev-Popov procedure shows to be efficient in that context having no restriction on the values of the gauge parameter α [20].
At this point, we showed that the Faddeev-Popov procedure allow us to treat the situation, but restricting the theory with the gauge fixing term to the Landau gauge. We will now see that taking this limit carefully we may impose the gauge condition in the CTP formalism.
Then, in the temporal gauge, the four-vector t µ will be taken as a time-like four-vector, being our particular choice the simplest one, t µ = (1, 0, 0, 0). The gauge condition reads F [A µ ] = A 0 . Note that over all the possible choices of the functions C, to which the gauge condition will be equal, the Landau gauge implies the gauge condition for C = 0, i. e., the field satisfies Landau gauge for the temporal gauge condition. This means that the solutions (and Green functions) obtained for an arbitrary value of the gauge parameter α must be taken to satisfy the Landau gauge. In this limit, as we can see from the retarded Green function equation [Eq. (31)], G 0ν Ret,α→0 must be identical to zero.
On the other hand, from the EM theory, the canonical momentum is given by Π µ ≡ F 0µ = δ µ j E j . In the temporal gauge, the electric field is given by E j = −∂ 0 A j , so the canonical momentum reads Π µ = −δ µ jȦ j . But in our theory with gauge fixing term (without imposing the gauge condition), the canonical momenta read as usual, Π µ = δ µ j ∂ j A 0 −Ȧ j . This makes that the initial time derivatives of each field component can be written as where the temporal component of the derivative is the only one that cannot be re-written in terms of the canonical momentum and the spatial derivative of the temporal component. Therefore, replacing this into the homogeneous solution A µ 0,α , we may write it in terms of the canonical momenta and, after integrating by parts the terms associated to ∂ j ΣA 0 i (and discard the boundary terms by convergence), the homogeneous solution reads: Then, the integration replacement can be done and then we can easily integrate on ∆A µ α,f and then on ΣA µ α,f thanks to the delta-functionals, and also functionally over ξ µ , clearly obtaining: Finally, we take the Landau gauge, α → 0, which implies G µν Ret,α → δ µ j δ ν k G jk Ret,LG on the retarded Green function and it also restricts the integrations because ΣA µ α,i → ΣA i and ΣȦ µ α,i → ΣΠ i = −ΣȦ i . However, this point encodes the action of imposing the temporal gauge on the field equations at every time. Moreover, if we consider that the field is free for times before the initial one, setting A 0 = 0 in the field equations, result in three equations for the the vector potential A i and an additional condition (which is a residual condition resulting from the equation for A 0 ) given by ∇ · Π i = 0. This completely defines the field, its components and canonical momenta for earlier times than the initial one, and also implies that the field has two independent components. Therefore, taking the limit α → 0 must includes, particularly, that the initial conditions check A 0 i = 0 and ∇ · Π i = 0 (being transverse to the direction of propagation of each field mode). Then, we can naturally write: where A i and Π i are taken as free initial conditions in the temporal gauge (and consequently perpendicular to the wave vector for each mode), having which is the homogeneous solution for the field equation after imposing the temporal gauge. It is worth noting that this is the natural and expected extension of the result obtained for the scalar field in Ref. [10]. Nevertheless, the gauge nature of the EM field causes that we have to choose a gauge condition in order to develop the calculations. In the temporal gauge, we see that the first factor in Eq. (35), involving the initial conditions, enforces the calculation by introducing the initial conditions in the chosen gauge.

Energy, Poynting Vector and Maxwell Tensor
Once we have calculated the CTP generating functional for the EM field satisfying a given gauge condition, we can proceed with the calculation of the field correlation as functional derivatives of the generating functional. Initially, we have introduced four-vector as classical CTPsources J µ , J ′ µ . In the temporal gauge the generating functional is functionally dependent on the spatial coordinates of the CTP-sources four-vectors J µ , J ′ µ , i. e., the generating functional depends on J, J ′ . Therefore, correlation functions involving the temporal coordinate of the field, which are constructed from functional derivatives of the generating functional with respect to the temporal coordinate of the source four-vectors J µ , J ′ µ , will vanish. This is clearly expected since we have chosen the temporal gauge, where A 0 ≡ 0. Thus, as it is well known, the field correlation can be written in general as: As in Ref. [10], since the generating functional has a simple form of Eq.(35), independently of the initial state of the field, we can easily compute its functional derivatives. Taking advantage of the symmetry properties of the noise kernel, we obtain: LG (x 1 , x 2 ) is the Jordan propagator [14] and A j 0 is the homogeneous solution of Eq.(35).
This correlation function corresponds to the Whightman function for the field in this open system and it is the EM field generalization of the results found in Refs. [15] and [10] for a quantum degree of freedom and a scalar field, respectively. In fact, considering that G jk Ret,LG is real, is clear that the correlation is a complex quantity, with the imaginary part given by G jk Jordan,LG , whereas the real part is formed by the others two terms. It is possible to show, using the typical relations between the propagators, that the Hadamard propagator is given by where this expression holds for every initial state of the field, and it is clear that this propagator is also gauge dependent. Note that the Hadamard propagator has two separated contributions. One is associated to the material degrees of freedom represented by the noise kernel N, which also splits into two contributions due to the composite nature of the material (polarization degrees of freedom plus bath in each point of space). The other contribution is entirely associated to the field's effective dynamics and the initial state.
With the field correlation at hand, we may calculate physical quantities of interest. We begin by giving a formal expression for the Poynting vector [21,22]. Considering that the real material considered is non-magnetic, we can define the Poynting vector as: where ǫ jkl is the Levi-Civita symbol and it is worth noting that the Poynting vector is a gauge invariant quantity, because the electric and magnetic fields are. Following Refs. [10] and [17], and also using the point splitting technique, the expectation value of the Poynting vector reads where ∂ t1 denotes time derivative on the point x 1 and so on for the other derivatives.
It is important to consider that to use the point splitting technique, the correlation function must be a regularized quantity, in order to have finite results. Moreover, all the expectation values of interest are expected to be real quantities, as it happens with the Poynting vector. However, this seems not to be the case due to the fact that the propagator is complex. But the coincidence limit combined with the symmetric definition of the Jordan propagator, make that the imaginary contribution vanishes at the end of the calculation. Then, the expectation value of the Poynting vector can be written only in terms of the (regularized) Hadamard propagator on Eq.(39) as follows: This last equation gives the full time evolution of the Poynting vector, which inherits the two contribution splitting from the Hadamard propagator.
Once we have worked out the expression for the Poynting vector in terms of the Hadamard propagator, we may do so with the EM energy and the Maxwell (or stress) tensor. However, this is not so straighforward since these quantities for the EM field in real materials have not unique definitions. This is related to the freedom on the arbitrary definitions of the mechanical and EM contributions where the matter is coupled to the EM field (see for example the discussion for the classical theory given in Ref. [21] for simple linear isotropic media and the general approach on Ref. [22]) due to the non-coincidence between the displacement and electric vector fields inside a macroscopic material.
However, since the energy density and the Maxwell tensor are locally defined at each point of space, we can avoid the discussion by calculating them in vacuum regions, independently if there are matter bodies in other points of space, i. e., whether or not there are material boundaries. As in that regions there is no distinction between electric and displacement vector fields, the definitions of the energy and the Maxwell tensor are unique. Therefore, the quantum definitions for both quantities are given by [21,22]: Using again the point splitting technique, we easily obtain for the expectation values of both quantities: where it is clear that both quantities also inherit the two contribution splitting of the Hadamard propagator.
It is important to remark that this expression indeed gives the key quantities to study the Casimir force between bodies separated by vacuum regions in a fully nonequilibrium situation for a EM field, generalizing the results given for the scalar field in Ref. [10].
As a final remark, we should note that, in a covariant formulation, these three quantities (Poynting vector, energy density and Maxwell tensor) are part of the (covariant) energy-momentum tensor for the EM field [22]. However, we have restricted our calculation for regions without material, avoiding the discussion for the definition of each quantity inside the material. The crucial point here is that, regardless of whether definition is considered, we can always write the expectation value in terms of the Hadamard propagator through the point splitting technique. Therefore, all the quantities will contain the contributions structure of the Hadamard propagator, given in Eq.(39). In a non-equilibrium situation, the transient dynamics of the EM field will have contributions from each part of the composite system. On the other hand, the long time regime (t 0 → −∞) is expected to be defined at most by the baths' contributions and the field's initial state contribution, having different steady situations depending on the chosen initial state, as it was shown for the scalar field in Ref. [10].

Open Electrodynamics in the Temporal Gauge
At this point we have fully developed the CTP formalism for the the EM field interacting with a linear, inhomogeneous and anisotropic material in a general context. In the last Section, we have written all the physical EM quantities in terms of the Hadamard propagator.
In this section, we are going to study a simple example in order to get a direct application of the developed formalism. We start by examining in a general way the EM equations of motion and then we focus on the dynamical aspects of the EM field in a infinite homogeneous and isotropic material.

EM Field's Equation in the Temporal Gauge
The EM retarded Green tensor is defined by the equations of motion (Eq.(30)) obtained from the EM CTP action after imposing the temporal gauge (A 0 ≡ 0). For this purpose, it is important to note that these equations include the gauge fixing term for the temporal gauge, corresponding to the term containing the gauge fixing parameter α. Considering the equation for the temporal coordinate (µ = 0), which is the only one which contain the gauge fixing term, we have: (46) Now, by choosing the Landau gauge, where α → 0, it naturally implies that A 0 ≡ 0 in order to do not have divergent terms. Then, the (vanishing) temporal gauge is naturally introduced by the choice of the Landau gauge. The equation in this case, remains: The dynamical equation for the temporal component A 0 in the vanishing gauge condition becomes a residual condition for the remaining components A m . From the definition of the EM dissipation kernel of Eq.(18), D 0m can be easily calculated and the condition reads: Ret,x (t − t ′ ) is a function of t − t ′ plus the distribution Θ(t − t ′ ) in order to write it derivative and that G Ret,x (0) = 0. By writing the first term as a integral: then, as for any function f we have that ∂ t f (t − t ′ ) = −∂ t ′ f (t − t ′ ) and the fact that the derivative of the Dirac delta function is an odd function, the condition can be written in general in its vectorial form as: where the permittivity tensor for the inhomogeneous and anisotropic material is given by: Ret,x (t − t ′ ) .
(51) It is worth noting that since in this basis the tensor is diagonal, it turns out that it is expressed in the Fresnel's principal axes basis. In fact, it turns out that any point of space having material we have the same Fresnel's basis. This is why, from the beginning, we considered the same three directions of oscillation for each polarization degree of freedom. Particularly, this also holds for disjoint material bodies, although in general is clear that is possible to have bodies with different Fresnel's bases. All these cases can be considered introducing changes of basis in order to relate them. However, in order to keep simplicity in the expressions these complications will be omitted in this work. Moreover, the residual condition expressed in Eq.(50) is, for the case of inhomogeneous and anisotropic material, close to conditions considered in the Literature. For example, on one hand, is close to the one considered in Ref. [23] (for dissipative isotropic materials) as a complementary condition for the temporal gauge, but for the case of anisotropic materials.
On the other hand, the present condition is also close to the generalized Coulomb gauge condition considered in Ref. [24]. For the case of non-dissipative and also nondispersive isotropic media, i. e., for the case of constant dielectric permittivity, a close condition can be obtained. First, isotropy implies that superscripts [m] are omitted. Now, for an arbitrary type of bath, the well-know QBM theory gives that Laplace transform of the retarded Green function for the material is given by (see Refs. [5,8,10]): Then, the constant dielectric permittivity case is given by setting s = 0 in the Laplace transform (Ref. [10]): By Mellin's transform, the associated retarded Green function results: All in all, the residual condition in Eq.(50) easily reads: where the inhomogeneous constant permittivity function results ε(x) = 1 + , which results clearly ensured by taking the generalized Coulomb gauge condition of Ref. [24] (∇·[ε(x) A(x, t)] = 0) for the case of inhomogeneous materials, since in the temporal gauge E = −∂ 0 A. Is clear that in this case, without dissipation, the permittivity function in the complex s−plane is real inside the material. Therefore, the refractive index n x = ε(x) is real in each point. Due to isotropy, the Fresnel ellipsoid's picture in each point is trivial and corresponds, in every point, to a sphere because all the ellipsoid axis are equal.
However, without assuming isotropy, a similar condition for anisotropic materials can be obtained. It is straightforward that for the anisotropic case, we would obtain the same equation as Eq.(55) by replacing the permittivity function ε(x) with the inhomogeneous constant permittivity tensor ← → ε (x) = I 1 + . As in the last case, since the material has no dissipation, the permittivity tensor inside the material in the complex s − plane is real. Therefore, the refractive indexes in each direction are n Then, in each point, we define the typical (real) Fresnel's ellipsoid.
We can conclude that, for the temporal gauge, the equation of motion for µ = 0 reduces to a residual condition given, for the general case, by Eq.(50). On the other hand, if we take the remaining equations of motion (µ = m), by imposing the temporal gauge, we clearly have: In this case, considering the components of the EM dissipation kernel from Eq. (18), it is straightforward that the equation reads: Ret, Considering twice the derivative property of a product between a function and a distribution plus the initial conditions for G [m] Ret,x , we finally obtain: Ret,x . Again, from the fact that we could write this diagonal tensor associated to the retarded Green functions, it is clear that the basis that we have chosen is the Fresnel's principal axes basis (in general, the tensor would be non-diagonal). It is also remarkable the appearance of the third term, which constitutes a finite renormalization position-dependent mass term for the EM field as the one found in the scalar case in Ref. [10].
In fact, Eq.(59), in some sense, can be considered as the electromagnetic or vectorial generalization of the equation of motion for the scalar field found in Ref. [10], including all the properties related to dissipation and inhomogeneity and also to anisotropy, which is a property entirely related to the vectorial nature of the EM field. However, the equations are not formally the same, because for the scalar case the second term of the l.h.s. of Eq.(59) is a Laplacian, while in the present case there is one more term related to the divergence of the field A.
In Ref. [25], a similar model for the interaction between the matter and the EM field is considered. For simplicity, an unidimensional problem is taken with an EM field in the Coulomb gauge from the start. The field equation is deduced by solving the Heisenberg equations of motion for the matter's degrees of freedom and a Laplacian is obtained in place of the second term of l.h.s. of Eq.(59) due to the Coulomb gauge condition. This is commented in the beginning in order to explain the context in which this toy model is inspired, but no formal treatment related to gauge invariance is given. This way, the field equation has obviously the same form as for the scalar field in Ref. [10]. However, the crucial point here is that the Coulomb gauge, unlike to what happens for the free field case, may not imply A 0 = 0. In fact, in that gauge, the four components of the EM field may still be coupled due to the term involving the EM dissipation kernel and the temporal component A 0 cannot be discarded as in Ref. [25]. On the other hand, by taking the temporal gauge, Eq.(50) is the residual condition that must be satisfied, which is not necessarily Coulomb condition for the case where there are variations in material (see next Section for the homogeneous case) and thus, the equation of motion is given by Eq.(59) differing from the one considered in Ref. [25]. All in all, we can say that the toy model of Ref. [25] shows not to represent one of the components of the EM field in the context commented (with spatial variations of the dielectric properties) and it is closer to the scalar field model considered in Ref. [10], which fully coincides for the one dimensional case.
Nevertheless, the temporal gauge shows to be adequate for the interaction with matter since it decouples the components of the EM field. Indeed, in the present model, a realistic EM field (with gauge and vectorial properties) interacting with matter must satisfies Eq.(59) and the residual condition given in general by Eq.(50).
Given the equation of motion, the 'temporal-gauged' EM retarded Green tensor ← → G Ret (x, x ′ , t) can be defined as (omitting the subscripts LG for Landau gauge because, from here and so on, they are useless since the vanishing temporal gauge was already introduced): With this analysis of the EM retarded Green tensor, we can study the time evolution of the Hadamard propagator of Eq.(39) (is clear that to study the energy, Poynting vector or Maxwell tensor in this case, we should consider it expressions in a material region, what we avoid specifically in Sec.V because the arbitrariness in their definitions). For the first term of the Hadamard propagator, by introducing the Fourier transform of the retarded Green tensor in the homogeneous solutions, we can easily write: where A j 0 (k, t) by direct construction is given by: where we explicitly write the integrals over x ′ since the notation ...

ΣAi,ΣΠi
implies functional integrations over The key point here is that, beyond the chosen initial state for the EM field, for the long-time limit (t i → −∞) we have that G jl Ret (k, t−t i ) → 0 and therefore A j 0 (x 1 ) A m 0 (x 2 ) ΣAi,ΣΠi −→ 0. In other words, the term associated to the initial conditions does not contribute to the steady regime.
On the other hand, the second term of Eq.(39), associated to the material's contribution present a simpler behaviour independent of the EM retarded Green tensor's properties. In this case, the initial time t i does not appears in the EM retarded Green tensors contained in the term, but it does in one of the parts of the EM noise kernel given in this gauge by ∂ 2 tt ′ N. As we have N = N B + N P , when each term is given by Eqs. (16) and (22) respectively, the former does not contain t i while the latter does. In fact, from its definition it is easy to see that only N P goes to 0 in the long-time limit (t i → −∞). Then, N → N B in the long-time limit, proving that the Hadamard propagator in the steady regime reduces to the baths' contribution: which is an analog result of the one obtained for the scalar field in Ref. [10] but for the EM field, and is also in agreement with the steady functional approach employed in Ref. [26].
This result is indeed physically expected because the dissipative dynamics of the EM field in every point of space. It can be also easily shown (see Ref. [10] for the scalar case) that if a non-dissipative material is considered, all the material dynamics should be erased by setting N ≡ 0 and the Hadamard propagator would be defined only by the initial conditions' term, which in this case does not vanishes. Moreover, as it is mentioned in Ref. [26] and shown for an unidimensional scalar field in Ref. [10], a caveat may be done around a scenario presenting vacuum regions (or at least, regions where the field has no damping dynamics), where more than one of the terms define the steady state. This will be the objective of forthcoming work [13].

Conclusions and Forthcoming Work
In this paper we have mainly extended a previous work [10] in order to calculate the CTP generating functional for the EM field in interaction with inhomogeneous anisotropic matter, through the open system framework. We have calculated a general expression for the EM field's influence action from the interaction of the field with a composite environment, consisting in quantum polarization degrees of freedom in each point of space, and connected to thermal baths (with arbitrary temperatures). Then, we have evaluated the CTP-EM-field generating functional in the temporal gauge by implementing the Faddeev-Popov procedure. Special care has been taken about how the gauge invariance must be treated in the CTP formalism when the EM field is interacting with inhomogeneous anisotropic matter.
In addition to previous works, we have also found closed expressions for the EM-energy, the Poynting vector and the Maxwell tensor, in vacuum regions, in terms of the Hadamard propagator, showing that all of these quantities present contributions from the field's initial conditions and also from the matter degrees of freedom in the material bodies.
We then study the open electrodynamics in the temporal gauge, obtaining the EM field's equation and a residual condition closely related to the gauge condition considered in Ref. [24]. Finally, we analyzed the dynamics and steady regime of the Hadamard propagator for the case of an infinite homogeneous and isotropic material interacting with the EM field, showing that in this case, in the long-time limit, the only contribution that survives is the one associated to the bath due to the damping dynamics of the EM field in every point of space.
Nevertheless, despite on the result obtained for this last situation, a few general considerations are in order about the contributions that we found. In the transient evolution, all the contributions will be present. In the long time regime, however, we expect that the baths' and the field's initial state contributions were the only that contribute at most depending if there are (or not) vacuum regions (or regions where the field has no damping dynamics) in the problem, as it happens for the scalar field in Ref. [10]. If there is such a dependence with the initial state contribution, it means that the steady situation is non unique. This is clearly caused by the fact that the field fluctuates freely and without damping in the vacuum regions, making the contribution of the initial conditions to reach the long time regime. Following this train of thought, we have also discussed in detail several implications of our results in relation to the nonequilibrium calculation of energies and forces in Casimir physics. As a forthcoming work is pending the study of the Casimir-Lifshitz problem in a fully nonequilibrium situation [13], exploiting all the features of the present quantum approach.