# Symmetric and antisymmetric constitutive tensors for
bi-isotropic and bi-anisotropic media

Pedro D. S. Silva<EMAIL_ADDRESS><EMAIL_ADDRESS>Rodolfo
Casana<EMAIL_ADDRESS><EMAIL_ADDRESS>Manoel M. Ferreira Jr
<EMAIL_ADDRESS><EMAIL_ADDRESS>Departamento de Física,
Universidade Federal do Maranhão, Campus Universitário do Bacanga, São Luís,
Maranhão, 65080-805, Brazil

###### Abstract

The Maxwell equations and the constitutive relations describe the classical
propagation of electromagnetic waves in continuous matter. Here, we
investigate the effects stemming from extended constitutive relations on the
propagation of waves in bi-isotropic and bi-anisotropic media using a
classical general approach based on the evaluation of dispersion relations and
refractive indices. For the bi-anisotropic media, we specify two classes of
magnetoelectric parameters represented by symmetric and antisymmetric tensors.
The three cases examined have provided real and distinct refractive indices
for two propagating modes, which implies birefringence. The propagating modes
were also carried out in all cases. The anisotropy or birefringence effect,
given by the rotatory power or phase difference, was evaluated in terms of the
magnetoelectric parameters of the theory in each case. The propagation
orthogonal to the vectors used to parametrize the symmetric and antisymmetric
magnetoelectric tensors is described by distinct modes, representing a route
to identify the kind of bi-anisotropic medium examined. The group velocity and
Poynting vector were also evaluated for all the cases examined to discuss the
energy propagation in these anisotropic media.

Electromagnetic wave propagation; Optical constants; Magneto-optical effects;
Birefringence

###### pacs:

41.20.Jb, 78.20.Ci, 78.20.Fm

##  I Introduction

As is well known, the propagation of electromagnetic waves in a continuous and
infinite medium is described by the standard Maxwell equations in the matter
Jackson ; Zangwill ,

$\displaystyle\mathbf{k}\cdot\mathbf{D}$
$\displaystyle=0\,,\quad\mathbf{k}\times\mathbf{H}+\omega\mathbf{D}=\mathbf{0},$
(1a) $\displaystyle\mathbf{k}\cdot\mathbf{B}$
$\displaystyle=0\,,\quad{\bf{k}}\times{\bf{E}}-\omega\mathbf{B}={\bf{0}},$
(1b)

where we use a typical plane wave ansatz,
$\mathbf{E}=\mathbf{E}_{0}e^{\mathrm{i}(\mathbf{k}\cdot\mathbf{r}-\omega t)}$
and $\mathbf{B}=\mathbf{B}_{0}e^{\mathrm{i}(\mathbf{k}\cdot\mathbf{r}-\omega
t)}$. The electric displacement and magnetic field, $\mathbf{D}$ and
$\mathbf{H}$, respectively, contain the electromagnetic response of the matter
in the form of electric polarization and magnetization, respectively. Besides,
these phenomena are related to the constitutive relations involving the fields
$(\mathbf{D},\mathbf{E})$ and $(\mathbf{H,B})$. For linear, isotropic, and
homogeneous dielectric matter, the constitutive relations take on the usual
form,

$\mathbf{D}=\epsilon\mathbf{E},\quad\mathbf{H}=\mu^{-1}\mathbf{B},$ (2)

where $\epsilon$ is the electric permittivity and $\mu$ is the magnetic
permeability of the medium (constant parameters, in this case), given by
$\epsilon=\epsilon_{0}(1+\chi^{E})$, $\mu=\mu_{0}(1+\chi^{M})$. Here,
$\chi^{E}$ and $\chi^{M}$ represent the electric and magnetic susceptibility,
respectively, contributing to the polarization,
$\mathbf{P}=\epsilon_{0}\chi^{E}\mathbf{E}$, and magnetization vector,
$\mathbf{M}=\chi^{M}\mathbf{H}$ Jackson ; Zangwill ; Landau . The simplest
configuration appearing in Eqs. (2) describes a medium, like water and glass,
where the physical properties do not depend on the direction of the wave
propagation.

The complexity and diversity of electromagnetic phenomena in matter are
addressed by general permittivity and permeability tensors, $\epsilon_{ij}$
and $\mu_{ij}$, written as $3\times 3$ matrices. These tensors are suitable to
describe interesting scenarios endowed with anisotropy, where the constitutive
relations (2) read

$D^{i}=\epsilon_{ij}E^{j},\quad H^{i}=(\mu^{-1})_{ij}B^{j},$ (3)

with

$\displaystyle\epsilon_{ij}$
$\displaystyle=\epsilon_{0}(\delta_{ij}+\chi_{ij}^{E}),$ (4a)
$\displaystyle\mu_{ij}$ $\displaystyle=\mu_{0}(\delta_{ij}+\chi_{ij}^{M}),$
(4b) and $\chi_{ij}^{E}$ and $\chi_{ij}^{M}$ representing the susceptibility
tensors. The expressions in Eqs. (4) include the polarization and
magnetization contributions, $P^{i}=\epsilon_{0}\chi_{ij}^{E}E^{j}$ and
$M^{i}=\chi_{ij}^{M}H^{j}$, which usually appears in the constitutive
relations as  $D^{i}=\epsilon_{0}E^{i}+P^{i},\quad
B^{i}=\mu_{0}H^{i}+\mu_{0}M^{i}.$ (5)

For anisotropic configurations, the tensor $\epsilon_{ij}$ describes uniaxial
and biaxial crystals Landau ; Bain ; Fowles ; Kurmanov , which display optical
activity (chirality) Condon  and birefringence Roth ; Hecht .

Effects of anisotropy may also appear in linear electrodynamics with linear
extended constitutive relations, envisaged as

$\displaystyle\mathbf{D}$
$\displaystyle=\hat{\epsilon}\,\mathbf{E}+\hat{\alpha}\,\mathbf{B},$ (6a)
$\displaystyle\mathbf{H}$
$\displaystyle=\hat{\beta}\,\mathbf{E}+\hat{\zeta}\,\mathbf{B},$ (6b)

where $\hat{\epsilon}=[\epsilon_{ij}]$, $\hat{\alpha}=[\alpha_{ij}]$,
$\hat{\beta}=[\beta_{ij}]$, and $\hat{\zeta}=[\zeta_{ij}]$ represent, in
principle, $3\times 3$ complex matrices. Such expressions above contain
additional magnetoelectric responses of the medium to electromagnetic fields:
$\hat{\alpha}$ measures the electric response to the magnetic field and
$\hat{\beta}$ represents the magnetic response to the electric field. In this
generalized context, if the relations (5) remain valid, the polarization and
magnetization vectors receive contributions from the magnetic and electric
sectors, respectively, being given as

$\displaystyle P^{i}$
$\displaystyle=\epsilon_{0}\chi_{ij}^{E}E^{j}+\alpha_{ij}B^{j},$ (7a)
$\displaystyle M^{i}$
$\displaystyle=\chi_{ij}^{M}H^{j}+\tilde{\beta}_{ij}E^{j},$ (7b)

where it holds that
$\hat{\zeta}=\hat{\mu}^{-1},\quad\hat{\beta}=-\mu_{0}\hat{\mu}^{-1}\tilde{\beta}$.

In order to ensure the energy conservation of the electromagnetic field in a
medium where the constitutive relations (6a) and (6b) prevail, the Poynting
theorem leads to the following set of relations for the complex matrices:

$\displaystyle\hat{\epsilon}$ $\displaystyle=\hat{\epsilon}^{\dagger},$ (8a)
$\displaystyle\hat{\mu}^{-1}$ $\displaystyle=(\hat{\mu}^{-1})^{\dagger},$ (8b)
$\displaystyle\hat{\alpha}$ $\displaystyle=-\hat{\beta}^{\dagger}.$ (8c)

For details, see Sec. II and Refs. Sihvola ; Kong ; Kamenetskii . Relation
(8c) will be crucially relevant in the analysis of the present work, as will
be clear in the next sections.

The simplest version of the relations in Eqs. (6),

$\begin{pmatrix}\mathbf{D}\\\
\mathbf{H}\end{pmatrix}=\begin{pmatrix}\epsilon&\alpha\\\
\beta&{\zeta}\end{pmatrix}\begin{pmatrix}\mathbf{E}\\\
\mathbf{B}\end{pmatrix}\,,$ (9)

includes $\epsilon$, $\alpha$, $\beta$, and $\zeta$ as single parameters and
describes the physics of bi-isotropic media (the most general linear,
homogenous and isotropic materials Sihvola ; Kong ), corresponding to the case
in which the matrices $\hat{\alpha}=[\alpha_{ij}]$, $\hat{\beta}=[\beta_{ij}]$
are diagonal and isotropic. In this case, in order to be consistent with
energy conservation, the relations (9) yield

$\alpha=-\beta^{*}.$ (10)

In the configuration the constitutive relations (9) have the form
$\mathbf{D}={\epsilon}\mathbf{E}+{\alpha}\mathbf{H}$ and
$\mathbf{B}={\beta}\mathbf{E}+{\zeta}\mathbf{H}$, the condition (10) becomes
$\alpha=\beta^{*}=\psi+i\kappa$, where $\psi$ is the Tellegen coefficient and
$\kappa$ is the chirality coefficient Aladadi . See Eqs. (212) for details.
The bi-isotropic relations (9) have been much studied in both theoretical
Sihvola1 ; Sihvola2 ; Sihvola3 ; Nieves ; Gauthier  and applied aspects Rado ;
Aladadi ; Jelinek , being also important to address optical properties Chang ;
Zou  and other properties of topological insulators Urrutia ; Urrutia2 ;
Lakhtakia ; Winder ; Li ; Li1 ; Tokura . Bi-isotropic relations are relevant
for axion electrodynamics citekey ; Sekine ; Tobar ; BorgesAxion ,
construction of optical isolators from chiral materials Silveirinha , the
Casimir effect in chiral media Casimir , and other applications Darinskii .
Furthermore, bi-anisotropic “chiral materials”, described by relations (6)
involving anisotropic tensors, were employed to investigate relativistic
electron gas Carvalho , time-dependent magnetoelectric parameters Lin , Weyl
semimetals Halterman ; Zu , magnetized materials Krupka1 ; Krupka2 , and
anisotropic dispersion relations Hillion ; Yakov ; Damaskos . It is also
worthy to mention some effects engendered by the anisotropic magnetoelectric
parameters, corresponding to the off-diagonal elements of the matrices
$\hat{\alpha}=[\alpha_{ij}]$ and $\hat{\beta}=[\beta_{ij}]$. Nondiagonal
terms, for instance, were examined in the discovery of electromagnons in
perovskites, which revealed an absorption difference of light propagating in
opposite directions (directional dichroism) Takahashi . Magnetoelectric
diagonal (and anisotropic) coefficients were investigated in the context of
multiferroic materials, where they induced a light polarization rotation angle
Kurumaji .

In extended scenarios, the constitutive tensors of relations (6) may also
depend on the space coordinates, standing for the description of
nonhomogeneous bi-isotropic and/or bi-anisotropic media Novitsky . These
tensors can present dependence on the magnitude of the electromagnetic fields
as well, $\hat{\epsilon}=\hat{\epsilon}({E,B})$ and
$\hat{\mu}=\hat{\mu}(E,B)$, a kind of approach which accounts for
birefringence in nonlinear electrodynamics Lorenci1 , allowing to recover the
Kerr and Cotton-Mouton effects in particular configurations Lorenci2 . A more
involved and general nonlinear construction, where the magnetoelectric
coefficients exhibit dependence at second order on the electromagnetic field
components, was recently examined Lorenci3 .

Generalized constitutive relations can also be envisaged for the current
density as an extension of the standard Ohm’s law. Such relations can be
written as $J^{i}=\sigma E^{i}+\sigma^{B}_{ij}B^{j}$, where $\sigma$ is the
usual Ohmic conductivity and $\sigma^{B}_{ij}$ is a general magnetic
conductivity tensor. An isotropic tensor, $\sigma^{B}_{ij}=\Sigma\delta_{ij}$,
stands for the chiral magnetic effect (CME) Kharzeev1 ; Qiu ; Fukushima .
Isotropic and anisotropic symmetric and antisymmetric conductivity tensors
were examined in Ref. Pedro1 . The antisymmetric parametrization of
$\sigma^{B}_{ij}$ also has found realization in some Weyl semimetals Kaushik .

Another possible extension occurs in the context of a Lorentz-violating
anisotropic electrodynamics Tobar1 ; Bailey , with constitutive relations
written as

$\begin{pmatrix}\mathbf{D}\\\
\mathbf{H}\end{pmatrix}=\begin{pmatrix}\epsilon\mathbb{1}+{\hat{\kappa}_{DE}}&{\hat{\kappa}_{DB}}\\\
&\\\
{\hat{\kappa}_{HE}}&{\mu}^{-1}\mathbb{1}+{\hat{\kappa}_{HB}}\end{pmatrix}\begin{pmatrix}\mathbf{E}\\\
\mathbf{B}\end{pmatrix}\,,$ (11)

where $\hat{\kappa}_{DE}$, $\hat{\kappa}_{DB}$, $\hat{\kappa}_{HE}$, and
$\hat{\kappa}_{HB}$ are dimensionless $3\times 3$ matrices composed of vacuum,
$\hat{\kappa}_{DE}^{vac}$, $\hat{\kappa}_{DB}^{vac}$,
$\hat{\kappa}_{HE}^{vac}$, and $\hat{\kappa}_{HB}^{vac},$ and matter pieces,
$\hat{\kappa}_{DE}^{matter}$, $\hat{\kappa}_{DB}^{matter}$,
$\hat{\kappa}_{HE}^{matter}$, and $\hat{\kappa}_{HB}^{matter}$. These
generalized scenarios lead to unusual electrodynamics where magnetoelectric
parameters stemming from Lorentz symmetry violation appear in matter or
vacuum, giving rise to interesting effects potentially related to the
phenomenology of new materials. A classical field theory approach to the
description of wave propagation in a continuous chiral medium supporting
higher-order derivative Lorentz-violating electrodynamics was recently
examined Pedro2 .

In this work, we investigate the possible effects stemming from extended
linear constitutive relations (6a) and (6b), assuming isotropic electric
permittivity and magnetic permeability, $\zeta_{ij}={\mu}^{-1}\delta_{ij}$,
$\epsilon_{ij}=\epsilon\delta_{ij}$, and that the tensors $\alpha_{ij}$ and
$\beta_{ij}$ may be described by symmetric and antisymmetric parametrizations.
For the three cases investigated, we have obtained general dispersive
equations which provide the refractive indices for any propagation direction.
For the anisotropic magnetoelectric tensors, we have worked out specific
solutions for the particular propagation axis in order to discuss the optical
repercussions. More specifically, the bi-anisotropic symmetric constitutive
relations are parametrized in terms of a 3-vector ${\bf{d}}$,

$\displaystyle{\bf{D}}=\epsilon{\bf{E}}+\tilde{\alpha}{\bf{d}}({\bf{d}}\cdot{\bf{B}}),\quad{\bf{H}}=\mu^{-1}{\bf{B}}+\tilde{\beta}{\bf{d}}({\bf{d}}\cdot{\bf{E}}),$
(12)

for which we discuss the dispersion relations, refractive indices, and
birefringence for special configurations where the propagation vector
direction is along and perpendicular to the vector ${\bf{d}}$. The
antisymmetric constitutive relations are parametrized in terms of two
3-vectors ${\bf{a}}$ and ${\bf{b}}$,

$\displaystyle{\bf{D}}=\epsilon{\bf{E}}+{\bf{a}}\times{\bf{B}},\quad{\bf{H}}=\mu^{-1}{\bf{B}}+{\bf{b}}\times{\bf{E}},$
(13)

satisfying ${\bf{b}}={\bf{a}}^{*}$ and used to describe the particular
scenarios where the propagation direction is longitudinal and orthogonal to
the vector ${\bf{a}}$.

The paper is outlined as follows: In Sec. II, we present the basic formalism
for obtaining the dispersion relations and refractive indices in a general
scenario of extended constitutive relations. In Sec. III, we discuss the
electromagnetic wave propagation in the bi-isotropic case. In the sequel, we
focus on the isotropic-anisotropic constitutive relations, examining symmetric
(see Sec. IV) and antisymmetric (see Sec. V) configurations for the tensors
$\alpha_{ka}$ and $\beta_{ka}$. Finally, in Sec. VI, we summarize our results.
Throughout the paper, we use natural units.

##  II  Dispersion relations for bi-isotropic and bi-anisotropic media
described by extended linear constitutive relations

In this section, we start from the Maxwell equations in a homogeneous
ponderable nonconducting medium endowed with general linear constitutive
relations in order to obtain the dispersion relations, which provide the
refractive index and the propagating modes. From Eq. (1), Ampère’s law reads

${\varepsilon_{ijk}}k^{j}H^{k}+\omega D^{i}=0,$ (14)

where $\varepsilon_{ijk}$ is the tridimensional Levi-Civita symbol. Replacing
the constitutive relations (6a) and (6b) in Eq. (14), one has

${\varepsilon_{ijk}}k^{j}\left(\beta_{ka}E^{a}+\zeta_{ka}B^{a}\right)+\omega\epsilon_{ij}E^{j}+\omega\alpha_{ij}B^{j}=0,$
(15)

Employing now Faraday’s law, $\omega\mathbf{B}=\mathbf{k\times E}$, one
obtains an equation totally in terms of the electric field,

$\displaystyle 0$ $\displaystyle=$
$\displaystyle\varepsilon_{ijk}\varepsilon_{amn}\zeta_{ka}k^{j}k^{m}E^{n}+\omega^{2}\epsilon_{ij}E^{j}+$
(16)
$\displaystyle+\omega\varepsilon_{jmn}\alpha_{ij}k^{m}E^{n}+\omega\varepsilon_{ijk}\beta_{ka}k^{j}E^{a}.$

Let us consider that the medium has isotropic both the electric permittivity
and magnetic permeability,

$\zeta_{ka}={\mu}^{-1}\delta_{ka},\text{ \ \
}\epsilon_{ij}=\epsilon\delta_{ij},$ (17)

