Birefringence

Crystal symmetries and Birefringence

Birefringence occurs when there is a varying dielectric constant due to an inherent crystal structure in a material. By default it is assumed that materials are amorphous and do not have any particular symmetry groups. But if one is to use a crystal within the program here are the following symmetry groups for permittivity. As a quick reminder, here is the relationship between permittivity and refractive index. Typically μ is close to 1 at optical frequencies.

$$n = \sqrt{\epsilon_r \mu_r}$$

Supported Crystal Structures:

Triclinic Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & b & c\\ b & d & e\\ c & e & f \end{bmatrix}$$

Monoclinic Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & b\\ 0 & c & 0\\ d & 0 & e \end{bmatrix}$$

Orthorhombic Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & 0\\ 0 & b & 0\\ 0 & 0 & c \end{bmatrix}$$

Tetragonal Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & 0\\ 0 & b & 0\\ 0 & 0 & b \end{bmatrix}$$

Trigonal Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & 0\\ 0 & a & 0\\ 0 & 0 & b \end{bmatrix}$$

Hexagonal Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & 0\\ 0 & a & 0\\ 0 & 0 & b \end{bmatrix}$$

Cubic Crystals have the following form:

$$\epsilon = \begin{bmatrix} a & 0 & 0\\ 0 & a & 0\\ 0 & 0 & a \end{bmatrix}$$

Crystal Axis rotations

Additionally any crystal axis can be rotated relative to the part's mechanical axis, which is the default crystal axis. We can do this by applying a rotation matrix to our permittivity. Where we define our rotation matrix similarly to how we define it for geometric rotations. (ZYX rotation order)

$$\epsilon' = \begin{bmatrix} a_{xx} & a_{xy} & a_{xz}\\ a_{yx} & a_{yy} & a_{yz}\\ a_{zx} & a_{zy} & a_{zz} \end{bmatrix} \epsilon \\ \quad \\ \begin{bmatrix} a_{xx} & a_{xy} & a_{xz}\\ a_{yx} & a_{yy} & a_{yz}\\ a_{zx} & a_{zy} & a_{zz} \end{bmatrix} = \begin{bmatrix} \cos\theta\cos\phi & \cos\theta\sin\phi\sin\psi-\sin\theta\sin\phi & \cos\theta\sin\phi\cos\psi+\sin\theta\sin\psi\\ \sin\theta\cos\phi& \sin\theta\sin\phi\sin\psi+\cos\theta\cos\psi& \sin\theta\sin\phi\cos\psi-\cos\theta\sin\psi\\ -\sin\phi & \cos\phi\sin\psi & \cos\phi\cos\psi \end{bmatrix}$$

Isotropic Media

For Isotropic media the usual approach is taken using P matrices and calculating the roots of the determinant for the refraction and reflection equation as we can see below. Then from there the Eigenvalues and vectors are found for calculating the change in the Polarization. Where the negative sign on the sqrt is the reflected vector and the positive is for the refracted.

$$0=\left\lfloor\epsilon+ n_2^2 \begin{bmatrix} 0& -\hat{k}_z' & \hat{k}_y'\\ \hat{k}_z' &0& -\hat{k}_x'\\ -\hat{k}_y' & \hat{k}_x' & 0 \end{bmatrix}^2 \right\rfloor,\quad \vec{k}' = n_1 \hat{k} \left[\pm\sqrt{n_2^2-n_1^2+n_1^2\hat{k}\cdot \hat{N} } - n_1 \hat{k}\cdot\hat{N} \right] \hat{N}$$

Non-Isotropic Media

Calculating the response in a non-isotropic medium can be tricky, but if the medium is smoothly varying (C1 continuous), then we can calculate ray trajectories using similar technique as what is done for the GRIN solver. But instead we have to be a little bit more careful in calculating this trajectory. Instead of 1 trajectory we will have 2. Because we do not have any discrete boundaries, we should be able to calculate how these two rays refract the standard GRIN tracing way, given that we keep the solutions straight. Essentially we find the inverse of the determinant to calculate what n2 is and calculate its derivative and perform GRIN calculations the usual way.