jonschlipf / RigorousCoupledWaveAnalysis.jl

Rigorous Coupled-Wave Analysis (RCWA) for nanophotonics simulations
GNU General Public License v3.0
33 stars 7 forks source link

axially anisotropic layer #2

Closed hadishamkhi closed 1 year ago

hadishamkhi commented 1 year ago

hi, would you please provide an example on setting anisotropic materials parameters for certain layer such as the liquid crystal?

jonschlipf commented 1 year ago

Hi, sorry for replying late. I would like to do so, but I never really worked with such materials, and lack the creativity to come up with an example. Do you have a favorite material (ideally non-dispersive) of which you could give me a permittivity matrix? Also recommendation for a layer stack would be nice. Maybe you have a reference publication. I can then implement an example, test and put it here.

hadishamkhi commented 1 year ago

Thanks for your response. I would recommend the structure presented in arxiv for testing. It is my group project. The liquid crystal fits the criteria, it is largely nondispersive and uniaxially anisotropic. I think you discussed this type of anisotropy in your paper supplementary. I would be happy to assist in implementing the Ph structure therein and testing it. Note, you can consider the liquid crystal to be isotropic near the interfaces due to the anchoring effect. The Ph also include a metasurface layer. The catch here is leverage the control over the Fabry-Perot resonance position of liquid crystal to tune its coupling strength with the metasurface Mie resonances.

jonschlipf commented 1 year ago

I wrote this code here for a plain LC layer from the paper. Unfortunately, discrimination between x and y polarized wave powers reflected and transmitted by the device is not fully implemented, so this one got lengthy and uses some RCWA internals.

using RigorousCoupledWaveAnalysis
using LinearAlgebra
no=1.52 #ordinary index
ne=1.81 #extraordinary index
LC=ConstantPermA([no^2,0,0,ne^2,no^2]) #sparse permittivity tensor (exx,exy,eyx,eyy,ezz)
#for other LC orientation, you would have to apply a rotation matrix to the first four elements
air=ConstantPerm(1.0)
wls=500:5:1500 #wavelength array
N=2
pitch=500
l1=AnisotropicLayer(1000,LC) #only 1um LC layer in air
mdl=RCWAModel([l1],air,air)
Rx=zeros(length(wls))
Ry=zeros(length(wls))
Tx=zeros(length(wls))
Ty=zeros(length(wls))
rx=zeros(length(wls))*1im
ry=zeros(length(wls))*1im
tx=zeros(length(wls))*1im
ty=zeros(length(wls))*1im
for i=1:length(wls) 
    grd=rcwagrid(N,N,500,500,1e-5,45,wls[i],air) #incoming wave is at 45 degrees polarization
    ste,stm=rcwasource(grd)
    kzin=grd.k0[3]#z component of incoming wave momentum
    ems=RigorousCoupledWaveAnalysis.eigenmodes(grd,wls[i],mdl.layers)#layer eigenmodes
    ref=RigorousCoupledWaveAnalysis.halfspace(grd.Kx,grd.Ky,mdl.εsup,wls[i])#superstrate eigenmodes
    tra=RigorousCoupledWaveAnalysis.halfspace(grd.Kx,grd.Ky,mdl.εsub,wls[i])#superstrate eigenmodes
    ro,to,r,t=etm_propagate(ref,tra,ems,stm,grd,false)#perform etm computation                                                        
    #need to create masks to select the 0-th order x and y
    maskx=zeros(length(ro))
    maskx[N*(2N+1)+N+1]=1     
    masky=zeros(length(ro))       
    masky[(2N+1)^2+N*(2N+1)+N+1]=1
    #compute powers
Tx[i]=-RigorousCoupledWaveAnalysis.a2p(to.*maskx,to*0,tra.V,I,kzin) #power in 0-th order x transmission
    Ty[i]=-RigorousCoupledWaveAnalysis.a2p(to.*masky,to*0,tra.V,I,kzin) #power in 0-th order y transmission
    Rx[i]=RigorousCoupledWaveAnalysis.a2p(0*ro,ro.*maskx,ref.V,I,kzin) #power in 0-th order x reflection
    Ry[i]=RigorousCoupledWaveAnalysis.a2p(0*ro,ro.*masky,ref.V,I,kzin) #power in 0-th order y reflection                    
    #for insight into these check src/ETM/ETM.jl
    #compute amplitudes (only valid in lossless media)
    ex,ey=RigorousCoupledWaveAnalysis.a2e2d(to,I) #identity matrix of correct size I             
    println(length(ex))
    tx[i]=ex[N*(2N+1)+N+1]
    ty[i]=ey[N*(2N+1)+N+1]
    ex,ey=RigorousCoupledWaveAnalysis.a2e2d(ro,I) #identity matrix of correct size
    rx[i]=ex[N*(2N+1)+N+1]
    ry[i]=ey[N*(2N+1)+N+1]
    print(kzin)
end
#now we can also compute polarization angle (the real part of this is the angle, the imaginary part is related to ellipticity:
phi_R=atand.(ry./rx).-45 #rotation in reflected wave
phi_T=atand.(ty./tx).-45 #rotation in transmitted wave

This is the polarization angle: image

These are the ref and tra powers in the respective polarization directions (Energy conservation Rx+Ry+Tx+Ty=1: image

I hope you can take it from here, and implement your designs. If you get something that resembles one of the major plots from the paper, I would like to include it as an example here.