in such a way that the anisotropy, typical of “chiral” media, is allowed to
exist in the magnetoelectric coefficients. Hence, Eq.(16) becomes

$\left[\mathbf{k}\times\left(\mathbf{k}\times\mathbf{E}\right)\right]^{i}+\omega^{2}\mu\bar{\epsilon}_{ij}E^{j}=0,$
(18)

where

$\bar{\epsilon}_{in}(\omega)=\epsilon\delta_{in}-\frac{1}{\omega}\left(\beta_{kn}{\varepsilon_{imk}}+\alpha_{ij}{\varepsilon_{jmn}}\right)k_{m},$
(19)

defines the frequency-dependent extended permittivity tensor, which carries
the electric and magnetic response of the medium. Equation (18) is also cast
in the form

$\left[k^{2}\delta_{ij}-k_{i}k_{j}-\omega^{2}\mu\bar{\epsilon}_{ij}\right]E^{j}=0.$
(20)

For a general anisotropic continuous scenario, we write
$\mathbf{k}=\omega\mathbf{n}$ where $\mathbf{n}$ is a vector pointing along
the direction of the wave vector and yields the refractive index:
$n=+\sqrt{\mathbf{n}^{2}}$. Here we consider that the index $n$ is nonnegative
and $\sqrt{\mathbf{n}^{2}}$ instead of $|\mathbf{n}|$, in order to permit
complex refractive indices. The refractive indices with negative real parts,
related to metamaterials, are not considered here. Hence, Eq.(20) becomes

$M_{ij}E^{j}=0,$ (21)

where the tensor $M_{ij}$ reads

$M_{ij}=n^{2}{\delta}_{ij}-n_{i}n_{j}-\mu\bar{{\epsilon}}_{ij},$ (22)

and $\bar{{\epsilon}}_{ij}$ is given by Eq.(19). This set of equations has a
nontrivial solution for the electric field if the determinant of the matrix
$M_{ij}$ vanishes. Such a condition provides the dispersion relations that
govern the wave propagation in the medium. In the case of standard media
described by anisotropic tensors $\epsilon_{ij}$ and $\zeta_{ij}$, and with no
extensions on the constitutive relations, $\alpha_{ij}=0$ and $\beta_{ij}=0$,
the dispersion relation can be found in Refs. Yakov ; Damaskos .

We next examine the propagation of electromagnetic waves in a dielectric
medium under the validity of anisotropic extended dispersion relations of the
form

$\displaystyle D^{i}$ $\displaystyle=$
$\displaystyle\epsilon\delta_{ij}E^{j}+\alpha_{ij}B^{j},$ (23a) $\displaystyle
H^{i}$ $\displaystyle=$
$\displaystyle\beta_{ij}E^{j}+\mu^{-1}\delta_{ij}B^{j},$ (23b) where Eqs. (17)
were considered. These relations may be considered isotropic-anisotropic since
they contain isotropic electric permittivity and magnetic permeability, but
anisotropic magnetoelectric tensors, $\alpha_{ij},\beta_{ij}$.

First, we present the general conditions on the constitutive tensors in order
to ensure energy conservation in the system. The Poynting theorem is given by
Jackson

$\displaystyle\nabla\cdot{\bf{S}}$
$\displaystyle=-\frac{\mathrm{i}\omega}{2}\left({\bf{E}}\cdot{\bf{D}}^{*}-{\bf{H}}^{*}\cdot{\bf{B}}\right)-\frac{({\bf{J}}^{*}\cdot{\bf{E}})}{2},$
(23c)

where

$\displaystyle{\bf{S}}=\frac{1}{2}\left({\bf{E}}\times{\bf{H}}^{*}\right),$
(23d)

is the Poynting vector. The real part of Eq. (23c) yields the energy
conservation law for the system. In the absence of sources and considering
there is no flux of energy density, the energy conservation condition is

$\displaystyle\mathrm{Re}\left[\mathrm{i}\omega\left({\bf{D}}^{\dagger}{\bf{E}}-{\bf{H}}^{\dagger}{\bf{B}}\right)\right]=0.$
(23e)

Depending on the form of the constitutive relations, constraints on the
parameters describing the medium (compatible with energy conservation) are
obtained. Indeed, replacing the constitutive relations in Eqs. (23) into Eq.
(23e) yields

$\displaystyle 0$
$\displaystyle=\mathbf{E}^{{\dagger}}\left(\hat{\epsilon}^{{\dagger}}-\hat{\epsilon}\right)\mathbf{E}-\mathbf{B{{}^{{\dagger}}}}\left(\hat{\zeta}\mathbf{{{}^{{\dagger}}}}-\hat{\zeta}\right)\mathbf{B}$
$\displaystyle+\mathbf{B}^{{\dagger}}\left(\hat{\beta}+\hat{\alpha}^{{\dagger}}\right)\mathbf{E}-\mathbf{E{{}^{{\dagger}}}}\left(\hat{\beta}\mathbf{{{}^{{\dagger}}}}+\hat{\alpha}\right)\mathbf{B},$
(23f)

which establishes a general relation involving all constitutive tensors with
the electromagnetic fields. A simple route to ensure energy conservation is to
set

$\displaystyle\hat{\epsilon}^{\dagger}=\hat{\epsilon},\quad\hat{\zeta}^{\dagger}=\zeta,\quad\hat{\beta}=-\hat{\alpha}^{\dagger}.$
(23g)

The last condition will be relevant in the discussions of the next sections,
as we will see. In the following, we write the dispersion relations from which
we obtain the refractive indices and the propagating modes for some special
configurations of $\alpha_{ij},\beta_{ij}$.

In the next sections, we study the dispersion relations, refractive indices,
propagating modes, group velocity, phase velocity, and Poynting vector of the
electromagnetic waves for the bi-isotropic and bi-anisotropic linear media.

##  III Bi-isotropic case

In the context of the constitutive relations (23), we begin considering the
total symmetric isotropic configuration, where the quantities $\hat{\alpha}$
and $\hat{\beta}$ are given by

$\alpha_{ij}=\alpha\delta_{ij},\quad\beta_{ij}=\beta\delta_{ij},$ (24)

with $\alpha$, $\beta\in\mathbb{C}$. The condition
$\alpha_{ij}=-\beta_{ij}^{\dagger}$, when applied on parametrization (24),
yields

$\displaystyle\beta^{*}=-\alpha.$ (25)

In this case, the constitutive relations take on the typical bi-isotropic
form,

$\displaystyle\mathbf{D}$ $\displaystyle=\epsilon\mathbf{E}+\alpha\mathbf{B},$
(26a) $\displaystyle\mathbf{H}$
$\displaystyle=\frac{1}{\mu}\mathbf{B}+\beta\mathbf{E},$ (26b) which represent
the simplest linear connection between $(\mathbf{D},\mathbf{H})$ and
$(\mathbf{E},\mathbf{B})$. As already mentioned, such relations play a
relevant role in topological insulators Chang ; Urrutia ; Urrutia2 ; Lakhtakia
; Winder ; Li ; Li1  and axion systems Sekine ; Tobar ; BorgesAxion .

Inserting Eq. (24) in Eq. (19), one obtains

$\bar{\epsilon}_{ij}=\epsilon\delta_{ij}+(\alpha+\beta){\varepsilon_{ijm}}n_{m},$
(27)

where the last term on the right-hand side represents the “magnetic-electric”
response of the medium. As we have started with isotropic tensors,
$\epsilon\delta_{ij}$, $\mu^{-1}\delta_{ij}$, $\alpha\delta_{ij}$, and
$\beta\delta_{ij}$, any effective arising anisotropy comes from the extended
structure of the constitutive relations (26). In this case, the tensor
$M_{ij}$ [Eq. (22)] has the form

$\displaystyle{M\equiv[M_{ij}]}$
$\displaystyle=\mathcal{N}-\mu(\alpha+\beta)\begin{pmatrix}0&n_{3}&-n_{2}\\\
-n_{3}&0&n_{1}\\\ n_{2}&-n_{1}&0\end{pmatrix},$ (28)

where

$\mathcal{N}=\begin{pmatrix}n_{2}^{2}+n_{3}^{2}-\mu\epsilon&-n_{1}n_{2}&-n_{1}n_{3}\\\
-n_{1}n_{2}&n_{1}^{2}+n_{3}^{2}-\mu\epsilon&-n_{2}n_{3}\\\
-n_{1}n_{3}&-n_{2}n_{3}&n_{1}^{2}+n_{2}^{2}-\mu\epsilon\end{pmatrix}.$ (29)

Requiring $\mathrm{det}[M_{ij}]=0$, one gets

$n^{4}-n^{2}\left[2\mu\epsilon-\mu^{2}(\alpha+\beta)^{2}\right]+\mu^{2}\epsilon^{2}=0.$
(30)

Solving for $n^{2}$, we obtain the following refractive indices

$\displaystyle n_{\pm}^{2}=$
$\displaystyle\mu{\epsilon}-2Z\pm\mathrm{i}\mu(\alpha+\beta)\sqrt{\mu\epsilon-Z},$
(31)

where

$Z=\frac{\mu^{2}(\alpha+\beta)^{2}}{4}.$ (32)

Thereby the corresponding $n_{\pm}$ read

$n_{\pm}=\sqrt{\mu\epsilon-Z}\pm\mathrm{i}\sqrt{Z},$ (33)

where we have considered only the indices with a positive real piece in order
to avoid metamaterial behavior. The refractive indices (33) are valid (and are
equal) for any propagation direction since the bi-isotropic case does not have
a preferred direction that could be represented by a constant vector. In spite
of that, the system may manifest anisotropic effect (circular birefringence)
due to the way the fields are coupled. Such an effect will be examined ahead.

The refractive indices can also be obtained by diagonalizing the electric
permittivity and setting each eigenvalue equal to $n^{2}/\mu$. The eigenvalues
$\epsilon_{a}$ ($a=1,2,3$) fulfill
$\bar{\epsilon}_{ij}e_{a}^{j}=\epsilon_{a}e_{a}^{i}$, where $\mathbf{e}_{a}$
represent the eigenvectors. Diagonalizing the matrix of the operator
${\bar{\epsilon}}$ [Eq. (27)], one finds the following eigenvalues:

$\displaystyle\epsilon_{1}$ $\displaystyle=\epsilon,$ (34)
$\displaystyle\epsilon_{2,3}$
$\displaystyle\equiv\epsilon_{\pm}=\epsilon\pm\mathrm{i}(\alpha+\beta)n,$ (35)

associated with the eigenvectors

$\displaystyle{\mathbf{e}}_{1}$ $\displaystyle=\frac{\mathbf{n}}{n},$ (36)
$\displaystyle\mathbf{e}_{2(3)}$
$\displaystyle=\frac{1}{n\sqrt{2(n^{2}_{1}+n_{3}^{2})}}\begin{pmatrix}n_{3}n\pm\mathrm{i}n_{1}n_{2}\\\
\mp\mathrm{i}(n_{1}^{2}+n^{2}_{3})\\\
\pm\mathrm{i}n_{2}n_{3}-n_{1}n\end{pmatrix}.$ (37)

Eigenvalues (34) and (35) are associated with the refractive indices

$\displaystyle n^{2}$ $\displaystyle=\mu\epsilon,$ (38) $\displaystyle
n_{\pm}^{2}$ $\displaystyle=\mu\epsilon\pm\mathrm{i}\mu(\alpha+\beta)n.$ (39)

We note that Eq.(38) represents the refractive index of an isotropic
dielectric medium. On the other hand, Eq.(39) recovers the result of Eq. (33),
meaning that only the eigenvalues $\epsilon_{2}$ and $\epsilon_{3}$ correspond
to the refractive indices of the medium $n_{+}$ and $n_{-}$, respectively.
This approach of finding the refractive indices $n$ via the relation
$n^{2}=\mu\epsilon_{a}(n)$, where $\epsilon_{a}$ stands for the eigenvalues of
the electric permittivity ${\bar{\epsilon}}_{ij}$, only works when the
electric field is orthogonal to the propagation direction. Here, such a
condition is guaranteed by the Gauss’ law ${\bf{k}}\cdot{\bf{D}}=0$, where the
electric displacement vector, given by

$\displaystyle{\bf{D}}$
$\displaystyle=\epsilon{\bf{E}}+\frac{\alpha}{\omega}{\bf{k}}\times{\bf{E}},$
(40)

provides ${\bf{k}}\cdot{\bf{E}}=0$. Then for a general vector
$\mathbf{n}=\mathbf{k}/\omega$, the related propagating electric field,
${\mathbf{E}}_{a}$, satisfies ${\mathbf{n}}\cdot{\mathbf{E}}_{a}=0$. This way,
Eqs.(21) and (22) simplify to

$\left[n^{2}\delta_{ij}-\mu\bar{\epsilon}_{ij}\right]E^{j}=0,$ (41)

or $n^{2}\delta_{ij}=\mu\bar{\epsilon}_{ij}$, creating the straightforward
correspondence between $n^{2}$ and $\bar{\epsilon}_{ij}$ eigenvalues, that is,
$n^{2}=\mu\epsilon_{a}(n)$. This is the reason by which the eigenvectors
($e_{a})$ represent the electric field modes,
$\mathbf{E}_{a}\sim\mathbf{e}_{a}$. Note that it also holds that
${\mathbf{n}}\cdot{\mathbf{e}}_{a}=0$. This situation is clearly illustrated
in the present case. In fact, the three normalized eigenvectors given in Eqs.
(36) and (37) are linearly independent and obey

$\mathbf{e}_{1}\cdot\mathbf{e}_{2}^{\ast}=\mathbf{e}_{1}\cdot\mathbf{e}_{3}^{\ast}=\mathbf{e}_{2}\cdot\mathbf{e}_{3}^{\ast}=0.$
(42)

In particular, ${\mathbf{e}}_{2}$ and ${\mathbf{e}}_{3}$ are orthogonal to
${\mathbf{e}}_{1}\sim{\mathbf{n}}$, thus indicating the transversality of the
propagating modes, $\mathbf{E}_{2}\sim\mathbf{e}_{2,}$
$\mathbf{E}_{3}\sim\mathbf{e}_{3}$, whose eigenvalues yield the correct
refractive indices $n_{\pm}$ [see Eq.(33)]. On the other hand,
${\mathbf{n}}\cdot{\mathbf{e}}_{1}$ is nonzero, meaning that the eigenvalue
$\epsilon_{1}$ does not yields a physical refractive index.

###  III.1 Propagation modes

As already explained, the electric field of the propagating modes is given by
solution (37). So, let us choose a a convenient coordinate system where

$\mathbf{n}=(0,0,n_{3}),$ (43)

with which the eigenvectors (37) are

$\mathbf{e}_{2(3)}=\frac{1}{\sqrt{2}}\begin{pmatrix}1\\\ \mp\mathrm{i}\\\
0\end{pmatrix},$ (44)

where $-\mathrm{i}$ represents a right-handed circular polarization (RCP) and
$+\mathrm{i}$ a left-handed circular polarization (LCP), respectively. We can
easily show that the same result stems directly from (21). By replacing the
simple choice of Eq.(43) in the matrix (28),

${M}=\begin{pmatrix}n_{3}^{2}-\mu\epsilon&&-\mu(\alpha+\beta)n_{3}&&0\\\
\mu(\alpha+\beta)n_{3}&&n_{3}^{2}-\mu\epsilon&&0\\\
0&&0&&-\mu\epsilon\end{pmatrix},$ (45)

and implementing the refractive indices (31), the condition $M_{ij}E^{j}=0$
provides the following normalized solutions of the electric field of the
propagating modes:

$\hat{\mathbf{E}}_{\pm}=\frac{1}{\sqrt{2}}\begin{pmatrix}1\\\ \pm\mathrm{i}\\\
0\end{pmatrix},$ (46)

where $\hat{\bf{E}}_{+}$ and $\hat{\bf{E}}_{-}$ represent the LCP and RCP
vectors, respectively. Solution (46) does not depend on the nature (real or
complex) of the parameters $\alpha$ and $\beta$, in such a way it will be
valid for all the cases examined in this section. The equality between
solution (46) and Eq. (44), $\mathbf{e}_{3(2)}\equiv\hat{\mathbf{E}}_{\pm}$,
confirms the approach here developed.

We point out that the circular polarization solution (46) is not exclusive of
the $z$-propagation direction. Indeed, taking on the propagation in the
$x$-axis, $\mathbf{n}=(n_{1},0,0)$, matrix (28) takes the form,

${M}=\begin{pmatrix}-\mu\epsilon&&0&&0\\\
0&&n_{1}^{2}-\mu\epsilon&&-\mu(\alpha+\beta)n_{1}\\\
0&&\mu(\alpha+\beta)n_{1}&&n_{1}^{2}-\mu\epsilon\end{pmatrix},$ (47)

whose associated modes,

$\hat{\mathbf{E}}_{\pm}=\frac{1}{\sqrt{2}}\begin{pmatrix}0\\\ \pm\mathrm{i}\\\
1\end{pmatrix},$ (48)

also correspond to transversal circularly polarized waves.

###  III.2 Optical effects of complex magnetoelectric parameters in
dielectrics

Since we have already found the refractive indices and the polarization of the
propagating modes, it is necessary to examine the physical behavior brought
about by the constitutive relations (26a) and (26b) on a conventional
dielectric substrate. In the limit $(\alpha+\beta)\rightarrow 0$, one recovers
the refractive index of an isotropic dielectric medium, given by Eq. (38),

$n_{\pm}^{2}=\mu\epsilon.$ (49)

Equation (33) provides

$n_{\pm}=\sqrt{\mu\epsilon-\frac{\mu^{2}(\alpha+\beta)^{2}}{4}}\pm\mathrm{i}\frac{\mu(\alpha+\beta)}{2}.$
(50)

Now we examine the refractive indices (50) in two cases: (a) $\alpha$,
$\beta\in\mathbb{C}$ and (b) $\alpha$, $\beta\in\mathbb{R}$. For $\alpha$ and
$\beta$ complex, one can write

$\alpha=\alpha^{\prime}+\mathrm{i}\alpha^{\prime\prime},\quad\beta=\beta^{\prime}+\mathrm{i}\beta^{\prime\prime},$
(51)

