# Choosing simulation time for Fabry Perot Cavity

Hello,
I am currently simulating a Gold Fabry perot cavity with the following geometry:
Mirror length (Y-axis) 150u , Mirror width (x-axis 5u) , Slit width in the middle of each mirror is 4u (in the y-axis) and 5u (in the x-axis) , and the separation between the 2 mirrors 50u.
I am using a gaussian source with a beam radius of 4.5u that is injected perpendicular to the plane of the mirror at the middle of the slit as indicated from the geometry in the attached photos. My objective is to measure the transmission at different wavelengths (1.85u-2.05u)
.
While choosing the FDTD solver region parameters I tried calculating the simulation time based on my theoretical knowledge using equations that describe the field decay inside the Fabry perot cavity which is mainly dependent on the Finesse and the length of the Cavity.
Td= FL/c*pi
where F is the Finesse , L is the separation between the 2 Mirrors and c is the speed of light.
I calculated the Finesse using the Gold reflectivity function and I calculated the reflictivty based on the refractive index of gold at a wavelength inside the range I am simulating.
Theoretically this Td is the time for the field to decay inside the cavity by a factor of e.
I tried interpolating the exponential decay function to calculate the simulation time when the field is at 10^-5 of itâ€™s initlal intensity which is the quantity used by Lumerical FDTD for the early shut off to operate. However, the simulation time I calculated is always faraway from the actual simulation time taken by the solver.
My conclusion is that there are several factors arenâ€™t taken into consideration in the equation however it affects the solverâ€™s algorithm. Forexample the mirror length (in the y-Axis) , the range of wavelength calculated and plotted against the transmission and the number of the frequency points.
I would like to ask if there is an approximate method to calculate the simulation time theoretically taking into consideration all the parameters affecting it before running the simulation then adding it in the solver properties?

Hello @amr.gamaleldin,

Thanks for posting and welcome to the community. I do not know of a way to determine the time for the fields to decays in a high finesse cavity a priori. A better way of estimating total sim time required would be running a shorter simulation and following the outline in this app gallery page on simple High Q cavity.

If you are simply interested in the transmission I would recommend a frequency domain solver which would be much more efficient. This could be accomplished by be defining the index and thickness as vectors and using the stackrt command or if you have the most recent version of FDTD (2019b r2) you could use the stack solver user-interface.

Regards,
Taylor

1 Like

Hello @trobertson
Thank you very much for your reply. Is the stackrt useful only for high finesse cavities or can I use it for low finesse cavities too such as silicon cavities? How can I use stackrt for gold slotted mirror in my original post while this method neglects the second dimension? I tried using stackrt for a Fabry perot cavity made of Silicon (for Near infrared range (Lamda : 1.4-1.6 Micro). I added the thickness and the refractive index of each layer including air but I obtained different result for the transmission spectrum than the FDTD graph that I obtained for the same simulation as follows.
I am not sure which one of them is more ideal and which one is more realistic . It seems for me that the stackrt method is less accurate as it neglect the 2nd dimension so this explains why the transmission reaches 100%. however,why is there a shift in the peak wavelengths between the two graphs. I also tried to use the fdtd method but with planewave source and periodic boundary conditions in the y-direction. I have taken an FDTD region with only 0.5micro meters length. and I obtained a transmission graph that is similar to the stackrt graph however, the shift in wavelengths is still there and the graph seems to be mirrored.
Note: for all the attached graph pictures , the y axis is the transmission percent and the x-axis is the wavelengths.
Here is the code I used for stackrt method

# start of code ____________________

f = linspace(c/1400e-9, c/1600e-9,20000); # frequency vector

d = [0; 2000e-9; 50000e-9; 2000e-9; 0]; # 5 layers (including air on top, bottom

# get RT using non-dispersive index data, and theta=0

#RT1 = stackrt(n1,d,f);

#field1 = stackfield(n1,d,f);

# alternate refractive index vector, for dispersive materials

nf = length(f);

nd = length(d);

n2 = matrix(nd,nf);

n2(1,1:nf) = 1; # air

n2(2,1:nf) = getfdtdindex(â€śSi (Silicon) - Palikâ€ť,f,min(f),max(f));

n2(3,1:nf) = 1; # air

n2(4,1:nf) = getfdtdindex(â€śSi (Silicon) - Palikâ€ť,f,min(f),max(f));

n2(5,1:nf) = 1; # air

# get RT using dispersive data (n2), and theta from 0 to 45 deg

RT2 = stackrt(n2,d,f);

#field2 = stackfield(n2,d,f);

# open data visualizer

#visualize(RT1);

visualize(RT2);

#visualize(field1);

#visualize(field2);

# make simple plots

#plot(RT1.lambda*1e6,RT1.Rp,RT1.Rs,RT1.Tp,RT1.Ts,â€śwavelength (um)â€ť,â€śPowerâ€ť,â€śnon-disperisive, theta=0â€ť);

#legend(â€śRpâ€ť,â€śRsâ€ť,â€śTpâ€ť,â€śTsâ€ť);

image(RT2.lambda*1e6,RT2.theta,RT2.Rp,â€śwavelength (um)â€ť,â€śtheta (deg)â€ť,â€śRp, dispersive exampleâ€ť);

# end of code _________________________

Stackrt transmission

Layout for fdtd

transmission for fdtd

Transmission for fdtd using planewavesource and very small FDTD Region

Hey @amr.gamaleldin,

Let me see if I can answer your questions above and give you some guidance on how to proceed with you application now that I understand it a bit better.

I believe that you are using two different fits for Silicon between stackRT and FDTD which is causing the shift in transmission. You can refer to getfdtdindex, for more information on this, and check the material explorer next to run button to see the index. The â€śdispersiveâ€ť index that you are returning is not in fact dispersion, see your stack below with frequency along x.

It looks like the sim time in the first simulation was not long enough for the fields to decay properly, causing the oscillations and incomplete transmission. This can be changed by editing the FDTD object. There are a few ways the simulation will terminate if the max sim time is reached or if the auto shutoff level is reached. There is an interesting discussion on the matter here. Look at the log files for information on your simulation.

Since you have a 2D structure you may have to perform 3D FDTD simulation to get the correct results. How long do these simulations usually take to run, and do they terminate due to auto-shutoff or max time reached? Maybe run a simplified version with unpatterned gold in 2D to get an an idea of the simulation time. You can set the sim time very high to make sure auto shutoff is reached, and reduce the simulation volume as much as possible.

Regards,
Taylor

1 Like

@trobertson
your reply is very useful, now everything is much more clear , thank you very much!
Update: I understand now the reason behind the difference between the resulting transmission from the two graphs and I have a better idea about choosing simulation time. However, I donâ€™t understand in my case specifically why will a 3D FDTD Simulation improves my results?
and did you mean by â€śthe first simulationâ€ť the gold slotted cavity or the Silicon Fabry perot cavity?

Hey @amr.gamaleldin,

Yes I was referring to the gold slotted cavity when I said 3D FDTD may be required. Is the Au cavity patterned or by slit width do you mean the separation of the metallic layers? If it is patterned then 3D FDTD will likely be required. If if it is unpatterned than 2D FDTD or, STACK frequency domain solver is appropriate as in the dielectric example you posted.

Does that make sense.

Regards,
Taylor

Hey @trobertson
yes there is a separation (etch material) between the metallic layers.