In this example we compare the results of our model with the analytical solution for a plane electromagnetic wave propagating in magnetic dielectric medium. The wave is linearly polarised and propagates in the direction of the constant applied magnetic field H (aligned with axis y). The electric conductivity of the medium  is zero and the relative permittivity
 is zero and the relative permittivity  is constant.  The analytical solution for the electric field as a function of spatial variable can be described as follows [1]:
 is constant.  The analytical solution for the electric field as a function of spatial variable can be described as follows [1]:
 (1)    ![Rendered by QuickLaTeX.com \begin{equation*} $\overline{E} = E_0 \Large{ \textnormal{[}}  \hat{z} \normalsize \textnormal{cos} \huge \textnormal{(}   \dfrac{\beta _+ - \beta _-}{2} \huge \textnormal{)}  y - \normalsize \hat{x} \textnormal{sin} \huge{(}\dfrac{\beta _+ - \beta _-}{2} \huge \textnormal{)}y \huge{ \textnormal{]}} e^{-i (\beta _+ + \beta _-) y/2} $ \end{equation*}](/wp-content/ql-cache/quicklatex.com-33cf1ce3713861467ae3834ae3b38c0f_l3.png)
where  are the propagation constants for right hand circular (subscript ‘+’) and left hand circular (subscript ‘-‘) modes, with
 are the propagation constants for right hand circular (subscript ‘+’) and left hand circular (subscript ‘-‘) modes, with  and
 and  being the diagonal and off-diagonal (gyrotropic) components of the permeability tensor [Pozar]. Due to different frequency dependence of the two constants, the polarisation of the electric field vector becomes a function of the spatial variable.  This is especially seen close to the resonance frequency
 being the diagonal and off-diagonal (gyrotropic) components of the permeability tensor [Pozar]. Due to different frequency dependence of the two constants, the polarisation of the electric field vector becomes a function of the spatial variable.  This is especially seen close to the resonance frequency  , where
, where  is the applied magnetic field.  Depending on the damping coefficient
 is the applied magnetic field.  Depending on the damping coefficient  ,
,  can significantly overcome
  can significantly overcome  leading to faster rotation of the polarisation. The angle of rotation, as a function of the spatial variable
 leading to faster rotation of the polarisation. The angle of rotation, as a function of the spatial variable  can be then determined from:
 can be then determined from:
 (2)    
For frequencies below  the sign of
 the sign of  is negative, indicating that the rotation of polarisation is counter-clockwise (the same as for the right-hand circular mode). For higher frequencies
 is negative, indicating that the rotation of polarisation is counter-clockwise (the same as for the right-hand circular mode). For higher frequencies  changes the sign to positive. In the special case of
 changes the sign to positive. In the special case of  and if no losses are present
 and if no losses are present  becomes purely imaginary, so only left circular mode can propagate in the media. If
 becomes purely imaginary, so only left circular mode can propagate in the media. If  is non-zero then RHCP mode can be present too, but only as an evanescent wave.
 is non-zero then RHCP mode can be present too, but only as an evanescent wave.
Results and Discussion
Figure 1 shows a field map of electric field component  produced by MaxLLG solver for a case of a magneto-dielectric medium with parameters:
 produced by MaxLLG solver for a case of a magneto-dielectric medium with parameters:  kOe,
 kOe,  kOe,
 kOe,  GHz,
 GHz,  ,
,  , cell dimensions
, cell dimensions  mm in all directions. The solution is produced for a case with Periodic Boundary Conditions (PBC) in transverse directions to propagation. The excitation of the electric field is achieved in the initial x-z plane (at pixel
 mm in all directions. The solution is produced for a case with Periodic Boundary Conditions (PBC) in transverse directions to propagation. The excitation of the electric field is achieved in the initial x-z plane (at pixel  ) with a uniform continuous source with the above frequency.
) with a uniform continuous source with the above frequency.

Fig. 1. a) Map of  component produced by MaxLLG (time step
 component produced by MaxLLG (time step  from x-y plane in the middle of the grid. Blue arrows show the position of the PBC boundaries. All fields within PBC boundaries are the same, as should be for a plane wave.
 from x-y plane in the middle of the grid. Blue arrows show the position of the PBC boundaries. All fields within PBC boundaries are the same, as should be for a plane wave.

Fig. 2. Electric field vector as function of propagation distance  for excitation frequency of
 for excitation frequency of  GHz (and other parameters described above). The field is linearly polarised but rotates with a constant phase according to equations 1 & 2. The inset shows the dependence of propagation constants
 GHz (and other parameters described above). The field is linearly polarised but rotates with a constant phase according to equations 1 & 2. The inset shows the dependence of propagation constants   and
 and  as function of the frequency.
 as function of the frequency.
Figure 2 demonstrate the rotation of the vector  as a function of propagation distance
 as a function of propagation distance  . The wave is initially linearly polarised in
. The wave is initially linearly polarised in  , but changes its polarisation with a constant phase
, but changes its polarisation with a constant phase  along the propagation direction.  It should be noted that the excitation frequency in this case is set much below the resonance condition (
 along the propagation direction.  It should be noted that the excitation frequency in this case is set much below the resonance condition ( GHz, at
 GHz, at  kOe), so the difference
 kOe), so the difference  is not significant (see figure #, inset) and the rotation of polarisation is not significant, keeping the electric vector linearly polarised within the cycles. However if the damping parameter is increased or the excitation frequency becomes close to resonance, the RHCP mode becomes damped leading to larger ellipticity and eventually only the LHCP mode remaining with circular polarisation (see figure 3).
 is not significant (see figure #, inset) and the rotation of polarisation is not significant, keeping the electric vector linearly polarised within the cycles. However if the damping parameter is increased or the excitation frequency becomes close to resonance, the RHCP mode becomes damped leading to larger ellipticity and eventually only the LHCP mode remaining with circular polarisation (see figure 3).

Fig. 3. Electric field vector as function of propagation distance  for excitation frequency of
 for excitation frequency of  GHz. The field vectors are not aligned within the cycle and rotate significantly because the of the larger difference between propagation constants
 GHz. The field vectors are not aligned within the cycle and rotate significantly because the of the larger difference between propagation constants  and
 and  .
.

Fig. 4. Electric field intensity  as function of propagation distance
 as function of propagation distance  as calculated using MaxLLG (red dots) and analytically derived using equation 1.
 as calculated using MaxLLG (red dots) and analytically derived using equation 1.
Figure 4 shows a comparison between the solution produced by MaxLLG solver (red dots) and the analytical solution calculated with equation 1. It is evident that within the first cycles of oscillations along the propagation direction the match is nearly ideal, indicating the validity of the modelling approach as well as the precision of the numerical solution. However the further away from the excitation centre the solutions deviate as the numerical solution is affected by the transient effects and the damping both of which are dependent of the frequency and the damping parameter alpha as discussed above.
[1] Pozar, David M. Microwave Engineering. Hoboken, NJ :Wiley, 2012.