where $\alpha^{\prime}=\mathrm{Re}[\alpha]$,
$\alpha^{\prime\prime}=\mathrm{Im}[\alpha]$,
$\beta^{\prime}=\mathrm{Re}[\beta]$ and
$\beta^{\prime\prime}=\mathrm{Im}[\beta]$. Condition (25) implies

$\alpha^{\prime}=-\beta^{\prime},\quad\alpha^{\prime\prime}=\beta^{\prime\prime},$
(52)

so that $\alpha+\beta=2{\mathrm{i}\alpha^{\prime\prime}}$. Therefore, Eq. (50)
is rewritten as

$n_{\pm}=\sqrt{\mu\epsilon+\mu^{2}{\alpha^{\prime\prime}}^{2}}\mp\mu\alpha^{\prime\prime},$
(53)

which are real, positive, and cause birefringence. Since the polarization
modes are circularly polarized vectors [see Eq. (46)], the birefringence
effect can be evaluated in terms of the rotatory power (see Appendix A),
defined as

$\delta=-\frac{[\mathrm{Re}(n_{+})-\mathrm{Re}(n_{-})]\omega}{2}\ $ (54)

Hence, using indices (53), the rotatory power is

$\delta=\mu\omega\alpha^{\prime\prime}.$ (55)

Such a birefringence effect [Eq. (55)] is a consequence of
$(\alpha+\beta)=2\mathrm{i}{\alpha}^{\prime\prime}$. Therefore, it only occurs
when the constitutive parameters possess an imaginary piece. On the other
hand, for $\alpha$, $\beta\in\mathbb{R}$, one has simply $\beta=-\alpha$,
$\alpha^{\prime\prime}=0$, and no birefringence takes place. This is the case
of the topological insulators bi-isotropic scenario Chang ; Urrutia ; Urrutia2
; Lakhtakia ; Winder ; Li ; Li1 , whose constitutive relations are

$\displaystyle\mathbf{D}$ $\displaystyle=$
$\displaystyle\epsilon\mathbf{E}-\alpha_{0}\mathbf{B},$ (56a)
$\displaystyle\mathbf{H}$ $\displaystyle=$
$\displaystyle\frac{\mathbf{B}}{\mu}+\alpha_{0}\mathbf{E},$ (56b)

with $\alpha_{0}=e^{2}/4\pi\hbar$ and $e$ being the elementary electric
charge. For Eq.(56), one has $(\alpha+\beta)=0$, so that no birefringence is
provided.

Concerning topological insulators, quantum effects of bulk interband
excitations may generate strong Faraday rotation associated with a type of
optical activity described by the Verdet constant Ohnoutek ; Liang . The
quantum origin of this effect does not represent a contradiction with the
classical absence of birefringence for topological insulators here remarked.

###  III.3 Group velocity, phase velocity, and Poynting vector

Using ${\bf{n}}={\bf{k}}/\omega$ in Eq. (33), one finds

$\displaystyle\omega_{\pm}$
$\displaystyle=\frac{k}{\sqrt{\mu\epsilon-Z}\pm\mathrm{i}\sqrt{Z}}.$ (57)

To assess the group and phase velocities, we need to consider the nature (real
or complex) of the magnetoelectric parameters.

  * •

For $\alpha$, $\beta\in\mathbb{C}$, it holds that
$(\alpha+\beta)=2\mathrm{i}\alpha^{\prime\prime}$ and
$Z=-\mu^{2}\alpha^{\prime\prime 2}$, so that

$\displaystyle\omega_{\pm}$
$\displaystyle=\frac{k}{\sqrt{\mu\epsilon+\mu\alpha^{\prime\prime
2}}\mp\mu\alpha^{\prime\prime}}.$ (58)

In this case, $\omega_{\pm}>0$ which guarantees propagation of physical modes
for all values of $k$. The phase and group velocities are equal,

$\displaystyle v_{\mathrm{ph}(\pm)}$
$\displaystyle\equiv\frac{\omega_{\pm}}{k}=\frac{1}{\sqrt{\mu\epsilon+\mu\alpha^{\prime\prime
2}}\mp\mu\alpha^{\prime\prime}},$ (59) $\displaystyle v_{\mathrm{g}(\pm)}$
$\displaystyle\equiv\left|\frac{\partial\omega_{\pm}}{\partial{\bf{k}}}\right|=\frac{1}{\sqrt{\mu\epsilon+\mu\alpha^{\prime\prime
2}}\mp\mu\alpha^{\prime\prime}}.$ (60)

  * •

For $\alpha$, $\beta\in\mathbb{R}$, $(\alpha+\beta)=Z=0$, one has

$\displaystyle\omega_{\pm}$ $\displaystyle=\frac{k}{\sqrt{\mu\epsilon}},$ (61)

which yields

$\displaystyle
v_{\mathrm{ph}(\pm)}=v_{\mathrm{g}(\pm)}=\frac{1}{\sqrt{\mu\epsilon}}.$ (62)

Since $v_{\mathrm{g}(\pm)}<1$ in both cases $i)$ and $ii)$, the classical
causality is ensured for all $k$ and any value of $\alpha^{\prime\prime}$.

To examine the energy flux propagation direction, we implement Faraday’s law
and the constitutive relation (26b) in the Poynting vector (23d), yielding

$\displaystyle{\bf{S}}$
$\displaystyle=\frac{1}{2\mu}\left[{\bf{n}}|{\bf{E}}|^{2}-({\bf{n}}\cdot{\bf{E}}){\bf{E}}^{*}\right]+\frac{\beta^{*}}{2}({\bf{E}}\times{\bf{E}}^{*}).$
(63)

In the absence of sources, the Gauss law (${\bf{k}}\cdot{\bf{D}}=0$), taking
into account Eq. (40), provides ${\bf{k}}\cdot{\bf{E}}=0$. Thus Eq. (63) is
rewritten as

$\displaystyle{\bf{S}}$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}+\frac{\beta^{*}}{2}({\bf{E}}\times{\bf{E}}^{*}).$
(64)

The real part of Eq. (64) provides the time-averaged Poynting vector, that is,

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}+\mathrm{Re}\left[\frac{\beta^{*}}{2}({\bf{E}}\times{\bf{E}}^{*})\right],$
(65)

where we have used $\mathrm{Re}[{\bf{n}}]={\bf{n}}=n\hat{\bf{n}}$, since the
refractive indices are real. Using the property $\mathrm{Re}[z]=(z+z^{*})/2$,
with a complex $z$,
${\bf{E}}={\bf{E}}^{\prime}+\mathrm{i}{\bf{E}}^{\prime\prime}$, and
$\beta=-\alpha^{*}$, the simplified Poynting vector takes the form

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}-\alpha^{\prime\prime}({\bf{E}}^{\prime}\times{\bf{E}}^{\prime\prime}).$
(66)

The Gauss law, ${\bf{n}}\cdot{\bf{E}}=0$, implies
${\bf{n}}\cdot{\bf{E}}^{\prime}=0$ and
${\bf{n}}~{}\cdot~{}{\bf{E}}^{\prime\prime}=0$, so both vectors
${\bf{E}}^{\prime}$ and ${\bf{E}}^{\prime\prime}$ are in the plane orthogonal
to $\bf{n}$. This way, the product
${\bf{E}}^{\prime}\times{\bf{E}}^{\prime\prime}$ is parallel or antiparallel
to $\bf{n}$. Therefore, in this bi-isotropic medium, the energy flux
propagates along the same direction of the electromagnetic wave, independently
of the value of the magnetoelectric parameter, $\alpha^{\prime\prime}$,
responsible for the birefringence.

##  IV Bi-anisotropic case with symmetric parameters

Now we explore the scenario where $\alpha_{ij}$ and $\beta_{ij}$ are
nondiagonal symmetric tensors, while the electric permittivity, $\epsilon$,
and the magnetic permeability, $\mu$, are simple numbers. They can be easily
parametrized by using a single 3-vector $\mathbf{d}$, that is,

$\alpha_{ij}=\tilde{\alpha}d_{i}d_{j},\quad\beta_{ij}=\tilde{\beta}d_{i}d_{j},$
(67)

in such a way that the constitutive relations take the form

$\displaystyle\mathbf{D}$
$\displaystyle=\epsilon\mathbf{E}+\tilde{\alpha}\mathbf{d}(\mathbf{d}\cdot\mathbf{B}),$
(68a) $\displaystyle\mathbf{H}$
$\displaystyle=\frac{1}{\mu}\mathbf{B}+\tilde{\beta}\mathbf{d}(\mathbf{d}\cdot\mathbf{E}).$
(68b)

The parameters $\alpha_{ij}$ and $\beta_{ij}$ in Eqs. (67) represent symmetric
matrices with trace given by ${\tilde{\alpha}\mathbf{d}^{2}}$ and
${\tilde{\beta}\mathbf{d}^{2}}$, respectively. These $3\times 3$ matrices
contain off-diagonal elements that could yield anisotropy. This is a bold
difference in relation to the bi-isotropic configuration (24), examined in
Sec. III. For a matter of generality, we suppose
$\alpha_{ij},\beta_{ij}\in\mathbb{C}$, which is compatible with
$\mathbf{d}\in\mathbb{R}^{3}$ and $\tilde{\alpha},\tilde{\beta}\in\mathbb{C}$.
This way, in accordance with condition (8c), the parameters (67) should obey

$\tilde{\beta}=-\tilde{\alpha}^{*},$ (69)

which implies

$\tilde{\alpha}^{\prime}=-\tilde{\beta}^{\prime},\quad\tilde{\alpha}^{\prime\prime}=\tilde{\beta}^{\prime\prime},$
(70)

for $\tilde{\alpha}^{\prime}=\mathrm{Re}[\tilde{\alpha}]$,
$\tilde{\alpha}^{\prime\prime}=\mathrm{Im}[\tilde{\alpha}]$,
$\tilde{\beta}^{\prime}=\mathrm{Re}[\tilde{\beta}]$ and
$\tilde{\beta}^{\prime\prime}=\mathrm{Im}[\tilde{\beta}]$.

The constitutive relation (68a) allows us to write the displacement vector in
the form

$\mathbf{D}=\epsilon\mathbf{E}+\frac{\tilde{\alpha}}{\omega}\mathbf{d}\left[\mathbf{d}\cdot(\mathbf{k}\times\mathbf{E})\right].$
(71)

The Gauss law, $\mathbf{k}\cdot\mathbf{D}=0$, requires electric field
configurations satisfying

$\left[\epsilon\mathbf{n}+\tilde{\alpha}(\mathbf{n}\cdot\mathbf{d})(\mathbf{d}\times\mathbf{n})\right]\cdot\mathbf{E}=0,$
(72)

for $\mathbf{k}=\omega\mathbf{n}$.

In the case $\mathbf{d}$ and $\mathbf{n}$ are parallel vectors, Eq.(72)
implies transversal electric field modes, $\mathbf{n}\cdot\mathbf{E}=0$. For
nonparallel vectors $\mathbf{d}$ and $\mathbf{n}$, the electric field may be
written as in Eq. (185) of Appendix A, which does not supply, in general,
transversal modes, that is, $\mathbf{n}\cdot\mathbf{E}\neq 0$; see Eq. (188).
Furthermore, for orthogonal vectors, ${\mathbf{n}\cdot\mathbf{d}=0}$, Eq. (72)
yields again transversal modes, ${\mathbf{n}}\cdot{\mathbf{E}}=0$, with the
electric field expressed as in Eq. (190).

Replacing relations (67) in the permittivity tensor (19), ones writes

$\displaystyle\bar{\epsilon}_{ij}=\epsilon\delta_{ij}-\frac{1}{\omega}\left(\tilde{\beta}{\varepsilon_{imn}}k_{m}d_{n}d_{j}+\tilde{\alpha}{\varepsilon_{amj}}d_{i}d_{a}k_{m}\right),$
(73)

in such a way the tensor $M_{ij}$, Eq. (22), provides

$\displaystyle{M}=\mathcal{N}-\mu(\mathcal{D}+\mathcal{E}),$ (74)

with $\mathcal{N}$ of Eq.(29), containing the usual constitutive elements,
while the magnetoelectric contributions are displayed in the following:

$\displaystyle\mathcal{D}$
$\displaystyle=-(\tilde{\alpha}-\tilde{\beta})\mathrm{diag}\left(D_{1},D_{2},D_{3}\right),$
(75a) $\displaystyle\mathcal{E}$
$\displaystyle=\begin{pmatrix}0&\epsilon_{12}&\epsilon_{13}\\\
\epsilon_{21}&0&\epsilon_{23}\\\ \epsilon_{31}&\epsilon_{32}&0\end{pmatrix},$
(75b) where $\displaystyle D_{1}$
$\displaystyle=d_{1}(d_{2}n_{3}-d_{3}n_{2}),$ (75c) $\displaystyle D_{2}$
$\displaystyle=d_{2}(d_{3}n_{1}-d_{1}n_{3}),$ (75d) $\displaystyle D_{3}$
$\displaystyle=d_{3}(d_{1}n_{2}-d_{2}n_{1}),$ (75e)

and

$\displaystyle\epsilon_{12}$
$\displaystyle=-\tilde{\beta}d_{2}(d_{3}n_{2}-d_{2}n_{3})+\tilde{\alpha}d_{1}(d_{1}n_{3}-d_{3}n_{1}),$
(76a) $\displaystyle\epsilon_{13}$
$\displaystyle=-\tilde{\beta}d_{3}(d_{3}n_{2}-d_{2}n_{3})+\tilde{\alpha}d_{1}(d_{2}n_{1}-d_{1}n_{2}),$
(76b) $\displaystyle\epsilon_{21}$
$\displaystyle=-\tilde{\beta}d_{1}(d_{1}n_{3}-d_{3}n_{1})+\tilde{\alpha}d_{2}(d_{3}n_{2}-d_{2}n_{3}),$
(76c) $\displaystyle\epsilon_{23}$
$\displaystyle=-\tilde{\beta}d_{3}(d_{1}n_{3}-d_{3}n_{1})+\tilde{\alpha}d_{2}(d_{2}n_{1}-d_{1}n_{2}),$
(76d) $\displaystyle\epsilon_{31}$
$\displaystyle=-\tilde{\beta}d_{1}(d_{2}n_{1}-d_{1}n_{2})+\tilde{\alpha}d_{3}(d_{3}n_{2}-d_{2}n_{3}),$
(76e) $\displaystyle\epsilon_{32}$
$\displaystyle=-\tilde{\beta}d_{2}(d_{2}n_{1}-d_{1}n_{2})+\tilde{\alpha}d_{3}(d_{1}n_{3}-d_{3}n_{1}).$
(76f) Evaluating $\mathrm{det}[M_{ij}]=0$, we obtain the dispersion relation,
$\displaystyle{\epsilon}\left(n^{2}-\mu\epsilon\right)^{2}+\tilde{\alpha}\tilde{\beta}\mu\left[\mu{\epsilon}d^{2}-(\mathbf{n}\cdot\mathbf{d})^{2}\right]{\left(\mathbf{d}\times\mathbf{n}\right)^{2}}=0,$
(77)

where $\left(\mathbf{d}\times\mathbf{n}\right)^{2}\equiv
d^{2}\mathbf{n}^{2}-\left(\mathbf{d}\cdot\mathbf{n}\right)^{2}$.

Relation (69) provides $\tilde{\alpha}\tilde{\beta}=-|\tilde{\alpha}|^{2}$.
Furthermore, implementing $\mathbf{n}\cdot\mathbf{d}=nd\cos\varphi$, Eq.(77)
yields

$n_{\pm}^{2}=\frac{1}{s}\left[N\pm\mu|\tilde{\alpha}|d^{2}\sin^{2}\varphi\sqrt{\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{4}}\right],$
(78)

or

$n_{\pm}=\sqrt{\frac{N+\mu\epsilon\sqrt{s}}{2s}}\pm\sqrt{\frac{N-\mu\epsilon\sqrt{s}}{2s}},$
(79)

where

$\displaystyle N$
$\displaystyle=\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{2}\sin^{2}\varphi.$
(80a) $\displaystyle s$
$\displaystyle=1+\frac{\mu}{\epsilon}|\tilde{\alpha}|^{2}d^{4}\sin^{2}\varphi\cos^{2}\varphi.$
(80b)

We notice in Eq.(79) two distinct refractive indices, both real and positive,
$n_{\pm}>0$, in such a way that birefringent electromagnetic propagation is
expected in this medium. We also highlight that Eq. (78) holds for any
propagation direction in relation to the vector $\mathbf{d}$, generally
expressed in terms of the angle $\varphi$, the relative angle between the
vector ${\bf{d}}$ and the propagation direction.

The behavior of $n_{\pm}$ [Eq. (79)] in terms of $\varphi\in[0,\pi]$ and the
dimensionless parameter $|\tilde{\alpha}|\in[0,1]$ is illustrated in Figs. 1
and 2. The anisotropy effect manifests itself by means of the angular
dependence of $n_{\pm}$ on $\varphi$. Notice the nonlinearity of $n_{\pm}$,
which behaves as a sinusoidal function, increasing with $|\tilde{\alpha}|$.

Figure 1:  Refractive index $n_{+}$ of Eq. (79) with $\mu=1$, $\epsilon=2$,
and $d=1$. The parameters $\epsilon$, $\mu$, and $|\tilde{\alpha}|d^{2}$ are
dimensionless. Figure 2:  Refractive index $n_{-}$ of Eq. (79) with $\mu=1$,
$\epsilon=2$, and $d=1$. The parameters $\epsilon$, $\mu$, and
$|\tilde{\alpha}|d^{2}$ are dimensionless.

In the following, we address the propagating modes and examine birefringence
effects for some specific propagation directions.

###  IV.1 Propagation modes

To examine the propagating modes, we first pay attention to the angle between
$\mathbf{n}$ and $\mathbf{d}$. For the configurations where $\mathbf{n}$ and
$\mathbf{d}$ are in the same direction, Eq.(78) yields
$n_{\pm}=\sqrt{\mu\epsilon}$, which is the usual scenario. Modified scenarios
arise when $\mathbf{n}$ and $\mathbf{d}$ are not aligned. To obtain the modes,
we choose a simplified coordinate system where the vector $\mathbf{n}$ is
along the $z$-axis, that is, $\mathbf{n}=(0,0,n)$. For such a choice, we
investigate a configuration where the background vector is longitudinal and
orthogonal to $\mathbf{n}$.

