Thank you for your response! I have been working on this in the meantime, considering both your suggestions and re-examining the index monitor results (I'll explain below).

First, I can see the use of encoding a custom material by it's refractive index, but in this situation it is possible and more efficient to code it directly.

When specifying custom material by its refractive index, the refractive index matrix must first be diagonalized so that at each point (x,y,z) only [n_{xx}, n_{yy}, n_{zz}] are specified (all off-diagonal elements are zero). My refractive index matrix is already diagonalized as seen in cylindrical coordinates, as it is:

n = [ n_{e}, 0, 0; 0, n_{o}, 0; 0, 0, n_{o} ] where n_{e} = 1.8 and n_{o} = 1.4

I have explored this method, and it gives the same results as specifying the anisotropy directly. When I state the refractive index at each point (x,y,z) in an 801x801x2 grid, the index monitor shows the same material properties as if I specify the refractive indices by material for the given structures -- though the latter is more exact and doesn't cause any points of intermeadiate refractive index along the borders of the structures, as in numerical calculation. (Note here that the environment is set to 1.34, the inner cylinder is set to 1.38, the outer layer radial direction is set to n_{e} = 1.49 and the outer layer phi and z directions is set to n_{0} = 1.48.)

Next, I appreciate your mentioning symmetry boundary conditions. I am concerned that this doesn't solve underlying asymmetry for the full-scale simulation. Symmetric boundary conditions are used for making the simulations more efficient but doesn't change any underlying behavior (when applied correctly). If I, for instance, only simulate the top-right quarter of the cross-section of the structure, and correctly apply boundary conditions, it should be the same as running the simulation with the whole cross-section. So the underlying asymmetry in these full simulations remains unexplained.

Due to the symmetry of the structure (including of the matrix transform, which is specifying radial and phi directions), anti-symmetric boundary conditions on the x and y axes should be appropriate. I am finding modes that have balanced intensity over the x and y axes, as shown in the image below (these show electric field intensities, the following behavior is also correspondingly very evident in the electric field vector plots of the modes). The position of the nodes seem to be affected by these boundary conditions -- again this is in a symmetric way, but I'm not confident that this describes the true behavior of the full cross-section (my point in the previous paragraph). These have resulted in rather "squarish" or "boxy" modes as opposed to smooth/round modes. My expectation (based off my prior simulations without the radial birefringence, as seen in the modes in the second picture in my original question) would be that the nodes and intensities are distributed more evenly, even to the point of forming a ring in the outer layer with uniform intensity.

You mention that mode degeneracy may be playing a role in the asymmetry. Can you expand on this a bit more? I have encountered the concept of mode degeneracy as being light of the same frequency comprising differing patterns -- usually related by some linear transformation (rotations and reflections)? -- but I don't see how this would affect the symmetry at all.

Lastly, I want to bring up another avneue of my investigation into this problem, and inquire if the software is performing calculations how I am expecting.

When I use the index monitor to view [n_{xx}, n_{yy}, n_{zz}] of the entire structure, I am finding constant values as shown in the picture below -- this is the index monitor image that I find with the software currently. On Lumerical's website (the link https://optics.ansys.com/hc/en-us/articles/360034915193-Tips-and-background-information-when-using-grid-attributes under section "Index monitor currently record only diagonal index") it mentions that the index monitor can only show these diagonal elements of the refractive index matrix, so encoding the anisotropy must be done carefully -- I would presume because the non-diagonal elements may not be zero. So I cannot view an element such as n_{xy} in the index monitor.

However, since I have the transform matrix between coordinate systems (the Jacobian J) and I know all the refractive indices, I can view the entire refractive index matrix in Cartesian coordinates by plotting this in MATLAB. I have included my calculation and the MATLAB plots for each matrix element below. Note here n_{rr} = 4, n_{phiphi} = 3, and n_{zz} = 2 were chosen to fully illustrate the behavior.

n_{Cartesian} = J^{-1} n_{cylindrical} J

n_{Cartesian} = [ cos(phi), sin(phi), 0; - sin(phi), cos(phi), 0; 0, 0, 1] [ n_{rr}, 0, 0; 0, n_{phiphi}, 0; 0, 0, n_{zz} ] [ cos(phi), - sin(phi), 0; sin(phi), cos(phi), 0; 0, 0, 1]

n_{Cartesian} = [ n_{rr} cos^2(phi) + n_{phiphi} sin^2(phi), ( n_{phiphi} - n_{rr}) cos(phi) sin(phi), 0; (n_{phiphi} - n_{rr}) cos(phi) sin(phi), n_{rr} sin^2(phi) + n_{phiphi} cos^2(phi), 0; 0, 0, n_{zz} ]

Even when restricting our view in my MATLAB plots to the diagonal elements n_{xx} and n_{yy}, there is clearly a gradient in the refractive index as observed in the xy plane.

If the matrix transformation is applied before the index monitor views the refractive index, then I should be seeing a gradient in the Lumerical index monitor as seen in these MATLAB plots. If the matrix transformation is applied after the index monitor views the refractive index, then I was correct in assuming that my current index monitor pictures are correct. Since the website notes that "index monitors do not record any information about the permittivity transformations" ( https://optics.ansys.com/hc/en-us/articles/360034915193-Tips-and-background-information-when-using-grid-attributes ) I assumed it was the latter, but I am reconsidering this and would appreciate clarification!

If my index monitor plots are correct, I am still unsure what the cause of the asymmetry is. If my index monitor plots are incorrect, there is clearly a problem encoding the anisotropy, and any guidance on how to fix this would be most appreciated!

Thank you again for your input, and I am open to further ideas from anyone!