  # Simulation of graphene antenna array

Hi,
I am trying to reproduce the results presented in the following paper:
https://pubs.acs.org/doi/pdf/10.1021/nl500602n

Which have been explained in the following post also:

But, I am unable to reproduce the results. I followed the instructions explained in the post and designed my structure using script, which is attached with the post, please check it and tell me, what is wrong with it? Apart from that, I have some queries regarding that post. such as;

1. How to check the auto shut off level, whether it is decreasing till 1e-5 or not?
2. It is suggested that PML should be at half of the maximum wavelength from the structure. How much distance from structure should be used for symmetric/asymmetric/periodic/Bloch boundaries?
3. At what height, plane wave should be kept from the structure? Its span should be larger that FDTD, is there any exact span size with respect to FDTD region?
4. In Fig. 1 of the paper following parameters are mentioned:
(i) Area of the square graphene sheet= 900 μm2
(ii) Antenna length= 0.95 μm ,
(iii) Antenna width= 0.24 μm
(iv) Antenna thickness= 40 nm thick (Pd thickness is 10 nm and Au thickness is 30 nm) (v)gap size = 60 nm.
The whole antenna array is composed of N = 24 rows (lateral
period Py = 1.2 μm), each comprising M = 30 end-to-end coupled antennas.
In the re implemented file Graphene10_2dcond.fsp, you didn’t use 0.95 μm length of antenna, infact x pan is taken 0.57 um, how did you select that value?
Authors have given several details, such as graphene area, array size, period etc., but in periodic boundary condition, we are not specifying these details. Results should be affected by period size and no. of antenna. Please clear this point, how to take care of period size and no. of antennas in x and y direction?.
updated capasso.lsf (6.3 KB)

Thanks and Regards
Prateeksha

Regards
Prateeksha

Hello @prateeksha.s1,

Thank you for the questions, and sorry for the late reply.

I am unable to reproduce the results. I followed the instructions explained in the post and designed my structure using script, which is attached with the post, please check it and tell me, what is wrong with it

It can be difficult to determine why your simulation results are different than a paper’s. It looks like this paper did not specify all of the details of their simulation, so you cannot expect your results to be exactly the same.

I would suggest you try using steep angle type PML boundaries for this type of device. You should also increase the number of PML layers. You could also try performing convergence testing to determine what mesh size you should be using.

1. How to check the auto shut off level, whether it is decreasing till 1e-5 or not?

You can see the auto shutoff level at various points during the simulation in the log file that is automatically generated and saved in the same directory as your simulation file when your simulation is run.

1. It is suggested that PML should be at half of the maximum wavelength from the structure. How much distance from structure should be used for symmetric/asymmetric/periodic/Bloch boundaries?

The placement of these BCs is determined by the symmetries or periodicities of the structure. The BCs should be placed at the planes of symmetry (for symmetric/antisymmetric BCs) or at the edges of the unit cell (halfway between two antennas, in your case) for periodic structures.

1. At what height, plane wave should be kept from the structure? Its span should be larger that FDTD, is there any exact span size with respect to FDTD region?

To be safe, the plane wave source should be placed about a half wavelength above the structure. If you are simulating an infinite plane wave in a periodic structure, as long as the plane wave source spans the FDTD region it should not matter how large it extends beyond the FDTD region.

1. In Fig. 1 of the paper following parameters are mentioned:
(i) Area of the square graphene sheet= 900 μm2
(ii) Antenna length= 0.95 μm ,
(iii) Antenna width= 0.24 μm
(iv) Antenna thickness= 40 nm thick (Pd thickness is 10 nm and Au thickness is 30 nm) (v)gap size = 60 nm.
The whole antenna array is composed of N = 24 rows (lateral
period Py = 1.2 μm), each comprising M = 30 end-to-end coupled antennas.
In the re implemented file Graphene10_2dcond.fsp, you didn’t use 0.95 μm length of antenna, infact x pan is taken 0.57 um, how did you select that value?

This topic is from several years ago, so I can’t know for sure why certain parameters were chosen. I would guess that this was chosen because the length of the antenna did not have a significant effect on the results, but I’m not sure.

Authors have given several details, such as graphene area, array size, period etc., but in periodic boundary condition, we are not specifying these details. Results should be affected by period size and no. of antenna. Please clear this point, how to take care of period size and no. of antennas in x and y direction?.