####  IV.1.1 $\mathbf{d}$-Longitudinal configuration:

We begin examining the $\mathbf{d}$-longitudinal configuration,

$\mathbf{d}=(0,0,d),$ (81)

for which Eq.(78) yields

$n^{2}=\mu\epsilon.$ (82)

For this case, matrix (74) is rewritten as

${M}=\begin{pmatrix}n^{2}-\mu\epsilon&0&0\\\ 0&n^{2}-\mu\epsilon&0\\\
0&0&-\mu\epsilon\end{pmatrix},$ (83)

in such a way that $M_{ij}E^{j}=0$ provides generic orthogonal modes

$\mathbf{E}=\begin{pmatrix}E_{x}\\\ E_{y}\\\ 0\end{pmatrix},$ (84)

representing a transversal wave with undefined polarization (linear, circular,
or elliptical). It is interesting to note that only one positive refractive
index was determined and it does not depend on the propagation direction,
which is a signal of isotropy. This means that the ${\bf{d}}$-direction
defines the optical axis of the medium.

####  IV.1.2 $\mathbf{d}$-transversal configuration

We follow writing the $\mathbf{d}$-transversal configuration,

$\mathbf{d}=(d_{1},d_{2},0),$ (85)

for which $s=1$ and $N=\mu{\epsilon}+|\tilde{\alpha}|^{2}\mu^{2}d^{4}/2$, so
that Eq. (78) yields

$n_{\pm}^{2}=\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{2}\pm\mu|\tilde{\alpha}|d^{2}\sqrt{\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{4}},$
(86)

and

$n_{\pm}=\sqrt{\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{4}}\pm\frac{\mu|\tilde{\alpha}|d^{2}}{2}.$
(87)

By using Eq. (87) we rewrite Eq. (86) as

$n_{\pm}^{2}=\mu\epsilon\pm\mu|\tilde{\alpha}|d^{2}n_{\pm}.$ (88)

Matrix (74) now is

${M}=\begin{pmatrix}n_{\pm}^{2}-\mu\epsilon+\Omega n_{\pm}&-\mu
n_{\pm}(\tilde{\beta}d_{2}^{2}+\tilde{\alpha}d_{1}^{2})&0\\\ \mu
n_{\pm}(\tilde{\beta}d_{1}^{2}+\tilde{\alpha}d_{2}^{2})&n_{\pm}^{2}-\mu\epsilon-\Omega
n_{\pm}&0\\\ 0&0&-\mu\epsilon\end{pmatrix},$ (89)

or better,

${M}=\begin{pmatrix}(\pm\mu|\tilde{\alpha}|d^{2}+\Omega)n_{\pm}&-\mu
n_{\pm}(\tilde{\beta}d_{2}^{2}+\tilde{\alpha}d_{1}^{2})&0\\\ \mu
n_{\pm}(\tilde{\beta}d_{1}^{2}+\tilde{\alpha}d_{2}^{2})&(\pm\mu|\tilde{\alpha}|d^{2}-\Omega)n_{\pm}&0\\\
0&0&-\mu\epsilon\end{pmatrix},$ (90)

where we have used Eqs.(70) and (88), writing

$\displaystyle\Omega$
$\displaystyle=\mu(\tilde{\alpha}-\tilde{\beta})d_{1}d_{2}=2\mu\tilde{\alpha}^{\prime}d_{1}d_{2}.$
(91)

The condition $M_{ij}E^{j}=0$ yields

$\mathbf{E}_{\pm}=E_{0}\begin{pmatrix}1\\\
\displaystyle\frac{\mu\left(\tilde{\beta}d_{1}^{2}+\tilde{\alpha}d_{2}^{2}\right)}{\Omega\mp\mu|\tilde{\alpha}|d^{2}}\\\
0\end{pmatrix},$ (92)

with an appropriately chosen amplitude $E_{0}$. For $d_{1}=0$, we achieve

$\mathbf{E}_{\pm}=E_{0}\begin{pmatrix}1\\\
\displaystyle\mp\frac{\tilde{\alpha}}{|\tilde{\alpha}|}\\\
0\end{pmatrix}=\frac{1}{\sqrt{2}}\begin{pmatrix}1\\\
\displaystyle\mp\frac{\tilde{\alpha}^{\prime}+\mathrm{i}\tilde{\alpha}^{\prime\prime}}{|\tilde{\alpha}|}\\\
0\end{pmatrix},$ (93)

which represents linear polarizations for $\tilde{\alpha}^{\prime\prime}=0$ or
circular polarizations for $\tilde{\alpha}^{\prime}=0$.

As Eq.(87) exhibits two real refractive indices, a scenario with birefringence
is set. The modes (93), however, do not represent RCP or LCP vectors, so the
birefringent propagation cannot be suitably described in terms of the rotatory
power (54). Rather, it can be characterized in terms of the phase shift
arising from the distinct phase velocities of the propagating modes, given by

$\Delta=\frac{2\pi}{\lambda_{0}}l(n_{+}-n_{-}),$ (94)

where $\lambda_{0}$ is the vacuum wavelength of incident light, $l$ is the
thickness of the medium or the distance traveled by the wave, and $n_{+}$ and
$n_{-}$ are the refractive indices of the medium. Note that this is the same
expression that controls the phase shift caused by “retarders” (for details,
see Chap. 8 of Ref. Hecht ). Using Eq.(87), one finds the corresponding phase
shift per unit length as

$\frac{\Delta}{l}=\frac{2\pi}{\lambda_{0}}\mu|\tilde{\alpha}|d^{2}.$ (95)

As the phase shift depends on the modulus of $\tilde{\alpha}$, the
birefringence now takes place for both real and imaginary magnetoelectric
parameters. This is a difference in relation to the bi-isotropic case (26), in
which the birefringence occurs only for imaginary parameters, as shown in
Eq.(55).

####  IV.1.3 $\mathbf{d}$-general configuration

Now, let us analyze the mixed case where the vector $\mathbf{d}$ has
orthogonal and longitudinal components relative to the propagation direction,
$\mathbf{n}$. In this sense, one can set

$\mathbf{d}=(0,d_{2},d_{3}).$ (96)

The refractive indices (78) are rewritten as

$n_{\pm}^{2}=\frac{1}{s}\left(\mu\epsilon+\Lambda_{\pm}\sin^{2}\varphi\right),$
(97)

where $s$ is given by Eq. (80b) and $\Lambda_{\pm}$ is defined as

$\Lambda_{\pm}=\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{2}\pm\mu|\tilde{\alpha}|d^{2}\sqrt{\mu\epsilon+\frac{\mu^{2}|\tilde{\alpha}|^{2}d^{4}}{4}}.$
(98)

For the coordinate system where $\mathbf{n}=(0,0,n_{3})$ and $\mathbf{d}$ is
given by Eq.(96), matrix (74) takes the form

${M}=\begin{pmatrix}n_{3}^{2}-\mu\epsilon&-\mu\tilde{\beta}d_{2}^{2}n_{3}&-\mu\tilde{\beta}d_{2}d_{3}n_{3}\\\
+\mu\tilde{\alpha}d_{2}^{2}n_{3}&n_{3}^{2}-\mu\epsilon&0\\\
+\mu\tilde{\alpha}d_{2}d_{3}n_{3}&0&-\mu\epsilon\end{pmatrix},$ (99)

which becomes

${M}=\begin{pmatrix}1/\Upsilon_{\pm}&-\mu\tilde{\beta}d_{2}^{2}n_{3}&-\mu\tilde{\beta}d_{2}d_{3}n_{3}\\\
&&\\\ +\mu\tilde{\alpha}d_{2}^{2}n_{3}&1/\Upsilon_{\pm}&0\\\ &&\\\
+\mu\tilde{\alpha}d_{2}d_{3}n_{3}&0&-\mu\epsilon\end{pmatrix},$ (100)

when Eq.(97) is taken into account and we have defined

$\Upsilon_{\pm}=\frac{s}{\mu\epsilon\left(1-s\right)+\Lambda_{\pm}\sin^{2}\varphi}.$
(101)

The condition $M_{ij}E^{j}=0$ yields

$\mathbf{E}_{\pm}=E_{0}\begin{pmatrix}1\\\
-\mu\left(\tilde{\alpha}^{\prime}+\mathrm{i}\tilde{\alpha}^{\prime\prime}\right)d_{2}^{2}n_{\pm}\Upsilon_{\pm}\\\
+\left(\tilde{\alpha}^{\prime}+\mathrm{i}\tilde{\alpha}^{\prime\prime}\right)d_{2}d_{3}n_{\pm}/\epsilon\end{pmatrix},$
(102)

where $n_{\pm}$ represents $n_{3\pm}$. Differently from the previous
$\mathbf{d}$-longitudinal or ${\bf{d}}$-orthogonal cases, modes (102) are
endowed with a longitudinal component, a feature of general solutions for
$\varphi\neq 0,\pi$ or $\varphi\neq\pi/2$.

The polarization of modes (102) can be read off their transversal sectors.
Since the transversal piece of Eq.(102) is neither RCP nor LCP, being linearly
(for $\tilde{\alpha}^{\prime\prime}=0$) or elliptically polarized (for
$\tilde{\alpha}^{\prime}=0$), the birefringence effect is expressed by means
of the phase shift per unit length defined in Eq.(94), here carried out as

$\frac{\Delta}{l}=\frac{4\pi}{\lambda_{0}}\sqrt{\frac{N-\mu\epsilon\sqrt{s}}{2s}},$
(103)

where we used Eq. (79). In Fig. 3 we have plotted the phase shift (103) per
unit length in terms of $\varphi$ and $|\tilde{\alpha}|$. We notice that the
birefringence effect is maximal for $\varphi=\pi/2$, which corresponds to
configurations where $n_{+}$ and $n_{-}$ have maxima and minima values,
respectively.

Figure 3:  Phase shift factor per unit length of Eq. (103). Here, we have used
$\mu=1$, $\epsilon=2$, and $d=1$.

As a final comment, note that for $d_{2}=0$ or $\mathbf{d}=(0,0,d_{3})$ the
field (102) simplifies as

$\mathbf{E}_{\pm}=\begin{pmatrix}1\\\ 0\\\ 0\\\ \end{pmatrix},$ (104)

which corresponds to the d-longitudinal case, whose mode is given by Eq.(84),
being compatible with the result (104). On the other hand, for the transversal
configuration, $d_{3}=0$ or $\mathbf{d}=(0,d_{2},0)$, the field (102) yields

$\mathbf{E}_{\pm}=E_{0}\begin{pmatrix}1\\\
-\mu\left(\tilde{\alpha}^{\prime}+\mathrm{i}\tilde{\alpha}^{\prime\prime}\right)d_{2}^{2}n_{\pm}\Upsilon_{\pm}\\\
0\end{pmatrix}.$ (105)

For this transversal configuration it holds that
$\Upsilon_{\pm}=1/\Lambda_{\pm}$ and
$\Lambda_{\pm}=\mu|\tilde{\alpha}|d^{2}n_{\pm}$. With that, solution (105)
recovers the one of Eq.(93).

###  IV.2 Group velocity, phase velocity, and Poynting vector

The dispersion equation (77) can be rewritten in the form

$\displaystyle\omega^{4}-2\omega^{2}\left\\{\frac{{\bf{k}}^{2}}{\mu\epsilon}+\frac{|\tilde{\alpha}|^{2}d^{2}}{2\epsilon^{2}}{({\bf{d}}\times{\bf{k}})^{2}}\right\\}+$
$\displaystyle+\frac{\mu}{\epsilon}\frac{|\tilde{\alpha}|^{2}}{\mu^{2}\epsilon^{2}}({\bf{k}}\cdot{\bf{d}})^{2}{({\bf{d}}\times{\bf{k}})^{2}}+\frac{{\bf{k}}^{4}}{\mu^{2}\epsilon^{2}}$
$\displaystyle=0,$ (106)

whose solutions for $\omega$,

$\displaystyle\omega_{\pm}^{2}$
$\displaystyle=\frac{{\bf{k}}^{2}}{\mu\epsilon}+\frac{|\tilde{\alpha}|}{\epsilon}{({\bf{d}}\times{\bf{k}})^{2}}\left[\frac{|\tilde{\alpha}|d^{2}}{2\epsilon}{\mp}\sqrt{\frac{1}{\mu\epsilon}+\frac{|\tilde{\alpha}|^{2}d^{4}}{4\epsilon^{2}}}\right].$
(107)

provide the phase and group velocities,

$\displaystyle
v_{\mathrm{ph}(\pm)}\equiv\frac{\omega_{\pm}}{k}=\sqrt{\frac{1}{\mu\epsilon}+\frac{|\tilde{\alpha}|}{\epsilon}{({\bf{d}}\times{\hat{\bf{k}}})^{2}}\left(\frac{|\tilde{\alpha}|d^{2}}{2\epsilon}{\mp}f_{\alpha}\right)},$
(108) $\displaystyle v_{g(\pm)}^{i}$
$\displaystyle=\frac{k^{i}}{\mu\epsilon\omega_{\pm}}+\frac{|\tilde{\alpha}|}{\epsilon\omega_{\pm}}\left[d^{2}k^{i}-({\bf{d}}\cdot{\bf{k}})d^{i}\right]\left(\frac{|\tilde{\alpha}|d^{2}}{2\epsilon}{\mp}f_{\alpha}\right),$
(109)

with

$\displaystyle
f_{\alpha}=\sqrt{\frac{1}{\mu\epsilon}+\frac{|\tilde{\alpha}|^{2}d^{4}}{4\epsilon^{2}}}.$
(110)

Both $v_{\mathrm{ph}(\pm)}$ and $v_{\mathrm{g}(\pm)}$ are valid for general
configurations, i. e., for any relative orientation between the vectors
${\bf{d}}$ and ${\bf{n}}$, for which the group velocity is no longer parallel
to ${\bf{n}}$, due to its component along the ${\bf{d}}$ vector. Considering
now the special cases discussed in Sec. IV.1, we state the following:

  * •

For the ${\bf{d}}$-longitudinal scenario where ${\bf{d}}\cdot\hat{\bf{k}}=d$,
Eqs. (108) and (109) provides

$\displaystyle
v_{\mathrm{ph}(\pm)}^{\mathrm{long.}}=\frac{1}{\sqrt{\mu\epsilon}},\quad{\bf{v}}_{\mathrm{g}(\pm)}^{\mathrm{long.}}=\frac{\hat{\bf{k}}}{\sqrt{\mu\epsilon}}.$
(111)

  * •

For the ${\bf{d}}$-transversal case where ${\bf{d}}\cdot\hat{\bf{k}}=0$, one
obtains

$\displaystyle v_{\mathrm{ph}(\pm)}^{\mathrm{trans.}}$
$\displaystyle=\sqrt{\frac{1}{\mu\epsilon}+\frac{|\tilde{\alpha}|^{2}d^{4}}{4\epsilon^{2}}}{\mp}\frac{|\tilde{\alpha}|d^{2}}{2\epsilon},$
(112) $\displaystyle{\bf{v}}_{\mathrm{g}(\pm)}^{\mathrm{trans.}}$
$\displaystyle=\left[\sqrt{\frac{1}{\mu\epsilon}+\frac{|\tilde{\alpha}|^{2}d^{4}}{4\epsilon^{2}}}{\mp}\frac{|\tilde{\alpha}|d^{2}}{2\epsilon}\right]\hat{\bf{k}}.$
(113)

The magnitude of the group velocity (113) is depicted in terms of the
magnetoelectric parameter in Fig. 4.

Figure 4:  Group velocity $v_{\mathrm{g}(\pm)}$ [Eq. (113)] of
${\bf{d}}$-transversal case. The blue dot-dashed curve indicates
$v_{\mathrm{g}(+)}$, while the red solid line represents $v_{\mathrm{g}(-)}$.
The vertical dashed line indicates $|\tilde{\alpha}|=1$, corresponding to the
value of the critical factor $\alpha_{c}$ [Eq. (114)] for the choices $\mu=1$,
$\epsilon=2$ and $d=1$, the ones adopted in this plot.

We note that $v_{g(-)}>1$ occurs for $|\tilde{\alpha}|>\alpha_{c}$, where the
critical value $\alpha_{c}$ is given by

$\displaystyle\alpha_{c}=\frac{\mu\epsilon-1}{\mu d^{2}}.$ (114)

Such a mode with $v_{\mathrm{g}(-)}>1$ is related to the refractive index
$n_{-}$ of Eq. (79), which can assume values $n<1$ (see Fig. 2).

Let us now evaluate the Poynting vector. Using the constitutive relation
(68b), we obtain

$\displaystyle{\bf{S}}$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}-\frac{1}{2\mu}({\bf{n}}\cdot{\bf{E}}){\bf{E}}^{*}+\frac{\beta^{*}}{2}({\bf{E}}\times{\bf{d}})({\bf{d}}\cdot{\bf{E}}^{*}).$
(115)

As the Gauss law yields
$\epsilon{\bf{k}}\cdot{\bf{E}}=-\tilde{\alpha}({\bf{k}}\cdot{\bf{d}})({\bf{d}}\cdot{\bf{B}})$,
Eq. (115) takes the form

$\displaystyle{\bf{S}}$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}+\frac{\tilde{\alpha}}{2\mu\epsilon}({\bf{n}}\cdot{\bf{d}})\left[{\bf{E}}\cdot({\bf{d}}\times{\bf{n}})\right]{\bf{E}}^{*}+$
$\displaystyle\phantom{=}+\frac{\tilde{\beta}^{*}}{2}({\bf{E}}\times{\bf{d}})({\bf{d}}\cdot{\bf{E}}^{*}).$
(116)

