{'fraction_non_alphanumeric': 0.03298009309695949, 'fraction_numerical': 0.016143408933346537, 'mean_word_length': 4.644214645053102, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Three Dimensional Lattice-Boltzmann Model for Electrodynamics

In this paper we introduce a novel 3D Lattice-Boltzmann model that recovers in the continuous limit the Maxwell equations in materials. In order to build conservation equations with antisymmetric tensors, like the Faraday law, the model assigns four auxiliary vectors to each velocity vector. These auxiliary vectors, when combined with the distribution functions, give the electromagnetic fields. The evolution is driven by the usual BGK collision rule, but with a different form for the equilibrium distribution functions. This LBGK model allows us to consider for both dielectrics and conductors with realistic parameters, and therefore it is adequate to simulate the most diverse electromagnetic problems, like the propagation of electromagnetic waves (both in dielectric media and in waveguides), the skin effect, the radiation pattern of a small dipole antenna and the natural frequencies of a resonant cavity, all with 2% accuracy. Actually, it shows to be one order of magnitude faster than the original FDTD formulation by Yee to reach the same accuracy. It is, therefore, a valuable alternative to simulate electromagnetic fields and opens lattice Boltzmann for a broad spectrum of new applications in electrodynamics. In this paper we develop a 3D Lattice-Boltzmann model that recovers in the continuous limit the Maxwell equations for macroscopic mediums. The model can sucessfully reproduces the propagation of the electromagnetic waves in dielectric mediums and waveguide, the skin effect, the electrical dipole radiation and the electromagnetic response of a resonant cavity.


I. INTRODUCTION
To simulate the electromagnetic fields inside materials is a fundamental tool in optics and electrodynamics, specially when boundary conditions are so complex that an analytical solution is impossible. Some examples of this phenomena are the design of antennas and the study of electrical discharges across media, among many others. For this reason, numerous numerical methods for the Maxwell equations have been implemented. One of these methods is the FDTD (Finite-difference time-domain) [1], which is based on finite differences and works pretty well, but consumes huge ammounts of computational resource. Another one is the ray tracing method [2], that is only useful when the electromagnetic signals are of a single frequency. In addition, one of the actual problems in computational electrodynamics is the simulation of electromagnetic pulses in heterogenous media, that is one of the main advantages of lattice-Boltzmann methods.
One of the numerical methods for simulating fluids is Lattice Boltzmann (LB) [3], which was developed from lattice-gas automata. Lattice Boltzmann simulations are performed on regular grids of many cells and a small number of velocity vectors per cell, each one associated to a set of density distribution functions, which evolve and spread together to the neighbohr cells according to the collisional Boltzmann equation. As far as we know, there is still not a lattice-Boltzmann model for electrodynamics, but the effects of Maxwell equations on plasmas has been included (in the form of the diffusion equation for the magnetic field) as part of LB models for magnetohydrodynamics. The first LB model for studying plasmas reproduces the resistive magnetohydrodynamic equations and was developed by Chen [4,5] as an extension of the Lattice-Gas model developed by Chen and Matthaeus [6] and Chen, Matthaeus and Klein [7]. This LB model uses 37 velocity vectors per cell on a square lattice and is developed for two dimensions. Thereafter, Martinez, Chen and Matthaeuss [8] decreased the number of velocity vec-tors from 37 to 13, which made easier a future 3D extension. Some of the first LB models for magnetohydrodynamics in 3D were developed by Bryan R. Osborn [9] and Breyiannis and Valougeorgi [10]. They used 19 vectors on a cubic lattice for the fluid, plus 7 vectors for the magnetic field, which makes a total number of 26 vectors per cell. By following a different path, Fogaccia, Benzi and Romanelli [11] introduced a 3D LB model for simulating turbulent plasmas in the electrostatic limit. All these models reproduce the resistive magnetohydrodynamc equations for a single fluid. Nevertheless, in a recent work [12], we have introduced a 3D LB model that reproduces the two fluids theory for plasmas. This model uses a cubic lattice and a set of 39 velocity vectors and includes explicitly the Maxwell equations in the vacuum.
In this paper, we introduce a 3D lattice-Boltzmann model that recovers the Maxwell equations in material media. The model is based on our previous develop for magnetohydrodynamics [12], but it is able to reproduce the behavior of electromagnetic fields inside dielectrics, magnets and conductors. The model uses 33 vectors per cell and 50 probability density functions, and performs well in a wide range of tests. Section II describes the model, with the evolution rules and the equilibrium expressions for the 50 density functions, plus the procedure to compute both the charge and current densities and the electric and magnetic fields. The Chapman-Enskog expansion showing how these rules recover the electrodynamic equations is developed in Appendix A. In order to validate the model, we simulate in section III the reflection of electromagnetic pulses on dielectric media, the propagation of electromagnetic waves along waveguide, the skin effect on conductors, the radiation pattern of an oscillating electrical dipole and the natural frequencies of a resonant cavity. The main results and conclusions are summarized in section IV. neighbors, and each vector is associated to a distribution function f i . The distribution function evolves according to the Boltzmann equation, where Ω i ( x, t) is a collision term, which is usually taken as a time relaxation to some equilibrium density, f eq i . This is known as the the Bhatnagar-Gross-Krook (BGK) operator [13], where τ is the relaxation time and f eq i ( x, t) is the equilibrium function. The equilibrium function is chosen in such a way, that (in the continuum limit) the model simulates the actual physics of the system. For our 3D model, we use a cubic regular grid of lattice constant δx= √ 2cδt, with c the light speed in vaccum (c≃3 × 10 8 m/s); in other words, c=1/ √ 2 in normalized lattice units (time unit = δt, spatial unit = δx). There are 13 velocity vectors (figure 1), 13 different vectors for This makes 12 vectors. The missing one is the rest vector v 0 , with components (0, 0, 0). Associated to each velocity vector v p i there are two electric vectors e p ij and two magnetic vectors b p ij (j=0, 1), which are used to compute the electromagnetic fields (figure 3). The electric vectors are perpendicular to v p i and lie on the same plane p. The magnetic vectors are perpendicular to v p i , too, but they are also perpendicular to the plane p. In terms of the velocity vectors (3), they are: with the rest vector e 0 =(0, 0, 0), and  The macroscopic variables are computed as follows: where D, E and J are subsidiary fields that represent the displacement field, the electric field and the total current density before external forcements, respectively (the actual mean fields, including external forcements, are described below). In addition, B is the induction field, H is the magnetic field, ρ c is the total charge density and ǫ r , µ r and σ are the dielectric constant, the permeability constant and the conductivity of the medium, repectively. For the evolution of these distribution functions we follow the proponsal of Zhaoli Guo, Chuguang Zheng and Baochang Shi [14] that includes external forcements in the lattice Boltzmann equations, as follows: The mean electric field, E ′ , is defined by and the mean density current vector J ′ is defined by In order to have the mean density vector J ′ in terms of the subsidiary fields, we replace Eq.(10) into Eq.(11) to obtain Next, we adopt BGK collision terms Ω p(r) ij and Ω where the relaxation time is set to τ = 1 2 . The equilibrium functions for the electromagnetic fields are This completes the definition of the lattice Boltzmann model. The proof that this lattice-Boltzmann model, via a Chapman-Enskog expansion, recovers the Maxwell equations is shown in Appendix A. This model reproduces the following equations: It is known [15] that we can derive from these equations the following expresions for the divergence of the electromagnetic fields: that is, the other two Maxwell equations hold if they are satisfied as initial conditions.

A. Dielectric Interface
First, we simulate the propagation of a electromagnetic gaussian pulse crossing an dielectric interface. For this purpose, we took an unidimensional array of L cells with periodic boundary conditions in the z coordinate and with each cell being its own neighbor in both directions of the x and y coordinates. One half of the simulation space, z<L/2, is vaccum (ǫ=1.0) and the other half, z>L/2, represents a dielectric medium with dielectric constant ǫ=2.5. For avoiding abrupts changes of the dielectric constant between two neigboring cells we choose the following distribution of the permitivity: The functional form of the incident gaussian electromagnetic pulse is given by the equations centered at z 0 . The constant α fixes the pulse width, E 0 is the pulse amplitude and the constant B 0 is related with E 0 according to the relation E 0 =cB 0 , with c the vacuum light speed. For the simulation we choose L=200, c=1/ √ 2, E 0 =0.001, α=0.01 y z 0 =40 (in normalizated units). The initial condition is shown in figure 4, and the electric field after 140 time steps is shown in figure 5.
The theoretical predictions of the amplitudes of the transmited and reflected pulses can be calculated from [15] For our case these predictions give a ratio between the transmited amplitude and the incident amplitude of 0.7751 and a ratio between the reflected and incident amplitudes of 0.2249. The values measured in the simulation are 0.7750 and 0.2249, respectively. This results correspond to errors smaller than 1% with the theoretical results. The results of this simulation show us that the model is able to reproduce the propagation of electromagnetic fields when dieletric media are present. This simulation takes less than 1 minute in a single Pentium IV de 3.0 GHz.

B. Skin Effect
For reproducing the skin effect, we simulate an electromagnetic plane wave that arrives perpendicular to a conductor material. For this purpose we construct a unidimensional space of L cells as before, with zero conductivity for z<L/4 and σ 0 conductivity for z>L/4 is σ 0 . In order to avoid abrupt changes between neigbohring cells we choose the following distribution for the conductivity: In order to generate an incoming plane wave we set at z=0 an electric oscillator by imposing there an electromagnetic field given by where E 0 is the wave amplitude, cB 0 =E 0 and ω is the angular frequency of the wave. Some simulation snapshots where E 0 is the amplitude of the electric field at the outside and ∆ is known as the skin thickness. For good conductors this thickness is given by [15] ∆ = 2 σµ 0 µ r ω . oscillation follows in excellent agreement the theoretical prediction. This simulation takes less than 1 minute in a single Pentium IV de 3.0 GHz.

C. Electric Dipole
To simulate an electric dipole we construct an array of 100 × 100 × 100 cells with free boundary conditions (each limit cell takes itself as his own missing neighbohr). In the center of this array we insert a small oscillating current density in the z direction, where J ′ 0 is the amplitude of the current density and T is the oscillation. For avoiding any abrupt change of the physical quantities between two neigboring cells we choose actually a gaussian functional form for the amplitude of the current density J ′ 0 , The periode was set to T =25.0 time steps and the amplitude to J 0 =0.0001 (in automaton units).
The simulation results are shown in Figures 9 and 10. Figure 9 draws the equipotential lines after 56 time steps. Figure 10 shows the amplitude of the radiated magnetic field as a function of the x component and the theoretical values given by [15] Here, k is the magnitude of the wave vector, k= ω c , with c the light speed in vaccum and ω the oscillation angular frequency. The magnitude P of the electric dipolar   momentum is computed as where Σ is the effective volume of the dipole. Table I shows the simulation peaks and the theoretical values for the magnetic field amplitude. Again, the simulation matches the theoretical predictions with an accuracy of less than 2%. This illustrates the possibilities of the method to simulate antennas. This simulation takes 40 minutes in a single Pentium IV de 3.0 GHz.

D. Waveguide. Microstrip
In this section we simulate a microstrip waveguide (see figure 11). For this purpose, we chose a 3D array of 100×50×50 cells with free boundary conditions and we insert two parallel metallic layers of conductivity σ, one wider than the other. The dimensions, depicted in figure  11, are (in normalizated units): t=5, d=20, W =50 and ǫ r =1. The signal enters into the waveguide by forcing the cells at y=0 to have the electric and magnetic fields of a plane wave (Eqs. (21)) with E 0 =0.01, c=1/ √ 2 and ω=π/100 (see figure 13). The conductivity is smoothly changed along three cells, as in Sec. (see figure 12). 4d 4d where With this equation we obtain a theoretical value of Z=72.6Ω. This gives us a 3% difference between the simulation and the theoretical predicition. This simulation takes 90 minutes in a single Pentium IV de 3.0 GHz.

E. Resonant Cavity
As a last test, we simulate a cubic resonant cavity and found its resonance frecuencies. The cavity we chose is an array of 50×50×50 cells, corresponding to a size of  5×5×5 cm, with periodic boundary conditions, but imposing a null electric field at the boundary (a perfect conductor). Inside the cavity we set an emmiter point as a single cell with an oscillating electric field in the x direction (the antenna) and a receptor point (the measurement instrument) as a sigle cell where we measure the electric field amplitude. The emiter point was set at (5, 5, 5) (in cell units) and the receptor point at (5,5,45). The oscillation frequency was scanned from 0.0182 to 0.027 oscillations per time step by multipying each previous value by a constant factor of 0.0072 (a total of 22 different frequencies). Each frequency starts a single run, consisting of 5000 time steps before equilibrium and four whole oscillations to estimate the amplitude of the electric field. Figure 16 shows amplitude of the electric field as a function of frequency. The resonant peaks are cleary identified. Table II compares the values these resonant frequencies with the theoretical predictions com-   puted from [15] ω nmp = c ǫ r µ r All differences are smaller than 1%. The whole simulation (22 runs of more than 5000 time steps each on 125000 cells) takes 18 days in a single Pentium IV of 3.0 GHz (actually, we ran the simulation on a cluster of five nodes on two days).

IV. CONCLUSION
In this paper we introduce a three dimensional lattice-Boltzmann model to simulate the Maxwell equations in macroscopic media (dielectric, magnetic and contuctor materials). The model actually reproduces the Maxwell equations in the continuous limit, as prooved in Appendix A. In order to test the model, we simulated some typical electromagnetic sets, namely: the amplitude of the transmitted and reflected pulses at a dielectric interface, the exponential decay in the amplitude of a incident plane wave on a conductor (skin effect), the amplitude (as a function of distance) of the magnetic field radiated by an electric dipole, the characteristic impedance of a microstrip waveguide and the resonant frequencies of a cubic resonant cavity. In all cases, the differences between the simulational results and the theoretical predictions lay between 1% and 3%. These five simulations proove that our LB model works very well in the more different situations.
The model let us to select the electromagnetic constants (magnetic permeability, electric permitivity and conductivity) for each cell at pleasure. Nevertheless, for avoiding numerical instabilities on the interfaces between two different materials we need to smooth the transition across three cells, approximately. For this purpose we have used tanh functions, but other functional forms can be also employed. This is a standard procedure in many numerical models [1,2].
The main advantage of our LB model is that it is about ten times faster than FDTD based models. Indeed, a the simulation with FTDT of the propagation of an electomagnetic pulse produced by a discharge antenna on a grid of tantos × tantas × tantas cells takes around Xh in a Pentium IV d. In comparison, our simulation of a Nos de que on a nose geu grid took xh. For this reason, the LB model seems to be promising in a variety of applications, including among others: the propagation of electromagnetic pulses produced by atmospheric rays or discharge antennas, the light diffraction in objects that have a complex geometry or the sputtering of electromagnetic signals across buidings. The velocity and precision of our lattice-Boltzmann model for electrodynamics are excellent for the design and optimization of antennas, which requires to simulate a large number of configurations, and very specially in the design of pulse antennas, where single-frequency numerical methods fail. All these applications can be thema of future works.
The developments in this paper also show how to introduce antisymmetrical tensors in conservation equations in Lattice-Boltzmann models. These equations are useful to simulate many interesting phenomena, like seismic or torsion waves. These are also topics of future development.
Hereby we have introduced a Lattice-Boltzmann model for electrodynamics that actually reproduces the Maxwell equations in media, with a plenty of future aplications. Moreover, this model shows how to construct LB models that fulfills conservation laws with antisymmetric tensors. We hope that this valuable theoretical develop will push the evolution of LB models further away in the horizon of even more exciting applications.

APPENDIX A: CHAPMAN-ENSKOG EXPANSION
The Boltzmann equations for each fluid, Eq. (7) and (8), determine the system evolution. This evolution rule gives in the continuum limit the macroscopic differential equations that the system satisfies. This is known as the Chapman-Enskog expansion.
To develop the Chapman-Enskog expansion, we start by taking the Taylor expansion of these equations until second order in spatial and temporal variables, where α, β=x, y, z denote the x, y and z component. Next, we expand the distribution functions, the spatial and temporal derivates in a power series of a small parameter, λ, It is assumed that only the 0th order terms in ǫ of the distribution functions contribute to the macroscopic variables. So, for n > 0 we have i,j,p G p(r)(n) ij e p ij = 0 , i,j,p G p(r)(n) ij b p ij = 0 .
The main current density J ′ is of the order λ [14], so we can write J ′ =λ J ′ 1 . Because G p(r)eq ij is now function of J ′ , we need to develop a Chapman-Enskog expansion of the equilibrium function, too: G p(r)eq ij = G p(r)(0)eq ij + λG p(r)(1)eq ij + λ 2 G p(r)(2)eq ij . (A7) So, if we replace these results into Eqs.(A1) and (A2), we obtain for the zeroth order in λ G p(r)(0)eq ij = G p(r)(0) ij . (A8a) For the first order in λ of the distribution function we have and for the second order in λ we obtain ) .
The first order and the second order terms for the equilibrium function of the electromagnetic fields are obtained by replacing the Eq. (10) into Eq. (14). From this, we can obtain, grouping in λ powers, and because the continuity equation Eq.(A14), we arrive to Finally, we obtain So, if the initial conditions of the electromagnetic fields satisfies the Maxwell equations ∇ · B = 0 . (A25) these equations will be reproduced correctly for all time.