## 6 Numerical Methods
In cylindrical coordinates, we end up with the scheme at the grid point \((r_{i},z_{j})\):  
\[(B_{\varphi})_{i,j}^{k+1}=(B_{\varphi})_{i,j}^{k}+\frac{\Delta t}{2h}\Big{(}- (E_{r})_{i,j+1}^{k}+(E_{r})_{i,j-1}^{k}+(E_{z})_{i+1,j}^{k}-(E_{z})_{i-1,j}^{k }\Big{)}\]  
\[(E_{r})_{i,j}^{k+1}=(E_{r})_{i,j}^{k}+\frac{c^{2}\Delta t}{2h}\Big{(}-(B_{ \varphi})_{i,j+1}^{k}+(B_{\varphi})_{i,j-1}^{k}\Big{)}-\Delta t\rho_{i,j}^{k}V _{r}\]  
\[(E_{z})_{i,j}^{k+1}=(E_{z})_{i,j}^{k}+\frac{c^{2}\Delta t}{2h}\Big{(}(B_{ \varphi})_{i+1,j}^{k}-(B_{\varphi})_{i-1,j}^{k}\Big{)}+\frac{c^{2}\Delta t}{r _{i}}(B_{\varphi})_{i,j}^{k}-\Delta t\rho_{i,j}^{k}V_{z} \tag{28}\]  
For the computation of the divergence of the electric field, we used the following approximation of (13), valid in cylindrical coordinates:  
\[\rho_{i,j}^{k}=\frac{(E_{r})_{i,j}^{k}}{r_{i}}+\frac{(E_{r})_{i+1,j}^{k}-(E_{ r})_{i-1,j}^{k}}{2h}+\frac{(E_{z})_{i,j+1}^{k}-(E_{z})_{i,j-1}^{k}}{2h} \tag{29}\]  
To ensure stability we must require a Courant-Friedrichs-Lewy (CFL) type condition--i.e., that the time-step \(\Delta t\) is proportional to the space discretization parameter \(h\). Such a condition is not particularly restrictive for our purposes. Equations (9) and (10) contain the time derivatives \(\partial p/\partial t\), \(\partial V/\partial t\). It is quite simple to come out with appropriate numerical schemes for their discretization. However, in this paper we are mainly concerned with the case of free-waves (see Section 5), so that \(p\) and the total derivative \(DV/Dt\) will be set to zero at any time. Since, in practice, \(\mathbf{E}\) and \(\mathbf{V}\) will stay orthogonal, Equation (9) becomes trivial, whereas (10) is substituted by (20). In the particular setting we are interested in, the vectors of the triplet \((\mathbf{E},\mathbf{B},\mathbf{V})\) will constitute an orthogonal frame and satisfy \(|\mathbf{E}|=|c\mathbf{B}|\) and \(|\mathbf{V}|=c\). What we just wrote is always true with the exception of some transition regions that will be described later on. Let us notice that the approximation of (10) in its general form is far from being trivial. First of all, there is no viscous term, a property that makes the theoretical analysis and the numerical implementation quite hard. Secondly, as already pointed out in the examples of Section 5, we also have \(\
abla\cdot\mathbf{V}\
eq 0\). This means that many numerical schemes available for the discretization of the equations of fluid dynamics may turn out to be unsuitable. Finally, we observe that there may be points in the computational domain where \(\rho=\
abla\cdot\mathbf{E}\) is zero (or relatively small). This is usually true in correspondence with points where \(\
abla p\) is also zero. As a consequence, the unreliability of the ratio \(\
abla p/\rho\) may give origins to severe instability problems in the treatment of (10). As mentioned above, here we will stay away from these situations; nevertheless, such problems exist and should be appropriately pondered within the general context. Later in this paper, we will study the transmission of an electromagnetic wave across the interface separating media displaying different conductive properties. We describe what happens inside each medium with the help of the equations relative to free-waves. In the instant of transition, the wave is not of free type, so that suitable corrections on \(\mathbf{V}\) are taken into account. For simplicity, we will continue however to neglect the role of the potential \(p\).