Considering $\tilde{\beta}=-\tilde{\alpha}^{*}$ and
${\bf{E}}={\bf{E}}^{\prime}+\mathrm{i}{\bf{E}}^{\prime\prime}$, the time-
averaged Poynting vector, obtained from the real part of Eq. (116), is given
by

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}+\frac{\tilde{\alpha}^{\prime}}{2\mu\epsilon}({\bf{n}}\cdot{\bf{d}})\left\\{\frac{}{}{\bf{E}}^{\prime}\left[{\bf{E}}^{\prime}\cdot({\bf{d}}\times{\bf{n}})\right]\right.+$
$\displaystyle\phantom{=}\left.+{\bf{E}}^{\prime\prime}\left[{\bf{E}}^{\prime\prime}\cdot({\bf{d}}\times{\bf{n}})\right]\right\\}-\frac{\tilde{\alpha}^{\prime\prime}}{2\mu\epsilon}({\bf{n}}\cdot{\bf{d}})\left\\{\frac{}{}{\bf{E}}^{\prime}\left[{\bf{E}}^{\prime\prime}\cdot({\bf{d}}\times{\bf{n}})\right]+\right.$
$\displaystyle\phantom{=}\left.-{\bf{E}}^{\prime\prime}\left[{\bf{E}}^{\prime}\cdot({\bf{d}}\times{\bf{n}})\right]\right\\}-\frac{\tilde{\alpha}^{\prime}}{2}\left[({\bf{d}}\cdot{\bf{E}}^{\prime})({\bf{E}}^{\prime}\times{\bf{d}})+\right.$
$\displaystyle\phantom{=}\left.+({\bf{d}}\cdot{\bf{E}}^{\prime\prime})({\bf{E}}^{\prime\prime}\times{\bf{d}})\right]+\frac{\tilde{\alpha}^{\prime\prime}}{2}\left[({\bf{d}}\cdot{\bf{E}}^{\prime})({\bf{E}}^{\prime\prime}\times{\bf{d}})+\right.$
$\displaystyle\phantom{=}\left.-({\bf{d}}\cdot{\bf{E}}^{\prime\prime})({\bf{E}}^{\prime}\times{\bf{d}})\right],$
(117)

or, equivalently,

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}+\frac{({\bf{n}}\cdot{\bf{d}})}{2\mu\epsilon}\left(\tilde{\alpha}^{\prime}{\bf{f}}^{\prime}-\tilde{\alpha}^{\prime\prime}{\bf{f}}^{\prime\prime}\right)+$
$\displaystyle\phantom{=}-\frac{1}{2}\left(\tilde{\alpha}^{\prime}{\bf{g}}^{\prime}-\tilde{\alpha}^{\prime\prime}{\bf{g}}^{\prime\prime}\right),$
(118)

where

$\displaystyle{\bf{f}}^{\prime,\prime\prime}$
$\displaystyle={\bf{E}}^{\prime}\left[{\bf{E}}^{\prime,\prime\prime}\cdot({\bf{d}}\times{\bf{n}})\right]\pm{\bf{E}}^{\prime\prime}\left[{\bf{E}}^{\prime\prime,\prime}\cdot({\bf{d}}\times{\bf{n}})\right],$
(119) $\displaystyle{\bf{g}}^{\prime,\prime\prime}$
$\displaystyle=({\bf{d}}\cdot{\bf{E}}^{\prime})({\bf{E}}^{\prime,\prime\prime}\times{\bf{d}})\pm({\bf{d}}\cdot{\bf{E}}^{\prime\prime})({\bf{E}}^{\prime\prime,\prime}\times{\bf{d}}).$
(120)

In the case $\mathbf{d}$ and $\mathbf{n}$ are parallel vectors, Eq. (72)
implies transversal electric field modes, $\mathbf{n}\cdot\mathbf{E}=0$, then
$\mathbf{d}\cdot\mathbf{E}=0$. In this case, we have
${\bf{f}}^{\prime,\prime\prime}=0$ and ${\bf{g}}^{\prime,\prime\prime}=0$
yielding simply
$\left\langle{\bf{S}}\right\rangle={\bf{n}}|{\bf{E}}|^{2}/{2\mu}$, and the
energy flux propagates along the same direction of $\mathbf{n}$.

For the case that the vectors $\mathbf{d}$ and $\mathbf{n}$ are mutually
orthogonal, $\mathbf{n}\cdot\mathbf{d}=0$, Eq. (72) also provides
$\mathbf{n}\cdot\mathbf{E}=0$, so that

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}{\bf{n}}|{\bf{E}}|^{2}-\frac{1}{2}\left(\tilde{\alpha}^{\prime}{\bf{g}}^{\prime}-\tilde{\alpha}^{\prime\prime}{\bf{g}}^{\prime\prime}\right).$
(121)

This scenario is such that the vectors $\mathbf{d}$, ${\bf{E}}^{\prime}$, and
${\bf{E}}^{\prime\prime}$ are in the same plane orthogonal to $\mathbf{n}$.
This way, the vectors ${\bf{E}}^{\prime}\times{\bf{d}}$ and
${\bf{E}}^{\prime\prime}\times{\bf{d}}$ are along the ${\bf{n}}$-direction or
are both null (when ${\bf{E}}$ and ${\bf{d}}$ are parallel vectors).
Consequently, the vectors ${\bf{g}}^{\prime,\prime\prime}$ come out parallel
(or antiparallel) to ${\bf{n}}$. Thus, in this case, the energy flux also
propagates along the same propagation direction of the electromagnetic wave,
whatever the $\tilde{\alpha}^{\prime}$, $\tilde{\alpha}^{\prime\prime}$
parameters values.

For the general case in which $\mathbf{n}$ and $\mathbf{d}$ are not collinear
or perpendicular, the energy flux is no longer parallel to the propagation
direction of the electromagnetic wave.

##  V Bi-anisotropic case with antisymmetric parameters

Now we analyze the case where the magnetoelectric parameters are described by
antisymmetric tensors, written as

$\displaystyle\alpha_{ij}$ $\displaystyle={\varepsilon_{ijk}}a_{k},$ (122a)
$\displaystyle\beta_{kn}$ $\displaystyle={\varepsilon_{knr}}b_{r},$ (122b)

where $\mathbf{a}=(a_{x},a_{y},a_{z})$ and $\mathbf{b}=(b_{x},b_{y},b_{z})$
are fixed and, in principle, complex 3-vectors, which induce preferred
direction in the system, while $\varepsilon_{ijk}$ represents the usual Levi-
Civita symbol in three dimensions. In order to satisfy Eq. (8c), the following
condition should hold:

$\mathbf{b}^{*}=\mathbf{a}.$ (123)

In the case the vectors $\mathbf{a}$ and $\mathbf{b}$ are real, this condition
reduces merely to $\mathbf{b}=\mathbf{a}$.

Under the validity of relations (122a) and (122b), the electric displacement
field and the magnetic field are given by

$\displaystyle\mathbf{D}$
$\displaystyle=\epsilon\mathbf{E}+\mathbf{a}\times\mathbf{B},$ (124)
$\displaystyle\mathbf{H}$
$\displaystyle=\frac{1}{\mu}\mathbf{B}+\mathbf{b}\times\mathbf{E}.$ (125)

Some analog antisymmetric constitutive relations have found application in the
description of electron gas systems Carvalho  and in the investigation of
electromagnetic propagation in time-dependent media with an antisymmetric
magnetoelectric coupling and an isotropic time-dependent permittivity Lin .

In the momentum space, Eq.(124) provides

$\displaystyle\mathbf{D}$
$\displaystyle=\epsilon\mathbf{E}+\frac{1}{\omega}\mathbf{a}\times\left(\mathbf{\mathbf{k}\times\mathbf{E}}\right),$
(126) $\displaystyle\mathbf{D}$
$\displaystyle=\left(\epsilon-\frac{\mathbf{a}\cdot\mathbf{\mathbf{k}}}{\omega}\right)\mathbf{\mathbf{E}}+\frac{\mathbf{a}\cdot\mathbf{\mathbf{E}}}{\omega}\mathbf{\mathbf{k}},$
(127)

where
$\mathbf{a}\times\left(\mathbf{\mathbf{k}\times\mathbf{E}}\right)=\left(\mathbf{a}\cdot\mathbf{\mathbf{E}}\right)\mathbf{\mathbf{k-}}\left(\mathbf{a}\cdot\mathbf{\mathbf{k}}\right)\mathbf{\mathbf{E}}$.
From the relation $\mathbf{k}\cdot\mathbf{D}=0$, for
$\mathbf{k}=\omega\mathbf{n}$, we obtain

$\left[\epsilon\mathbf{n}+\mathbf{n\times}\left(\mathbf{a\times
n}\right)\right]\cdot\mathbf{\mathbf{E}}=0.$ (128)

In general, the propagating modes are no longer transversal. Yet, when
$\mathbf{a}$ and $\mathbf{n}$ are parallel vectors, Eq. (128) becomes

$\epsilon\mathbf{n}\cdot\mathbf{\mathbf{E}}=0,$ (129)

recovering a transversal electric field.

Replacing Eqs. (122a) and (122b) in Eq. (19), one obtains the following
extended electric permittivity tensor:

$\bar{\epsilon}_{ij}=\left(\epsilon-\frac{\mathbf{a}\cdot\mathbf{k}}{\omega}-\frac{\mathbf{b}\cdot\mathbf{k}}{\omega}\right)\delta_{ij}+\frac{b_{i}k_{j}}{\omega}+\frac{a_{j}k_{i}}{\omega},$
(130)

where there appear direction dependent terms: $(\mathbf{a}\cdot\mathbf{k})$,
$(\mathbf{b}\cdot\mathbf{k})$, $a_{i}k_{j}$, $a_{j}k_{i}$. Now the matrix
$M_{ij}$ is written as

$\displaystyle{M}$
$\displaystyle=\mathcal{N}+\mu[({\bf{a}}+{\bf{b}})\cdot{\bf{n}}]\mathbb{1}_{3\times
3}-\mu\left(\mathcal{A}+\mathcal{B}\right),$ (131)

where $\mathcal{N}$ is given by Eq.(29), and

$\displaystyle\mathcal{A}$
$\displaystyle=\mathrm{diag}\left(\frac{}{}(b_{1}+a_{1})n_{1},(b_{2}+a_{2})n_{2},(b_{3}+a_{3})n_{3}\right),$
(132a) $\displaystyle\mathcal{B}$
$\displaystyle=\begin{pmatrix}0&b_{1}n_{2}+a_{2}n_{1}&b_{1}n_{3}+a_{3}n_{1}\\\
b_{2}n_{1}+a_{1}n_{2}&0&b_{2}n_{3}+a_{3}n_{2}\\\
b_{3}n_{1}+a_{1}n_{3}&b_{3}n_{2}+a_{2}n_{3}&0\end{pmatrix}.$ (132b)

For $\mathrm{det}[M_{ij}]=0$, the following dispersion equation is attained:

$\displaystyle 0$
$\displaystyle=\left[n^{2}-\mu\epsilon+\mu(\mathbf{c}\cdot\mathbf{n})\right]\left\\{\left[n^{2}-\mu\epsilon+\mu(\mathbf{c}\cdot\mathbf{n})\right]\frac{}{}\right.$
$\displaystyle\left.\hskip
28.45274pt+\frac{\mu}{\epsilon}\left[(\mathbf{a}\cdot\mathbf{b})n^{2}-(\mathbf{a}\cdot\mathbf{n})(\mathbf{b}\cdot\mathbf{n})\right]\right\\},$
(133)

which can be cast in the form:

$\left[\frac{{}}{{}}n^{2}-\mu\epsilon+\mu(\mathbf{c}\cdot\mathbf{n})\right]\left[\mathbf{n}^{T}\mathbb{\tilde{Q}}\mathbf{n}+\mu(\mathbf{c}\cdot\mathbf{n})-\mu\epsilon\right]=0,$
(134)

where $\mathbf{c}=\mathbf{a}+\mathbf{b}$ is always a real vector [due to Eq.
(123)] and $\mathbb{Q}$ is a complex $\left(3\times 3\right)$ self-adjoint
matrix defined by

$\mathbb{Q}=1+\frac{\mu\left(\mathbf{a}\cdot\mathbf{b}\right)}{\epsilon}-\frac{\mu\left(\mathbf{ab}^{T}+\mathbf{ba}^{T}\right)}{2\epsilon}.$
(135)

The first dispersion relation, expressed as

$\left(\mathbf{n}+\frac{\mu{\mathbf{c}}}{2}\right)^{2}=\mu\epsilon+\frac{\mu^{2}{\mathbf{c}}^{2}}{4},$
(136)

describes a sphere centered in $\mathbf{n}_{0}=-\mu{\mathbf{c}}/2$ with radius
$\sqrt{\mu\epsilon+\mu^{2}{\mathbf{c}}^{2}/4}$. The second dispersion relation
is

$\mathbf{n}^{T}\mathbb{Q}\mathbf{n}+\mu(\mathbf{a}+\mathbf{b})\cdot\mathbf{n}=\mu\epsilon.$
(137)

The eigenvalues of $\mathbb{Q}$ determine the surface described by the
dispersion relation. If all eigenvalues are positive, the surface becomes an
ellipsoid whose center is not at the origin [if
$(\mathbf{a}+\mathbf{b)}\cdot\mathbf{n}=0$, it is centered at the origin],
with the principal axes oriented along the respective eigenvectors. When there
are at least two distinct eigenvalues, the medium produces birefringence.

###  V.1 Propagation properties of magnetoelectric parameters in dielectrics

In order to investigate the electromagnetic propagation in a dielectric medium
governed by the constitutive relations (124) and (125), we suppose that the
magnetoelectric parameters are constrained by relation (123). In the case the
vectors $\mathbf{a}$ and $\mathbf{b}$ have a complex piece, that is

$\mathbf{a}=\mathbf{a}^{\prime}+\mathrm{i}\mathbf{a}^{\prime\prime},\quad\mathbf{b}=\mathbf{b}^{\prime}+\mathrm{i}\mathbf{b}^{\prime\prime},$
(138)

one writes $\mathbf{a}^{\prime}=\mathbf{b}^{\prime}$ and
$\mathbf{a}^{\prime\prime}=-\mathbf{b}^{\prime\prime}$, as a consequence of
Eq. (123). In this case, we have

$\displaystyle(\mathbf{a}+\mathbf{b})\cdot\mathbf{n}=2(\mathbf{a}^{\prime}\cdot\mathbf{n}),$
(139)
$\displaystyle(\mathbf{a}\cdot\mathbf{n})(\mathbf{b}\cdot\mathbf{n})=(\mathbf{{a^{\prime}}}\cdot\mathbf{n})^{2}+(\mathbf{{a^{\prime\prime}}}\cdot\mathbf{n})^{2},$
(140)
$\displaystyle\mathbf{a}\cdot\mathbf{b}={a^{\prime}}^{2}+{a^{\prime\prime}}^{2}=|a|^{2},$
(141)

for a real $\mathbf{n}$ vector. We also suppose that the 3-vector $\mathbf{a}$
fulfills $\mathbf{a}^{\prime}\cdot\mathbf{n}=a^{\prime}n\cos\varphi$,
$\mathbf{a}^{\prime\prime}\cdot\mathbf{n}=a^{\prime\prime}n\cos\varphi$, so
that the involved relation (133) provides,

$\displaystyle\left[n^{2}-\mu\epsilon+2\mu a^{\prime}n\cos\varphi\right]=0$
(142a)
$\displaystyle\left[n^{2}\left(\epsilon+\mu|a|^{2}\sin^{2}\varphi\right)+2\mu\epsilon
a^{\prime}n\cos\varphi-\mu\epsilon^{2}\right]=0,$ (142b)

from which the following (positive) indices are achieved:

$\displaystyle n_{(1)}$
$\displaystyle=\mu\sqrt{{a^{\prime}}^{2}\cos^{2}\varphi+\epsilon/\mu}-\mu
a^{\prime}\cos\varphi,$ (143) $\displaystyle n_{(2)}$
$\displaystyle=\frac{1}{r}\left(\sqrt{\mu\epsilon+\mu^{2}{a^{\prime}}^{2}+\mu^{2}{a^{\prime\prime}}^{2}\sin^{2}\varphi}-\mu
a^{\prime}\cos\varphi\right),$ (144)

where

$\displaystyle r$
$\displaystyle=1+\frac{\mu}{\epsilon}|a|^{2}\sin^{2}\varphi.$ (145)

Above, we have retained only the roots corresponding to positive refractive
indices, since we are not addressing metamaterials. The general behavior of
$n_{(1,2)}$ is illustrated in Figs. 5 and 6 in terms of $\varphi\in[0,\pi]$
and the dimensionless parameter $a^{\prime}\in[0,1]$.

Figure 5: Refractive index $n_{(1)}$ of Eq. (143) in terms of $\varphi$ and
$a^{\prime}$. Here we have set $\mu=1$ and $\epsilon=2$. The parameters
$\epsilon$, $\mu$, and $a^{\prime}$ are dimensionless. Figure 6: Refractive
index $n_{(2)}$ of Eq. (144) in terms of $\varphi$ and $a^{\prime}$. We have
used $\mu=1$, $\epsilon=2$ and $a^{\prime\prime}=1$. The parameters
$\epsilon$, $\mu$, and $a^{\prime}$ are dimensionless.

The anisotropy effect is described by the angular dependence on $\varphi$, not
exactly equal for $n_{(1)}$ and $n_{(2)}$. In fact, note that by setting
$a^{\prime}\mapsto 0$, one obtains $n_{1}\mapsto\sqrt{\mu\epsilon}$, a
constant value which corresponds to the straight border line of Fig. 5 (for
$a^{\prime}=0$). Under such a limit,
$n_{2}\mapsto\mu\epsilon\left(\mu\epsilon+\mu^{2}a^{\prime\prime
2}\sin^{2}\varphi\right)^{-1/2}$, which is represented by the sinusoidal
border line at $a^{\prime}=0$ in Fig. 6. Furthermore, the border behavior at
$a^{\prime}=1$ is also distinct for the two plots. Both $n_{(1)}$ and
$n_{(2)}$ have maximal values at $\varphi=\pi$.

It is important to mention that both Eq. (133) as the refractive indices (143)
and (144) hold for any propagation direction in relation to the vector
$\mathbf{a}$, parametrized by the angle $\varphi$. Now, we investigate
propagation for some special angles between $\mathbf{a}$ and $\mathbf{n}$.

####  V.1.1 Particular case 1: $\mathbf{a}$-orthogonal configuration.