By using periodic BCs, we are assuming that there are an infinite number of periods in the x and y directions. If you want to include the effects of the number of periods, you would have to simulate the entire device with all of the antennas. This would result in a very long simulation, however.

I hope this helps. Let me know if you have any questions.

I have few more queries regarding absorption calculation in graphene sheet, which is explained in the same post Simulation of graphene photodetector

I calculated graphene absorption in free standing graphene sheet for Ef= 0 eV and I am able to reproduce the results obtained in supplementry material of paper https://pubs.acs.org/doi/pdf/10.1021/nl500602n (Figure S3), the 2.2% absorption for broad range of wavelength. In that post three methods are explained for graphene absorption:

1. Conductivity method
2. Difference method
3. Analytical method based on transfer matrix method

I am able to reproduce the results using analytical method. Conductivity method follows the same trend as analytical method, but it shows lower values compared to the analytical model results. Moreover, the absorption graph by conductivity method shows so many fluctuations. I tried with coarse meshing also, but it didn’t work.

How to improve the absorption results provided by conductivity method, since analytical method can not be applied in case of antenna structures?

Difference method is not providing any results. Please check my code.

I am a bit skeptical about analytical method. It works on TMM method, which needs thickness of individual layer of stack. Absorption should increase with the increasing thickness of graphene, however in my calculations, I am getting the same absorption values for different thicknesses( at constant fermi level and scattering rate). Please clear this point.

Is there any other method for absorption calculation in 2D graphene layer?

I am attaching my .lsf and .fsp files along with the supplementary file of the above mentioned paper.

Thanks and Regards
Prateeksha

How to improve the absorption results provided by conductivity method, since analytical method can not be applied in case of antenna structures?

I noticed a couple of issues with your simulation file. The first is that you are using a plane wave with PML boundary conditions. This will lead to some issues, as described on this page. Because you are modelling an infinite graphene plane, you should use periodic boundary conditions instead. I would also recommend you try a 2D FDTD simulation rather than 3D, which will be faster and easier to set up.

I would recommend you base your simulation on the 2Dtest.fsp simulation file from this post on the topic you posted:

To simulate the graphene layer by itself you can just remove the other layers from this simulation.

Difference method is not providing any results. Please check my code.

In this line:

Pgraphene_diff = abs(Tbelow.T-Tabove.T)/sourcepower(fi);

You don’t need to normalize the monitor transmission results by the source power, they are already normalized. This is most likely the cause of the issue.

I am a bit skeptical about analytical method. It works on TMM method, which needs thickness of individual layer of stack. Absorption should increase with the increasing thickness of graphene, however in my calculations, I am getting the same absorption values for different thicknesses( at constant fermi level and scattering rate). Please clear this point.

This is because the graphene conductivity is divided by the graphene thickness in this line:
graphene_eps = 1 + 1i*sigmax/(eps0*2*pi*f*graphene_t);

I believe that this is done to convert the surface conductivity sigmax to a bulk conductivity.

Is there any other method for absorption calculation in 2D graphene layer?

Not that I know of, but I am not an expert in this field. The normal methods for calculating absorption in FDTD will not work because the graphene is a 2D layer.

I hope this helps. Let me know if you have any questions.

Thank you very much for your valuable suggestions!!

I simulated graphene using periodic boundary condition and I got the same results for “analytical” and “conductivity” methods, but still “difference method” was giving ‘0’ constant line. But when I used 1 nm meshing in z-direction, then I could get the same graph but with some ripples. So my questions regarding that are:

1. Is 2D-graphene very sensitive for fine meshing? Even though, we do not define any thickness in 2D graphene, then how did it understand 1 nm meshing? This finer meshing has been used in [2Dtest.fsp] also.

2. I got almost same flow in case of difference method, but small ripples are still there (even I used 24 PML layers and 1e-11 auto shut off level), If I reduced these no., then ripples are more evident. How to improve this graph?