We begin bychoosing a scenario where the vector $\mathbf{a}$ is orthogonal to
the propagation axis, that is, $\mathbf{a}\cdot\mathbf{n}=0$ (or
$\varphi=\pi/2$). Thus, the dispersion relation (142) provides two solutions
for the refractive index,

$\displaystyle n_{(1)}$ $\displaystyle=\sqrt{\mu\epsilon},$ (146)
$\displaystyle n_{(2)}$
$\displaystyle=\frac{\sqrt{\mu\epsilon}}{\sqrt{1+(\mu/\epsilon){|a|}^{2}}}.$
(147)

The two refractive indices are real. Therefore, there occurs birefringence but
not absorption. To obtain the propagating modes, we take the vector
$\mathbf{n}$ along the $z$-axis, that is, $\mathbf{n}=(0,0,n_{3})$. Orthogonal
to it, we set $\mathbf{{a^{\prime}}}=(a_{1},a_{2},0)$,
$\mathbf{{a^{\prime\prime}}}=(a_{1},a_{2},0)$. This choice leads to a very
simple expression for the matrix of Eq. (131), that is,

${M}=\begin{pmatrix}n_{3}^{2}-\mu\epsilon&0&-\mu a_{1}^{*}n_{3}\\\
0&n_{3}^{2}-\mu\epsilon&-\mu{a_{2}^{*}}n_{3}\\\
-\mu{a_{1}}n_{3}&-\mu{a_{2}}n_{3}&-\mu\epsilon\end{pmatrix}.$ (148)

Replacing index (146) in matrix (148), the condition $M_{ij}E^{j}=0$ provides
as a solution a transversal mode,

$\mathbf{E}_{(1)}=\frac{1}{|\mathbf{a}|}\begin{pmatrix}{a_{2}}\\\ -{a_{1}}\\\
0\end{pmatrix},$ (149)

with $|\mathbf{a}|=\sqrt{|{a_{1}}|^{2}+|{a_{2}}|^{2}}$. Now, replacing
relation (147) in matrix (148), one has

$\begin{pmatrix}-\mu\epsilon f_{a}&0&-\mu{a_{1}^{*}}n_{3}\\\ 0&-\mu\epsilon
f_{a}&-\mu{a_{2}^{*}}n_{3}\\\
-\mu{a_{1}}n_{3}&-\mu{a_{2}}n_{3}&-\mu\epsilon\end{pmatrix}\begin{pmatrix}E_{x}\\\
E_{y}\\\ E_{z}\end{pmatrix}=0,$ (150)

whose solution is

$\mathbf{E}_{(2)}=\frac{1}{|a|\sqrt{1+f_{a}}}\begin{pmatrix}{a_{1}^{*}}\\\
{a_{2}^{*}}\\\ -|a|\sqrt{f_{a}}\end{pmatrix},$ (151)

with

$\displaystyle f_{a}$
$\displaystyle=\frac{(\mu/\epsilon)|a|^{2}}{1+(\mu/\epsilon){|a|^{2}}}.$ (152)

We note that Eq.(151) represents a mixed mode, endowed with a longitudinal
component. In this case, it is not possible to find a pure transversal mode
for the field $\mathbf{E}_{(2)}$. The transversal mode $\mathbf{E}_{(1)}$ and
the transversal sector of the mixed mode $\mathbf{E}_{(2)}$ could exhibit
linear polarization, or circular or elliptical polarization. Indeed, in
principle one writes the $\mathbf{{a}}$-vector components as
$a_{1}=(a_{1}^{\prime}+\mathrm{i}a_{1}^{\prime\prime})$,
$a_{2}=(a_{2}^{\prime}+\mathrm{i}a_{2}^{\prime\prime})$. One notices that for
either $a_{1}^{\prime\prime}=a_{2}^{\prime\prime}=0$ or
$a_{1}^{\prime}=a_{2}^{\prime}=0$, Eq. (149) and the transversal part of Eq.
(151) yield linearly polarized modes. On the other hand, for either
$a_{1}^{\prime}=a_{2}^{\prime\prime}=0$ or
$a_{1}^{\prime\prime}=a_{2}^{\prime}=0$, that is, for $a_{1}=a_{1}^{\prime}$,
$a_{2}=\mathrm{i}a_{2}^{\prime\prime}$, or
$a_{1}=\mathrm{i}a_{1}^{\prime\prime}$, $a_{2}=a_{2}^{\prime}$, the
polarization is elliptical. Circularly polarized modes only occur when either
$a_{1}^{\prime}=a_{2}^{\prime\prime}=0$ and
$a_{2}^{\prime}=a_{1}^{\prime\prime}$ or
$a_{1}^{\prime\prime}=a_{2}^{\prime}=0$ and
$a_{1}^{\prime}=a_{2}^{\prime\prime}$.

After finding the refractive indices (146) and (147), and the corresponding
modes (149) and (151), we need to discuss the physical effects on wave
propagation. In the case the associated modes are linearly or elliptically
polarized, the implied birefringence is expressed in terms of the phase shift
(94), namely,

$\Delta=\frac{2\pi}{\lambda_{0}}l[n_{(1)}-n_{(2)}]\,,$ (153)

where $n_{1}$ and $n_{2}$ are the refractive indices (146) and (147),
respectively. The phase shift per unit length is

$\frac{\Delta}{l}=\frac{2\pi}{\lambda_{0}}\sqrt{\mu\epsilon}\left[1-\frac{1}{\sqrt{1+(\mu/\epsilon)|a|^{2}}}\right]\,,$
(154)

which, for $(\mu/\epsilon)|a|^{2}\ll 1$, simplifies as

$\frac{\Delta}{l}=\frac{\pi\mu|a|^{2}}{\lambda_{0}}.$ (155)

####  V.1.2 Particular case 2: $\mathbf{a}$-longitudinal configuration.

Let us now consider the case where the vectors $\mathbf{a}$ and $\mathbf{n}$
point along the same direction, $\mathbf{a}\cdot\mathbf{n}=an$, for which Eq.
(133) is written as

$\displaystyle\left(n^{2}+2\mu a^{\prime}n-\mu\epsilon\right)^{2}=0,$ (156)

which involves the square of a quadratic polynomial in $n$. Thus, there is a
doubly degenerate refractive index, namely,

$\displaystyle n=\sqrt{\mu\epsilon+\mu^{2}{a^{\prime}}^{2}}-\mu a^{\prime}.$
(157)

The latter corresponds exactly to the solutions of Eqs. (143) and (144) for
$\varphi=0$, as expected.

To obtain the propagating modes, we take the vector $\mathbf{n}$ along the
$z$-axis, that is, $\mathbf{n}=(0,0,n_{3})$, in such a way that
$\mathbf{a}=\mathbf{b^{*}}=(0,0,a_{3}^{\prime}+\mathrm{i}a_{3}^{\prime\prime})$.
It leads to a simple form of the matrix of Eq.(131), that is,

${M}=\begin{pmatrix}n_{3}^{2}-\mu\epsilon+A&0&0\\\
0&n_{3}^{2}-\mu\epsilon+A&0\\\ 0&0&-\mu\epsilon\end{pmatrix},$ (158)

where $A=2\mu na_{3}^{\prime}$. With index (157), matrix (158) reads simply as

${M}=\begin{pmatrix}0&0&0\\\ 0&0&0\\\ 0&0&-\mu\epsilon\end{pmatrix},$ (159)

and provides generic transversal modes,

$\mathbf{E}=\begin{pmatrix}E_{x}\\\ E_{y}\\\ 0\end{pmatrix},$ (160)

with arbitrary $E_{x}$ and $E_{y}$. Note this transversality occurs in
accordance with Eq. (129). The solution (160) may represent a linear,
circular, or elliptic polarization mode, depending on the nature and relation
between $E_{x}$ and $E_{y}$. As this relation is not supplied by the features
examined so far, we conclude that the $\mathbf{a}$-longitudinal configuration
allows any polarization, in principle. Furthermore, only one refractive index
expression was achieved, with no signal of anisotropy. Thus, we conclude that
the ${\bf{a}}$ direction defines the optical axis of the medium.

It is worthwhile to observe that the complex vectors $\mathbf{a}$ and
$\mathbf{b}$, given in Eq. (138), yield real refractive indices as it happens
in the case of real vectors. For the $\mathbf{a}$-orthogonal configuration,
the complex vectors of Eq. (138) may provide circular or elliptical
polarization (in its transversal sector), besides the linear one.

###  V.2 Group velocity, phase velocity, and Poynting vector

From Eq. (133), we can write down two dispersion equations,

$\displaystyle 0$
$\displaystyle=\omega^{2}-2\frac{\omega}{\epsilon}({\bf{a}}^{\prime}\cdot{\bf{k}})-\frac{{\bf{k}}^{2}}{\mu\epsilon},$
(161) $\displaystyle 0$
$\displaystyle=\omega^{2}-2\frac{\omega}{\epsilon}({\bf{a}}^{\prime}\cdot{\bf{k}})-\left(\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}\right){\bf{k}}^{2}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}}+$
$\displaystyle\phantom{=}+\frac{({\bf{a}}^{\prime\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}},$
(162)

which provide, respectively, the following solutions for $\omega$ (with
$\omega>0$):

$\displaystyle\omega_{(1)}$
$\displaystyle=\sqrt{\frac{{\bf{k}}^{2}}{\mu\epsilon}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}}}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})}{\epsilon},$
(163) $\displaystyle\omega_{(2)}$
$\displaystyle=\sqrt{\left(\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}\right){\bf{k}}^{2}-\frac{({\bf{a}}^{\prime\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}}}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})}{\epsilon}.$
(164)

The phase and group velocities are given by

$\displaystyle v_{\mathrm{ph}(1)}$
$\displaystyle=\sqrt{\frac{1}{\mu\epsilon}+\frac{({\bf{a}}^{\prime}\cdot\hat{\bf{k}})^{2}}{\epsilon^{2}}}+\frac{({\bf{a}}^{\prime}\cdot\hat{\bf{k}})}{\epsilon},$
(165) $\displaystyle v_{\mathrm{ph}(2)}$
$\displaystyle=\sqrt{\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}-\frac{({\bf{a}}^{\prime\prime}\cdot\hat{\bf{k}})^{2}}{\epsilon^{2}}}+\frac{({\bf{a}}^{\prime}\cdot\hat{\bf{k}})}{\epsilon},$
(166) $\displaystyle{\bf{v}}_{\mathrm{g}(1)}$
$\displaystyle=f_{a^{\prime}}\left[\frac{{\bf{k}}}{\mu\epsilon}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})}{\epsilon^{2}}{\bf{a}}^{\prime}\right]+\frac{{\bf{a}}^{\prime}}{\epsilon},$
(167) $\displaystyle{\bf{v}}_{\mathrm{g}(2)}$
$\displaystyle=f_{a^{\prime\prime}}\left[\left(\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}\right){\bf{k}}-\frac{({\bf{a}}^{\prime\prime}\cdot{\bf{k}})}{\epsilon^{2}}{\bf{a}}^{\prime\prime}\right]+\frac{{\bf{a}}^{\prime}}{\epsilon},$
(168)

with

$\displaystyle f_{a^{\prime}}$
$\displaystyle=\left[\frac{{\bf{k}}^{2}}{\mu\epsilon}+\frac{({\bf{a}}^{\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}}\right]^{-1/2},$
(169) $\displaystyle f_{a^{\prime\prime}}$
$\displaystyle=\left[\left(\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}\right){\bf{k}}^{2}-\frac{({\bf{a}}^{\prime\prime}\cdot{\bf{k}})^{2}}{\epsilon^{2}}\right]^{-1/2}.$
(170)

Relations (167) and (168) reveal that the group velocity is in general not
parallel to ${\bf{n}}$. Starting from Eqs. (163) and (164), we now
particularize the phase and group velocities for propagation orthogonal and
longitudinal to ${\bf{a}}^{\prime}$.

  * •

For the ${\bf{a}}$-orthogonal case where
${\bf{a}}^{\prime}\cdot{\bf{k}}={\bf{a}}^{\prime\prime}\cdot{\bf{k}}=0$, one
finds

$\displaystyle v_{\mathrm{ph}(1)}^{\mathrm{orth}}$
$\displaystyle=\sqrt{\frac{1}{\mu\epsilon}},\quad\mathbf{v}_{\mathrm{g}(1)}^{\mathrm{orth}}=\hat{\bf{k}}\sqrt{\frac{1}{\mu\epsilon}},$
(171) $\displaystyle v_{\mathrm{ph}(2)}^{\mathrm{orth}}$
$\displaystyle=\sqrt{\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}},\quad{\bf{v}}_{\mathrm{g}(2)}^{\mathrm{orth}}=\hat{\bf{k}}\sqrt{\frac{1}{\mu\epsilon}+\frac{|a|^{2}}{\epsilon^{2}}}.$
(172)

  * •

For the ${\bf{a}}$-longitudinal scenario
${\bf{a}}^{\prime,\prime\prime}\cdot{\bf{k}}=a^{\prime,\prime\prime}k$, we
obtain

$\displaystyle v_{\mathrm{ph}(1)}^{\mathrm{long}}$
$\displaystyle=v_{\mathrm{g}(1)}^{\mathrm{long}}=\sqrt{\frac{1}{\mu\epsilon}+\frac{a^{\prime
2}}{\epsilon^{2}}}+\frac{a^{\prime}}{\epsilon},$ (173) $\displaystyle
v_{\mathrm{ph}(2)}^{\mathrm{long}}$
$\displaystyle=v_{\mathrm{g}(2)}^{\mathrm{long}}=\sqrt{\frac{1}{\mu\epsilon}+\frac{a^{\prime
2}}{\epsilon^{2}}}+\frac{a^{\prime}}{\epsilon},$ (174)

with

$\quad{\bf{v}}_{\mathrm{g}(1)}^{\mathrm{orth}}=v_{\mathrm{g}(1)}^{\mathrm{orth}}\hat{\bf{k}},\quad{\bf{v}}_{\mathrm{g}(2)}^{\mathrm{orth}}=v_{\mathrm{g}(2)}^{\mathrm{orth}}\hat{\bf{k}}.$
(175)

In Fig. 7 we plot the group velocities $v_{\mathrm{g}(1,2)}$ of both
${\bf{a}}$-longitudinal and ${\bf{a}}$-orthogonal cases in terms of the
parameter $a^{\prime}$.

Figure 7:  Group velocities for ${\bf{a}}$-longitudinal and
${\bf{a}}$-orthogonal cases. The blue solid line indicates
$v_{\mathrm{g}(1,2)}^{\mathrm{long.}}$ of Eqs. (173) and (174). The red curve
illustrates $v_{\mathrm{g}(1)}^{\mathrm{orth.}}$ of Eq. (171), and the green
line represents $v_{\mathrm{g}(2)}^{\mathrm{orth.}}$ of Eq. (172). The gray
dashed vertical lines indicate the values of $a^{\prime}\in\\{1/2,1\\}$ above
which the group velocities become greater than 1, in agreement with Eqs.
(176), and (177), respectively. Here we have used $\mu=1$, $\epsilon=2$ and
$a^{\prime\prime}=1$.

We notice the group velocities can be greater than 1 depending on the values
of $a^{\prime}$.

  * •

For the ${\bf{a}}$-longitudinal case, in order to ensure
$v_{\mathrm{g}}^{\mathrm{long}}<1$, the following must occur:

$\displaystyle a^{\prime}<\frac{\mu\epsilon-1}{2\mu}.$ (176)

  * •

For the ${\bf{a}}$-orthogonal configuration, one finds that
$v_{\mathrm{g}(2)}^{\mathrm{orth}}<1$ only when

$\displaystyle|a|<\sqrt{\frac{\epsilon}{\mu}(\mu\epsilon-1)}.$ (177)

These results are important to constrain the magnetoelectric parameters in
ranges suitable to ensure group velocity smaller than 1.

With the constitutive relation (125), the Poynting vector is

$\displaystyle{\bf{S}}$
$\displaystyle=\frac{1}{2\mu}({\bf{n}}+\mu{\bf{a}})|{\bf{E}}|^{2}-\frac{1}{2\mu}({\bf{n}}\cdot{\bf{E}}){\bf{E}}^{*}-\frac{1}{2}({\bf{a}}\cdot{\bf{E}}){\bf{E}}^{*}.$
(178)

By considering ${\bf{a}}={\bf{a}}^{\prime}+\mathrm{i}{\bf{a}}^{\prime\prime}$
and ${\bf{E}}={\bf{E}}^{\prime}+\mathrm{i}{\bf{E}}^{\prime\prime}$, we firstly
analyze the case $\bf{a}^{\prime},\bf{a}^{\prime\prime}$ parallel to
$\mathbf{n}$. From Eq. (128) we obtain $\mathbf{n}\cdot{\bf{E}}=0$, i.e.,
$\mathbf{n}\cdot\mathbf{E}^{\prime}=0=\mathbf{n}\cdot\mathbf{E}^{\prime\prime}$
implying
$\mathbf{a}^{\prime}\cdot\mathbf{E}=0=\mathbf{a}^{\prime\prime}\cdot\mathbf{E}$
or $\mathbf{a}\cdot\mathbf{E}=0$. This way Eq. (178) simplifies to

$\displaystyle{\bf{S}}^{\mathrm{long.}}$
$\displaystyle=\frac{1}{2\mu}({\bf{n}}+\mu{\bf{a}})|{\bf{E}}|^{2},$ (179)

whose time-averaged form becomes

$\displaystyle\left\langle{\bf{S}}\right\rangle^{\mathrm{long.}}$
$\displaystyle=\frac{1}{2\mu}({\bf{n}}+\mu{\bf{a}^{\prime}})|{\bf{E}}|^{2}.$
(180)

As ${\bf{a}^{\prime}}\parallel{\bf{n}}$, the energy flux propagates along the
wave propagation direction in this situation.

For a general case, the time-averaged Poynting vector is

$\displaystyle\left\langle{\bf{S}}\right\rangle$
$\displaystyle=\frac{1}{2\mu}({\bf{n}}+\mu{\bf{a}}^{\prime})|{\bf{E}}|^{2}+$
$\displaystyle\phantom{=}+\frac{1}{2}(F_{1}{\bf{E}}^{\prime}-F_{2}{\bf{E}}^{\prime\prime})\left[({\bf{a}}^{\prime}\cdot{\bf{E}}^{\prime})+({\bf{a}}^{\prime\prime}\cdot{\bf{E}}^{\prime\prime})\right]+$
$\displaystyle\phantom{=}+\frac{1}{2}(F_{1}{\bf{E}}^{\prime\prime}+F_{2}{\bf{E}}^{\prime})\left[({\bf{a}}^{\prime}\cdot{\bf{E}}^{\prime\prime})-({\bf{a}}^{\prime\prime}\cdot{\bf{E}}^{\prime})\right],$
(181)

where

$\displaystyle
F_{1}=\frac{n^{2}}{\mu}\frac{\left[\epsilon-\left(\mathbf{n\cdot
a}^{\prime}\right)\right]}{\Delta}-1,\quad
F_{2}=\frac{n^{2}}{\mu}\frac{\left(\mathbf{n\cdot
a}^{\prime\prime}\right)}{\Delta},$ (182)

with

$\displaystyle\Delta=\left[\epsilon-\left(\mathbf{n\cdot
a}^{\prime}\right)\right]^{2}+\left(\mathbf{n\cdot
a}^{\prime\prime}\right)^{2}.$ (183)

In particular, by considering $\bf{a}^{\prime},\bf{a}^{\prime\prime}$
orthogonal to $\mathbf{n}$, we obtain

$\displaystyle\left\langle{\bf{S}}\right\rangle^{\mathrm{orth.}}$
$\displaystyle=\frac{1}{2\mu}({\bf{n}}+\mu{\bf{a}}^{\prime})|{\bf{E}}|^{2}+\frac{(n^{2}-\mu\epsilon)}{2\mu\epsilon}\left[\frac{}{}({\bf{a}}^{\prime}\cdot{\bf{E}}^{\prime}){\bf{E}}^{\prime}+\right.$
$\displaystyle\left.\phantom{=}+({\bf{a}}^{\prime\prime}\cdot{\bf{E}}^{\prime\prime}){\bf{E}}^{\prime}+({\bf{a}}^{\prime}\cdot{\bf{E}}^{\prime\prime}){\bf{E}}^{\prime\prime}-({\bf{a}}^{\prime\prime}\cdot{\bf{E}}^{\prime}){\bf{E}}^{\prime\prime}\frac{}{}\right],$
(184)

which yields a propagation not aligned to $\bf{n}$. This result is consistent
with the group velocities (167) and (168), which contain a piece not belonging
to the $\bf{n}$ axis.

Finally, we note that for $\bf{a}^{\prime},\bf{a}^{\prime\prime}$ nonparallel
to ${\bf{n}}$, the energy flow does not propagate along the propagation axis
because there is a contribution along the electric field direction.

##  VI Final Remarks

In this work, we have examined the propagation of electromagnetic waves in bi-
isotropic and bi-anisotropic matter. As for the bi-anisotropic scenarios, we
have taken isotropic electric permittivity and magnetic permeability tensors,
while the magnetoelectric parameters were supposed as symmetric and
antisymmetric anisotropic complex tensors. As an initial action, a modified
permittivity tensor was written as part of the matrix equation, which allowed
us to obtain the dispersion equations and refractive indices.

The symmetric magnetoelectric tensors were parametrized in terms of a single
3-vector, $\mathbf{d}$, and two complex scalars,
$\tilde{\alpha},\tilde{\beta}$. For propagation along the $z$-axis,
$\mathbf{n}=(0,0,n)$, arbitrary transversal propagating modes were obtained
for the $\mathbf{d}$-longitudinal configuration. The $\mathbf{d}$ axis
coincides with the optical axis of the medium. On the other hand, when the
$\mathbf{d}$-vector is orthogonal to the propagation axis, the associated
modes are also transversal, with polarization linear or elliptical.

The antisymmetric magnetoelectric tensors were parametrized in terms of two
3-vectors, $\mathbf{a}$ and $\mathbf{b}$, related by
$\mathbf{b}^{*}=\mathbf{a}$. We have considered complex vectors and obtained
real and positive refractive indices [see Eqs. (143) and (144)]. For the
$\mathbf{a}$-longitudinal configuration, arbitrary transversal propagating
modes, associated with one unique refractive index, were obtained for the
propagation along the $z$ axis, $\mathbf{n}=(0,0,n)$. Therefore, in this case,
the vector $\mathbf{a}$ determines the optical axis of the medium. For the
$\mathbf{a}$-orthogonal configuration, there appears a transversal and a mixed
mode, composed of a longitudinal and a transversal piece, associated with two
refractive indices. The polarization may be linear, elliptical, or circular.
The birefringence effect was evaluated in terms of the $\mathbf{a}$-vector
magnitude.

In both symmetric and antisymmetric cases, the propagation along the
$\mathbf{d}$ or $\mathbf{a}$ directions is isotropic, since they define the
optical axis of the medium. On the other hand, the propagation orthogonal to
the magnetoelectric vectors $\mathbf{d}$ or $\mathbf{a}$ provides a route of
phenomenological distinction between the symmetric [Eq. (67)] or antisymmetric
tensor [Eq. (122)], due to the observed difference between the associated
propagating modes. See Eqs. (93) and (151). We have also evaluated the group
velocities for the bi-isotropic and bi-anisotropic cases examined. In the bi-
isotropic case, the phase and group velocities turned out equal. For the bi-
anisotropic configurations, these velocities are different. The general group
velocities present a component along the 3-vectors $\mathbf{d}$ or
$\mathbf{a}$, being no longer necessarily parallel to $\mathbf{n}$ in these
cases. We have also carried out the Poynting vector, observing that the
electromagnetic energy flux does not occur along the wave propagation
direction for general configurations. The energy flux direction coincides with
the $\mathbf{n}$ axis in the following situations: when $\mathbf{d}$ is
parallel or orthogonal to $\mathbf{n}$ (in the symmetric case) and when
$\mathbf{a}^{\prime}$ and $\mathbf{a}^{\prime\prime}$ are parallel to
$\mathbf{n}$ (in the antisymmetric case).

As a final remark, we may try to state a parallel between our results and the
physics of anisotropic media described by constitutive relations (3), with
permittivity and permeability described by general tensors $\epsilon_{ij}$ and
$\mu_{ij}$. This is the case examined in Refs. Yakov ; Markel . The general
anisotropic dispersion relations of these references are different from our
relations, which are based on isotropic permittivity and permeability,
$\epsilon_{ij}=\epsilon\delta_{ij}$ and $\mu_{ij}=\mu\delta_{ij}$, and non-
null magnetoelectric parameters, $\alpha_{ij}$ and $\beta_{ij}$. These basic
distinctions, consequently, do not favor straightforward comparisons between
such anisotropic systems.

### Acknowledgments

The authors P.D.S.S., R.C., and M.M.F. express their gratitude to FAPEMA,
CNPq, and CAPES (Brazilian research agencies) for invaluable financial
support. R.C. acknowledges the support from the Grants No. CNPq/306724/2019-7,
No. FAPEMA/Universal-01131/17, and No. FAPEMA/Universal-00812/19. M.M.F. is
supported by No. FAPEMA/Universal/01187/18, No. CNPq/Produtividade
311220/2019-3, and No. CNPq/Universal/422527/2021-1. Furthermore, we are
indebted to CAPES/Finance Code 001 and FAPEMA/POS-GRAD-02575/21.

##  Appendix A  Vector evaluation for the electric field of the modes

In this Appendix, we present an alternative route for carrying out the
electric field of the propagating modes in the symmetric and antisymmetric
constitutive parameters.

###  A.1 Symmetric parameters case

For non-parallel vectors $\mathbf{d}$ and $\mathbf{n}$, the electric field may
be written as

$\displaystyle\mathbf{E}$ $\displaystyle=$
$\displaystyle\frac{E_{1}}{\left|\mathbf{d}\times\mathbf{n}\right|}\mathbf{d}\times\mathbf{n}+\frac{E_{2}}{n\left|\mathbf{d}\times\mathbf{n}\right|}\mathbf{n\times}\left(\mathbf{d}\times\mathbf{n}\right)$
(185) $\displaystyle-\frac{E_{1}}{\epsilon
n^{2}}\tilde{\alpha}(\mathbf{n}\cdot\mathbf{d})\left|\mathbf{d}\times\mathbf{n}\right|\mathbf{n,}$

which does not supply, in general, transversal modes, that is,
$\mathbf{n}\cdot\mathbf{E}\neq 0$. However, for orthogonal vectors,
$\mathbf{n}\cdot\mathbf{d=0}$, Eq. (72) yields transversal modes,
$\mathbf{n}\cdot\mathbf{E}=0$, whose electric field can be written as

$\mathbf{E}={E_{1}}(\hat{\mathbf{d}}\times\hat{\mathbf{n}})+{E_{2}}\hat{\mathbf{d}},$
(186)

where $\hat{\mathbf{d}}={\mathbf{d}}/|{\mathbf{d}}|$ and
$\hat{\mathbf{n}}={\mathbf{n}}/|{\mathbf{n}}|$. In general, $E_{1}$ and
$E_{2}$ are arbitrary constants, and for a normalized electric field, they
satisfy $|E_{1}|^{2}+|E_{2}|^{2}=1$.

For $\mathbf{n}$ and $\mathbf{d}$ nonparallel, the equation
$\mathbb{M}\mathbf{E}=0$, with the electric field given in Eq. (185), sets

${E}_{2}=E_{1}\frac{\mu\tilde{\alpha}(d^{2}\sin^{2}\varphi)n_{\pm}}{n_{\pm}^{2}-\mu\epsilon}=E_{1}\mu\tilde{\alpha}(d^{2}\sin^{2}\varphi)\Upsilon_{\pm}n_{\pm},$
(187)

where we have used Eqs. (97) and (101), and $n_{\pm}$ is given by Eq. (79). So
the electric field [Eq. (185)] reads

$\displaystyle\mathbf{E}_{\pm}$ $\displaystyle=$ $\displaystyle
E_{0}\frac{\mathbf{d}\times\mathbf{n}_{\pm}}{\left|\mathbf{d}\times\mathbf{n}_{\pm}\right|}+E_{0}\frac{\mu\tilde{\alpha}(d^{2}\sin^{2}\varphi)\Upsilon_{\pm}}{\left|\mathbf{d}\times\mathbf{n}_{\pm}\right|}\mathbf{n_{\pm}\times}\left(\mathbf{d}\times\mathbf{n}_{\pm}\right)$
(188) $\displaystyle-
E_{0}\frac{\tilde{\alpha}d^{2}\cos\varphi\sin\varphi}{\epsilon}\mathbf{n}_{\pm},$

setting $E_{1}=E_{0}$. It provides the same conclusions obtained in Sec.
IV.1.3 for the general $\mathbf{d}$ configuration.

For $\mathbf{n}$ and $\mathbf{d}$ orthogonal vectors, Eq. (187) simplifies as

$E_{2}=E_{1}\frac{\mu\tilde{\alpha}d^{2}n_{\pm}}{n^{2}_{\pm}-\mu\epsilon}=\pm\frac{\tilde{\alpha}}{|\tilde{\alpha}|}E_{1}.$
(189)

where $n_{\pm}$ is now given by Eq. (88). So the electric field [Eq. (188)]
becomes written as

$\mathbf{E}_{\pm}=E_{0}\frac{\mathbf{d}\times\mathbf{n}_{\pm}}{|\mathbf{d}\times\mathbf{n}_{\pm}|}\pm
E_{0}\frac{\tilde{\alpha}}{|\tilde{\alpha}|}\frac{\mathbf{d}}{d},$ (190)

providing the same conclusions obtained in the particular case of Sec. IV.1.2.

The vector formalism can also be used to alternatively express the
permittivity tensor and the dispersion relation. In fact, replacing relations
(67) in the permittivity tensor (19), one writes

$\hat{\bar{\epsilon}}=\epsilon\mathbb{1}{-\tilde{\alpha}\mathbf{d}\left(\mathbf{d}\times\mathbf{n}\right)^{T}+\tilde{\beta}\left(\mathbf{d}\times\mathbf{n}\right)\mathbf{d}^{T}},$
(191)

with $\hat{\bar{\epsilon}}=[{\bar{\epsilon}}_{ij}]$. In such a way the tensor
$\mathbb{M}=[M_{ij}]$ [Eq. (22)], takes the form

$\mathbb{M}=(n^{2}-\mu\epsilon)\mathbb{1}-\mathbf{n}\mathbf{n}^{T}{+\mu\tilde{\alpha}\mathbf{d}(\mathbf{d}\times\mathbf{n})^{T}-\mu\tilde{\beta}(\mathbf{d}\times\mathbf{n})\mathbf{d}^{T}.}$
(192)

Evaluating $\mathrm{det}\mathbb{M}=0$, we read the dispersion relation (77) in
the form,

${\epsilon}\left(n^{2}-\mu\epsilon\right)^{2}+\tilde{\alpha}\tilde{\beta}\mu\left[\mu{\epsilon}d^{2}-(\mathbf{n}\cdot\mathbf{d})^{2}\right]\left|\mathbf{d}\times\mathbf{n}\right|^{2}=0.$
(193)

###  A.2 Antisymmetric parameters case

For nonparallel $\mathbf{a}$ and $\mathbf{n}$, the electric field is expressed
as

$\mathbf{\mathbf{E}}=\frac{E_{1}}{\left|\mathbf{a\times
n}\right|}\mathbf{\left(\mathbf{a\times n}\right)-}\frac{\epsilon
E_{3}}{n\left|\mathbf{a\times
n}\right|^{2}}\mathbf{\mathbf{n\times}\left(\mathbf{a\times
n}\right)}+\frac{E_{3}}{n}\mathbf{n}.$ (194)

The equation $\mathbb{M}\mathbf{E}=0$ with the matrix (131) and the electric
field given in Eq. (194) provides the equations

$E_{1}\left[\frac{{}}{{}}n^{2}-\mu\epsilon+2\mu\mathbf{a}\cdot\mathbf{n}\right]=0,$
(195)
$E_{3}\left[\epsilon\left(n^{2}-\mu\epsilon+2\mu\mathbf{a}\cdot\mathbf{n}\right)+\mu\left|\mathbf{a\times
n}\right|^{2}\right]=0.$ (196)

where the expressions in the brackets are the dispersion relations obtained in
Eq. (133) when $\mathbf{b}=\mathbf{a}$.

If the dispersion relation in Eq. (195),
$n^{2}-\mu\epsilon+2\mu\mathbf{a}\cdot\mathbf{n}=0$, is valid, then $E_{1}$
remains arbitrary and $E_{3}=0$. Thus, the electric field mode (194) becomes
transversal to the propagation direction, that is,

$\mathbf{\mathbf{E}}=\frac{E_{1}}{\left|\mathbf{a\times
n}\right|}(\mathbf{a}\times\mathbf{n}),$ (197)

being in accordance with the result (149) of the case analyzed in Sec. V.1.1.
On the other hand, if the dispersion relation in Eq. (196) holds,
$\epsilon\left(n^{2}-\mu\epsilon+2\mu\mathbf{a}\cdot\mathbf{n}\right)+\mu\left|\mathbf{a\times
n}\right|^{2}=0$, then $E_{1}=0$ and $E_{3}$ remains arbitrary, so that the
electric field modes (194) are no longer transverse,

$\mathbf{\mathbf{E}}=\mathbf{-}\frac{\epsilon E_{3}}{n\left|\mathbf{a\times
n}\right|^{2}}\mathbf{\mathbf{n\times}\left(\mathbf{a\times
n}\right)}+\frac{E_{3}}{n}\mathbf{n}.$ (198)

It is in accordance with the result (151) of the case analyzed in Sec. V.1.1.

Replacing Eqs. (122a) and (122b) in Eq. (19), one obtains the following
extended electric permittivity tensor:

$\bar{\epsilon}=\left(\epsilon-\mathbf{c}\cdot\mathbf{n}\right)\mathbb{1}+\mathbf{b}\mathbf{n}^{T}+\mathbf{n}\mathbf{a}^{T},$
(199)

where $\mathbf{c}=\mathbf{a}+\mathbf{b}$ is a real vector [see Eq. (123)]. All
the pieces of the new electric permittivity coming from antisymmetric
parameters contain direction-dependent terms: $(\mathbf{c}\cdot\mathbf{n})$,
$b_{i}n_{j}$, $n_{i}a_{j}$. Now, the explicit form of the matrix
$\mathbb{M}=[M_{ij}]$ is

$\mathbb{M}=[n^{2}-\mu\epsilon+\mu(\mathbf{c}\cdot\mathbf{n})]\mathbb{1}-\mathbf{n}\mathbf{n}^{T}-\mu\mathbf{b}\mathbf{n}^{T}-\mu\mathbf{n}\mathbf{a}^{T}.$
(200)

##  Appendix B  Rotatory power and dichroism coefficient

Relation (54) holds for the situation that the propagating modes are given by
circularly polarized waves (LCP and RCP). To examine the effect of the optical
activity on the modes, we start from a linearly polarized wave propagating
through a medium along the $z$ axis. As is well known, a wave with linear
polarization,

$\mathbf{E}_{i}=\mathbf{E}_{0i}\mathrm{e}^{\mathrm{i}(kz-\omega t)}\,,$ (201a)
can be split into two circularly polarized waves,
$\displaystyle\mathbf{E}_{0i}=\begin{pmatrix}1\\\ 0\\\
0\end{pmatrix}=\frac{1}{2}\begin{pmatrix}1\\\ -\mathrm{i}\\\
0\end{pmatrix}+\frac{1}{2}\begin{pmatrix}1\\\ \mathrm{i}\\\ 0\end{pmatrix}\,,$
(201b) corresponding to the sum of RCP and LCP waves, respectively. After the
initial wave passes through a distance $z$ in the medium, the final electric
field can be obtained as the combination of two components, $\mathbf{E}_{+}$
and $\mathbf{E}_{-}$, with the wave vectors $\mathbf{k}_{+}$ and
$\mathbf{k}_{-}$, respectively. One then has  $\displaystyle\mathbf{E}_{f}$
$\displaystyle=\mathbf{E}_{+}\mathrm{e}^{\mathrm{i}(k_{+}z-\omega
t)}+\mathbf{E}_{-}\mathrm{e}^{\mathrm{i}(k_{-}z-\omega t)}$
$\displaystyle=\frac{1}{2}\begin{pmatrix}1\\\ \mathrm{i}\\\
0\end{pmatrix}\mathrm{e}^{\mathrm{i}k_{+}z}\mathrm{e}^{-\mathrm{i}\omega
t}+\frac{1}{2}\begin{pmatrix}1\\\ -\mathrm{i}\\\
0\end{pmatrix}\mathrm{e}^{\mathrm{i}k_{-}z}\mathrm{e}^{-\mathrm{i}\omega
t}\,,$ (202)