3. I am interested in calculating graphene absorption at 1550 nm wavelength, hence I am using 2D model of graphene. But, I am totally confused seeing the post regarding graphene’s material properties (Scattering rate in surface conductivity of graphene) It is mentioned on (https://support.lumerical.com/hc/en-us/articles/360042244874-Graphene-surface-conductivity-material-model) that relaxation time and scattering rate are related as ;ζ = 1/2Г for volumetric model. I assume, it must be same in case of 2D model also. But what is the input value, which lumerical understand for scattering rate? For example, if tau= 1Ps, then from above relation, Г = 0.5×10^12 S^-1, lumerical should take h_bar*Г= 6.58×10^-16×0.5×10^12=0.000329 eV.
Is it correct or it will consider something else? Please rectify me, if I am doing wrong calculations.

4. In the post (Scattering rate in surface conductivity of graphene), they have done same calculations with Drude approximations. I could not match those values. However, this approximation is valid for mid IR range, that implies that chemical potential and scattering rate won’t be interdependent at high frequency region, is it true? Otherwise, if chemical potential is varied, then scattering rate also needs to varied.

5. How did you calculate fermi level and chemical potential in (https://support.lumerical.com/hc/en-us/articles/360042243634-Electro-optical-modulator-based-on-a-graphene-coated-waveguide), since they are also using 1.53 um wavelength, which is similar to my operating region.

Regarding, (Simulation of graphene photodetector), I am able to reproduce the result (Figure S3) using scattering rate = 0.001645 eV and fermi level Ef= 0.4 eV, however scattering rate is different according to my calculations.(https://pubs.acs.org/doi/suppl/10.1021/nl500602n/suppl_file/nl500602n_si_001.pdf.
I have following two questions regarding that post:

1. At first, I used 8 PML layers and auto shut off level as 1e-5, then I was getting some ripples. For improving the results, I increased PML layers to 24 and auto shut off level as 1e-11, then also I got similar flow and ripples, but absorption level is highly improved. I have attached both results and .lsf files for your reference. Am I using correct simulation set up? How to recognize, that which simulation result is correct?

2. They have given mobility µ=1,000 cm2/Vs. By assuming Ef = 0.4 eV, if I calculate using formula (µ=(eζv_f^2 A)/E_f ) (the formula is mention in eq (27) https://arxiv.org/ftp/arxiv/papers/1712/1712.08965.pdf), since it is valid for mid IR then I get scattering rate as 0.008225 eV. Is it correct calculation, if not then please tell me, where am I making mistake? I have attached the snapshot of calculation for your reference.  1.lsf (6.4 KB) capasso_analytical_absorption.lsf (2.3 KB)

Thanks and Regards
Prateeksha

Hi,

Thanks and Regards
Prateeksha

Sorry for the slow response @prateeksha.s1 .

Yes, the graphene layer is not affected by the mesh size but the fields around the graphene need to be properly resolved, so a fine mesh is required. As with most meshing issues, this requires convergence testing to find the right settings.

You are correct that ripples in frequency domain results often mean the simulation is too short. You will have to increase the simulation time as well in addition to reducing the autoshutoff level: Generally increasing the PML layers and decreasing the autoshutoff level will improve results, so I would guess that result is more accurate. As I mentioned above, you will probably have to increase simulation time as well as decreasing the autoshutoff levels.

Regarding your other questions, unfortunately I am not an expert on graphene so I can’t be very helpful with these types of theory-related questions. I can help you with questions related to setting up your simulation, but I would recommend you consult the literature referenced on our graphene pages for details on the models used in the graphene materials.

Let me know if you have any more questions.

Dear Sir,
Regarding the material properties of graphene, I have checked the graphene page and related reference also, but I didn’t find any relation between relaxation time and scattering rate which is clearly mentioned in case of volumetric model ;ζ = 1/2Г. This is important to know that what value of scattering rate does Lumerical understand otherwise it might create some ambiguity during fabrication.

Thanks and Regards
Prateeksha

Hello @prateeksha.s1,

Sorry for the slow reply. I agree that the definitions on the graphene material page and in the references are a bit confusing.

There is a definition of \tau in reference  from this page (which I will call \tau_{ref}). By comparing eqn. 3 in the reference with the second to last equation in this page, we can see that the Lumerical definition of \tau is different than the reference definition of \tau_{ref} by a factor of \hbar:

\tau_{ref} = \tau/\hbar

From this, and the relation \Gamma = 1/(2\tau), we can see that

\Gamma = \frac{1}{2 \tau_{ref} \hbar}

This is the definition that Lumerical is using. Note that the “scattering rate” used in the material creation is in eV (equal to Γ/ℏ, as discussed on this KX post).

You should compare the definitions used by the source of the material parameters that you are using in your simulations with Lumerical’s definitions to make sure they are correct. In particular, make sure that the units are correct.

I hope this helps. Let me know if you have any questions.