which can be cast into the form

$\displaystyle\mathbf{E}_{f}$
$\displaystyle=\frac{1}{2}\mathrm{e}^{\mathrm{i}\psi}\mathrm{e}^{-\mathrm{i}\omega
t}\left[\mathrm{e}^{-\mathrm{i}\theta}\begin{pmatrix}1\\\ \mathrm{i}\\\
0\end{pmatrix}+\mathrm{e}^{\mathrm{i}\theta}\begin{pmatrix}1\\\ -\mathrm{i}\\\
0\end{pmatrix}\right]$
$\displaystyle=\mathrm{e}^{\mathrm{i}\psi}\mathrm{e}^{-\mathrm{i}\omega
t}\begin{pmatrix}\cos\theta\\\ \sin\theta\\\ 0\end{pmatrix}\,,$ (203a) with
the quantities $\displaystyle\theta$
$\displaystyle=-\frac{(k_{+}-k_{-})z}{2}\,,$ (203b) $\displaystyle\psi$
$\displaystyle=\frac{(k_{+}+k_{-})z}{2}\,.$ (203c) Notice that Eq. (203a)
describes a linearly polarized wave whose polarization vector is rotated by an
angle $\theta$. From Eq. (203b), one obtains
$\theta=-\frac{(n_{+}-n_{-})z\omega}{2}\,,$ (204)

where we have used $\mathbf{k}=\omega\mathbf{n}$. In general, the refractive
indices can be complex quantities. Because of this, one can infer from Eq.
(204)

$\frac{\theta}{z}=-\frac{\omega}{2}\left[\mathrm{Re}(n_{+})+\mathrm{i}\mathrm{Im}(n_{+})-\mathrm{Re}(n_{-})-\mathrm{i}\mathrm{Im}(n_{-})\right]\,,$
(205)

from which we define the specific rotatory power,

$\delta=-\frac{[\mathrm{Re}(n_{+})-\mathrm{Re}(n_{-})]\omega}{2}\,,$ (206)

as well as the dichroism coefficient,

$\delta_{\mathrm{d}}=-\frac{[\mathrm{Im}(n_{+})-\mathrm{Im}(n_{-})]\omega}{2}\,.$
(207)

Notice that when the medium is nonbirefringent, $\theta=0$ and $\psi=kz$.
Then, the form (201a) is recovered from Eq.(203a).

##  Appendix C  Relations for constitutive parameters

The energy conservation in electromagnetic systems is established by the
Poynting theorem, presented in Eq. (23e) as

$\mathrm{Re}\left[\mathrm{i}\omega\left(\mathbf{E}\cdot\mathbf{D}^{\ast}-\mathbf{H}^{\ast}\cdot\mathbf{B}\right)\right]=0,$
(208)

or

$\mathrm{Re}\left[\mathrm{i}\omega\left(\mathbf{D}^{{\dagger}}\mathbf{E}-\mathbf{{H^{{\dagger}}B}}\right)\right]=0.$
(209)

We can now consider a medium described by the following constitutive relations

$\displaystyle\mathbf{D}$
$\displaystyle=\hat{\epsilon}\mathbf{E}+\hat{\lambda}\mathbf{H},$ (210a)
$\displaystyle\mathbf{B}$
$\displaystyle=\hat{\mu}\mathbf{H}+\hat{\gamma}\mathbf{E}.$ (210b)

with $\hat{\tilde{\epsilon}}$, $\hat{\lambda}$, $\hat{\mu}$ and $\hat{\gamma}$
being nonsingular $3\times 3$ complex matrices. Following the same previous
procedure, Eq. (209) yields

$\displaystyle\mathbf{E}^{\dagger}\left(\hat{\epsilon}^{\dagger}-\hat{\epsilon}\right)\mathbf{E}+\mathbf{H}^{\dagger}\left(\hat{\mu}^{\dagger}-\hat{\mu}\right)\mathbf{H}$
$\displaystyle+$
$\displaystyle+\mathbf{H}^{\dagger}\left(\hat{\lambda}^{\dagger}-\hat{\gamma}\right)\mathbf{E}-\mathbf{E}^{\dagger}\left(\hat{\lambda}-\hat{\gamma}^{\dagger}\right)\mathbf{H}$
$\displaystyle=0.$ (211)

In this case, for arbitrary fields, in order to be consistent with Eq. (211),
one finds

$\hat{\epsilon}^{\dagger}=\hat{\epsilon},\quad\hat{\mu}^{\dagger}=\hat{\mu},\quad\hat{\lambda}=\hat{\gamma}^{\dagger},$
(212)

conditions compatible with energy conservation.

## References

  * (1) J.D. Jackson, Classical Electrodynamics, 3rd ed. (John Wiley & Sons, New York, 1999).
  * (2) A. Zangwill, Modern Electrodynamics (Cambridge University Press, New York, 2012).
  * (3) L.D. Landau and E.M. Lifshitz, Electrodynamics of Continuous Media, Course of Theoretical Physics, Vol. 8, 2nd ed. (Pergamon Press, New York, 1984).
  * (4) A. K. Bain, Crystal optics: properties and applications (Wiley-VCH Verlag GmbH & Co. KGaA, Germany, 2019).
  * (5) G. R. Fowles, Introduction to modern optics, 2nd ed. (Dover Publications, INC., New York, 1975).
  * (6) A. A. Kurmanov, N. A. Ispulov , A. Qadir, A. Zh Zhumabekov, S. N. Sarymova, and K. R. Dossumbekov, Propagation of electromagnetic waves in stationary anisotropic media, Phys. Scr. 96, 085505 (2021).
  * (7) E. U.Condon, Theories of Optical Rotatory Power, Rev. Mod. Phys. 9, 432 (1937).
  * (8) T. Roth and G. L. J. A. Rikken, Observation of Magnetoelectric Linear Birefringence, Phys. Rev. Lett. 88, 063001 (2002).
  * (9) E. Hecht, Optics, 4nd ed. (Addison Wesley, San Francisco, 2002).
  * (10) I. V. Lindell, A. H. Sihvola, S. A. Tretyakov, and A. J. Viitanen, Electromagnetic Waves in Chiral and Bi-Isotropic Media (Artech House, Boston, 1993).
  * (11) J. A. Kong, Electromagnetic Wave Theory (Wiley, New York, 1986).
  * (12) E. O. Kamenetskii, Energy balance equation for electromagnetic waves in bianisotropic media, Phys. Rev. E 54, 4359 (1996).
  * (13) A. H. Sihvola and I. V. Lindell, Bi-isotropic constitutive relations, Microw. Opt. Technol. Lett., 4 (8), 295-297 (1991).
  * (14) A. H. Sihvola and I. V. Lindell, Properties of bi-isotropic Fresnel reflection coefficients, Optics Communications 89, 1 (1992).
  * (15) S. Ougier, I. Chenerie, A. Sihvola, and A. Priou, Propagation in bi-isotropic media: effect of different formalisms on the propagation analysis, Progress In Electromagnetics Research 09, 19 (1994).
  * (16) J. F. Nieves and P. B. Pal, Third electromagnetic constant of an isotropic medium, Am. J. Phys. 62, 207 (1994).
  * (17) R. C. Gauthier, Bi-anisotropic resonators analyzed using Fourier–Bessel numerical formulation; Sagnac effect application, Optics Communications 435, 413 (2019).
  * (18) G. T. Rado and V. J. Folen, Observation of the Magnetically Induced Magnetoelectric Effect and Evidence for Antiferromagnetic Domains, Phys. Rev. Lett. 7, 310 (1961); D. N. Astrov, Magnetoelectric effect in chromium oxide, Soviet Physics JETP 13, 729 (1961).
  * (19) Y. T. Aladadi and M. A. S. Alkanhal, Classification and characterization of electromagnetic materials, Sci. Rep. 10, 11406 (2020).
  * (20) L. Jelinek, R. Marqués, F. Mesa, and J. D. Baena, Periodic arrangements of chiral scatterers providing negative refractive index bi-isotropic media, Phys. Rev. B 77, 205110 (2008).
  * (21) Ming-Che Chang and Min-Fong Yang, Optical signature of topological insulators, Phys. Rev. B 80, 113304 (2009).
  * (22) Zheng-Wei Zuo, Dong-Bo Ling, L.Sheng, D.Y.Xing, Optical properties for topological insulators with metamaterials, Phys. Lett. A 377, 2909 (2013).
  * (23) A. Martín-Ruiz, M. Cambiaso, and L. F. Urrutia, The magnetoelectric coupling in Electrodynamics.  Int. J. Mod. Phys. A 34, 1941002 (2019).
  * (24) A. Martín-Ruiz, M. Cambiaso, and L.F. Urrutia, Electro- and magnetostatics of topological insulators as modeled by planar, spherical, and cylindrical $\theta$ boundaries: Green’s function approach, Phys. Rev. D 93, 045022 (2016).
  * (25) A. Lakhtakia and T. G. Mackay, Classical electromagnetic model of surface states in topological insulators, J. Nanophoton. 10 (3), 033004 (2016).
  * (26) T. M. Melo, D. R. Viana, W. A. Moura-Melo, J. M. Fonseca, A. R. Pereira, Topological cutoff frequency in a slab waveguide: Penetration length in topological insulator walls, Phys, Lett. A 380, 973 (2016).
  * (27) Z.-X. Li, Yunshan Cao, Peng Yan, Topological insulators and semimetals in classical magnetic systems, Phys. Report 915, 1 (2021).
  * (28) R. Li, J. Wang, Xiao-Liang Qi and S.-C. Zhang, Dynamical axion field in topological magnetic insulators, Nature Phys. 6, 284 (2010).
  * (29) Y.Tokura, K.Yasuda, A. Tsukazaki, Magnetic topological insulators. Nat. Rev. Phys. 1, 126 (2019).
  * (30) A. Sekine and K. Nomura, Axion electrodynamics in topological materials, J. Appl. Phys. 129, 141101 (2021).
  * (31) D.M. Nenno, C.A.C. Garcia, J. Gooth, C. Felser, and P. Narang, Axion physics in condensed-matter systems. Nat. Rev. Phys. 2, 682–696 (2020).
  * (32) M. E. Tobar, B. T. McAllister, and M. Goryachev, Modified axion electrodynamics as impressed electromagnetic sources through oscillating background polarization and magnetization, Phys. Dark Universe 26, 100339 (2019).
  * (33) L. H. C. Borges, A. G. Dias, A. F. Ferrari, J. R. Nascimento, A. Yu. Petrov, Generation of Axion-Like Couplings via Quantum Corrections in a Lorentz Violating Background,  Phys. Rev. D 89, 045005 (2014).
  * (34) F. R. Prudêncio and M. G. Silveirinha, Optical isolation of circularly polarized light with a spontaneous magnetoelectric effect, Phys. Rev. A 93, 043846 (2016).
  * (35) R. Zhao, J. Zhou, Th. Koschny, E. N. Economou, and C. M. Soukoulis, Repulsive Casimir Force in Chiral Metamaterials, Phys. Rev. Lett. 103, 103602 (2009); M. G. Silveirinha and S. I. Maslovski, Comment on Repulsive Casimir Force in Chiral Metamaterials, Phys. Rev. Lett. 105, 189301 (2010).
  * (36) A. N. Darinskii, Surface plasmon polaritons in metal films on anisotropic and bianisotropic substrates, Phys. Rev. A 104, 023507 (2021).
  * (37) C. A. A. de Carvalho, The relativistic electron gas: a candidate for nature’s left-handed material,  Phys. Rev. D 93, 105005 (2016); E. Reyes-Gómez, L. E. Oliveira and C. A. A. de Carvalho, The electromagnetic response of a relativistic Fermi gas at finite temperatures: Applications to condensed-matter systems,  EPL 114, 17009 (2016).
  * (38) R.-Y. Zhang, Y.-W. Zhai, S.-R. Lin, Q. Zhao, W. Wen, M.-L. Ge, Time Circular Birefringence in Time-Dependent Magnetoelectric Media, Sci. Rep. 5, 13673 (2015); S.-R. Lin, R.-Y. Zhang, Y.-R. Ma, W. Jia, Q. Zhao, Electromagnetic wave propagation in time-dependent media with antisymmetric magnetoelectric coupling, Phys. Lett. A 380, 2582 (2016).
  * (39) K. Halterman, M. Alidoust and A. Zyuzin, Epsilon-near-zero response and tunable perfect absorption in Weyl semimetals, Phys. Rev. B 98, 085109 (2018).
  * (40) R. Zu, M. Gu, L. Min, C. Hu, N. Ni, Z. Mao, J. M. Rondinelli and V. Gopalan, Comprehensive anisotropic linear optical properties of Weyl semimetals, TaAs and NbAs, [arXiv:2011.09036] (2021).
  * (41) J. Krupka, Measurement of the complex permittivity, initial permeability, permeability tensor and ferromagnetic linewidth of gyromagnetic materials, Meas. Sci. Technol. 29, 092001 (2018).
  * (42) J. Krupka, A. Pacewicz, B. Salski, P. Kopyt, J. Bourhill, M. Goryachev and M. Tobar, Electrodynamic improvements to the theory of magnetostatic modes in ferrimagnetic spheres and their applications to saturation magnetization measurements, J. Magnetism and Magnetic Materials, 487, 165331 (2019).
  * (43) P. Hillion, Manifestly covariant formalism for electromagnetism in chiral media, Phys. Rev. E 47, 1365 (1993).
  * (44) Y. Itin, Dispersion relation for electromagnetic waves in anisotropic media, Phys. Lett. A 374, 1113 (2010).
  * (45) V. A. Markel and J. C. Schotland, Homogenization of Maxwell’s equations in periodic composites: Boundary effects and dispersion relations, Phys. Rev. E 85, 066603 (2012).
  * (46) N.J. Damaskos, A.L. Maffett and P.L.E. Uslenghi, Dispersion relation for general anisotropic media, IEEE Trans. Antennas Propagat. AP-30, 991 (1982).
  * (47) Y. Takahashi, R. Shimano, Y. Kaneko, H. Murakawa, and Y. Tokura, Magnetoelectric resonance with electromagnons in a perovskite helimagnet, Nature Phys. 8, 121 (2012); S. Iguchi, R. Masuda, S. Seki, Y. Tokura, Y. Takahashi, Enhanced gyrotropic birefringence and natural optical activity on electromagnon resonance in a helimagnet, Nature Communications 12, 6674 (2021).
  * (48) T. Kurumaji, Y. Takahashi, J. Fujioka, R. Masuda, H. Shishikura, S. Ishiwata, and Y. Tokura, Optical Magnetoelectric Resonance in a Polar Magnet $(Fe,Zn)_{2}2Mo_{3}O_{8}$ with Axion-Type Coupling, Phys. Rev. Lett. 119, 077206 (2017).
  * (49) A. Novitsky, A. S. Shalin, and A. V. Lavrinenko, Spherically symmetric inhomogeneous bianisotropic media: Wave propagation and light scattering,  Phys. Rev. A 95, 053818 (2017).
  * (50) V. A. De Lorenci and D. D. Pereira, Magnetoelectric birefringence as unique effect in isotropic media, Phys. Rev. E 82, 036605 (2010).
  * (51) V. A. De Lorenci and G. P. Goulart, Magnetoelectric birefringence revisited, Phys. Rev. D 78, 045015 (2008) .
  * (52) V. A. De Lorenci, Aspects of wave propagation in nonlinear medium: birefringence and the second-order magnetoelectric coefficients, Phys. Rev. A 105, 023530 (2022).
  * (53) D.E. Kharzeev, The chiral magnetic effect and anomaly-induced transport, Prog. Part. Nucl. Phys. 75, 133 (2014).
  * (54) K. Fukushima, D.E. Kharzeev, and H.J. Warringa, Chiral magnetic effect, Phys. Rev. D 78, 074033 (2008).
  * (55) Z. Qiu, G. Cao and X.-G. Huang, Electrodynamics of chiral matter, Phys. Rev. D 95, 036002 (2017).
  * (56) P.D.S. Silva, M.M. Ferreira Jr., M. Schreck, and L.F. Urrutia, Magnetic-conductivity effects on electromagnetic propagation in dispersive matter, Phys. Rev. D 102, 076001 (2020).
  * (57) S. Kaushik, D.E. Kharzeev, and E.J. Philip, Transverse chiral magnetic photocurrent induced by linearly polarized light in symmetric Weyl semimetals, Phys. Rev. Research 2, 042011(R) (2020).
  * (58) M . E. Tobar, P. Wolf, A. Fowler, and J. G. Hartnett, New methods of testing Lorentz violation in electrodynamics,  Phys. Rev. D 71, 025004 (2005).
  * (59) Q.G. Bailey and V.A. Kostelecký, Lorentz-violating electrostatics and magnetostatics, Phys. Rev. D 70, 076006 (2004).
  * (60) P.D.S. Silva, L. Lisboa-Santos, M. M. Ferreira Jr., and M. Schreck, Effects of CPT-odd terms of dimensions three and five on electromagnetic propagation in continuous matter, Phys. Rev. D 104, 116023 (2021).
  * (61) L. Ohnoutek et. al., Strong interband Faraday rotation in 3D topological insulator $Bi_{2}Se{3}$, Sci. Rep.6, 19087 (2016).
  * (62) L. Wu, M. Salehi, N. Koirala, J. Moon, S. Oh, N. P. Armitage, Quantized Faraday and Kerr rotation and axion electrodynamics of a 3D topological insulator, Science 354, 6316 (2016).