sisl to extract RHO file issue

[el-abed]

Describe the issue Good evening to whom it may concern, I am using sisl generally and specifically the following command to plot the charge transfer in 3D

sgrid DeltaRho.grid.nc --geometry zeroPG.fdf Rho.cube

I would expect 2 colors for zero twist phospherene/ graphene. One color to represent excessive number of electrons on Phospherene and another color for deprive electrons on graphene. However i got the following figure of mixture of colors and I am not sure why.

I tried to send the DeltaRHO.nc and RHO.cube files but they are very very large. Any other way of sending files here. I tried zipping them still did not work. Thank you very much and looking forward to your reply

[zerothi]

Colors in VMD (It seems this is VMD) is determined using the coloring method. This has nothing to do with sisl :(

If you want two colors, you need to do something like this:

1. Add a "Representation", select "Drawing Method" as "Isosurface", select "Coloring Method" as "ColorID", then select an isovalue and a color value.
2. Repeat step 1 but select the opposite isovalue and another color.

[el-abed]

I am using VESTA not VMD. Is that a wrong choice to plot cube? Should I be using VMD because i made sure i had a color for each case but yet the isosurfaces are not on the atom. They are around it as shown:

Does that mean I have to choose another software ? which one is best for such studies? Thank you for sisl! Really enjoying it as a beginner!

[zerothi]

I have no clue about vesta... :(

You can use whatever to plot cube files, XCrySDen, VMD or VESTA. However, the different colors does not originate from sisl. It is a pure graphical thing from VESTA. About the colors you should ask VESTA.

As for the mis-alignment of the geometry and the grid. Could you please do (in your work directory):

ls -l
sgrid --version

Also, were you running an MD simulation?

[el-abed]

To answer your first question:

Second question I was doing a siesta run to find the charge transfer with the input file attached. fdf-file.txt

What do you think? I will ask VESTA and get back to you just for people to know! Thank you alot for the speedy replies

[zerothi]

The ls -l should be in the directory where you did your siesta calculation.

Anyways, I can see you did an MD calculation. This means that the fdf coordinates are not necessarily the same as the final structure.

Try this instead:

sgrid DeltaRho.grid.nc --geometry zeroPG.XV Rho.cube

[el-abed]

For ls -l :

[eh8122@gadi-login-03 V_0]$ ls -l total 8332480 -rw-r--r-- 1 eh8122 g46 30 May 5 09:34 0_NORMAL_EXIT -rw-r--r-- 1 eh8122 g46 249 May 5 09:31 BASIS_ENTHALPY -rw-r--r-- 1 eh8122 g46 258 May 5 09:31 BASIS_HARRIS_ENTHALPY -rw-r--r-- 1 eh8122 g46 88473940 May 5 09:31 BaderCharge.grid.nc -rw-r--r-- 1 eh8122 g46 298069 May 1 12:49 C.dat -rw-r--r-- 1 eh8122 g46 321638 May 4 13:47 C.ion -rw-r--r-- 1 eh8122 g46 53740 May 4 13:47 C.ion.nc -rw-r--r-- 1 eh8122 g46 298027 May 4 13:47 C.ion.xml -rw-r--r-- 1 eh8122 g46 143097 May 1 13:44 C.psf -rw-r--r-- 1 eh8122 g46 18385 May 5 09:34 CLOCK -rw-r--r-- 1 eh8122 g46 176947540 May 5 09:31 DeltaRho.grid.nc -rw-r--r-- 1 eh8122 g46 88473952 May 5 09:31 ElectrostaticPotential.grid.nc -rw-r--r-- 1 eh8122 g46 25678 May 5 08:07 FORCE_STRESS -rw-r--r-- 1 eh8122 g46 20319 May 4 13:47 INPUT_TMP.17575 -rw-r--r-- 1 eh8122 g46 20319 Apr 30 14:30 INPUT_TMP.47185 -rw-r--r-- 1 eh8122 g46 20 May 5 09:34 MESSAGES -rw-r--r-- 1 eh8122 g46 471 May 4 13:47 NON_TRIMMED_KP_LIST -rw-r--r-- 1 eh8122 g46 480 May 5 09:34 OCCS -rw-r--r-- 1 eh8122 g46 1299 May 1 11:08 OUTVARS.yml -rw-r--r-- 1 eh8122 g46 247681 May 1 12:50 P.dat -rw-r--r-- 1 eh8122 g46 321638 May 4 13:47 P.ion -rw-r--r-- 1 eh8122 g46 53740 May 4 13:47 P.ion.nc -rw-r--r-- 1 eh8122 g46 298027 May 4 13:47 P.ion.xml -rw-r--r-- 1 eh8122 g46 153443 May 1 13:44 P.psf -rw-r--r-- 1 eh8122 g46 2496 May 4 13:47 PARALLEL_DIST -rw-r--r-- 1 eh8122 g46 272295218 Jun 9 21:16 Rho.cube -rw-r--r-- 1 eh8122 g46 176947532 May 5 09:34 Rho.grid.nc -rw-r--r-- 1 eh8122 g46 413666 May 1 12:46 T.dat -rw-r--r-- 1 eh8122 g46 5403 May 5 09:34 TIMES -rw-r--r-- 1 eh8122 g46 88473940 May 5 09:31 TotalCharge.grid.nc -rw-r--r-- 1 eh8122 g46 176947544 May 5 09:31 TotalPotential.grid.nc -rw-r--r-- 1 eh8122 g46 4705 May 5 17:05 V_0.csv -rw-r--r-- 1 eh8122 g46 60524 May 5 09:31 fdf-17782.log -rw-r--r-- 1 eh8122 g46 60524 May 1 11:05 fdf-47370.log -rw-r--r-- 1 eh8122 g46 1299 May 5 09:34 fort.14 -rw-r--r-- 1 eh8122 g46 230 May 3 15:07 run_qsub -rw------- 1 eh8122 g46 14 May 1 11:08 run_qsub.e4899052 -rw------- 1 eh8122 g46 86 May 3 15:07 run_qsub.e4989678 -rw------- 1 eh8122 g46 14 May 5 09:34 run_qsub.e5024055 -rw------- 1 eh8122 g46 757 May 1 11:09 run_qsub.o4899052 -rw------- 1 eh8122 g46 754 May 3 15:08 run_qsub.o4989678 -rw------- 1 eh8122 g46 757 May 5 09:34 run_qsub.o5024055 -rw-r--r-- 1 eh8122 g46 10615 May 5 17:05 trial.csv -rw-r--r-- 1 eh8122 g46 28943 May 5 17:05 trial.txt -rw-r--r-- 1 eh8122 g46 51388 Apr 30 14:30 vdw_kernel.table -rw-r--r-- 1 eh8122 g46 1247040 May 5 07:50 zeroPG.ANI -rw-r--r-- 1 eh8122 g46 90112104 May 5 09:31 zeroPG.BADER -rw-r--r-- 1 eh8122 g46 262392 May 4 13:47 zeroPG.BONDS -rw-r--r-- 1 eh8122 g46 262392 May 5 08:07 zeroPG.BONDS_FINAL -rw-r--r-- 1 eh8122 g46 18440 May 5 08:07 zeroPG.CG -rw-r--r-- 1 eh8122 g46 94528796 May 5 08:06 zeroPG.DM -rw-r--r-- 1 eh8122 g46 183000 May 5 09:30 zeroPG.DOS -rw-r--r-- 1 eh8122 g46 180224104 May 5 09:31 zeroPG.DRHO -rw-r--r-- 1 eh8122 g46 1509707 May 5 08:08 zeroPG.EIG -rw-r--r-- 1 eh8122 g46 23187 May 5 08:07 zeroPG.FA -rw-r--r-- 1 eh8122 g46 132583564 May 5 08:07 zeroPG.HSX -rw-r--r-- 1 eh8122 g46 555 May 4 13:47 zeroPG.KP -rw-r--r-- 1 eh8122 g46 180224104 May 5 09:34 zeroPG.LDOS -rw-r--r-- 1 eh8122 g46 1460160 May 5 07:50 zeroPG.MD -rw-r--r-- 1 eh8122 g46 1478080 May 5 07:50 zeroPG.MD.nc -rw-r--r-- 1 eh8122 g46 5658 May 5 07:50 zeroPG.MDE -rw-r--r-- 1 eh8122 g46 1504880 May 5 07:50 zeroPG.MD_CAR -rw-r--r-- 1 eh8122 g46 3468694 May 4 13:47 zeroPG.ORB_INDX -rw-r--r-- 1 eh8122 g46 550111382 May 5 09:31 zeroPG.PDOS -rw-r--r-- 1 eh8122 g46 601085122 May 5 09:31 zeroPG.PDOS.xml -rw-r--r-- 1 eh8122 g46 180224104 May 5 09:31 zeroPG.RHO -rw-r--r-- 1 eh8122 g46 24507 May 5 07:50 zeroPG.STRUCT_NEXT_ITER -rw-r--r-- 1 eh8122 g46 24507 May 5 08:07 zeroPG.STRUCT_OUT -rw-r--r-- 1 eh8122 g46 90112104 May 5 09:31 zeroPG.TOCH -rw-r--r-- 1 eh8122 g46 90112104 May 5 09:31 zeroPG.VH -rw-r--r-- 1 eh8122 g46 180224104 May 5 09:31 zeroPG.VT -rw-r--r-- 1 eh8122 g46 46338 May 5 07:50 zeroPG.XV -rw-r--r-- 1 eh8122 g46 414 May 4 13:47 zeroPG.bib -rw-r--r-- 1 eh8122 g46 20827 May 1 13:44 zeroPG.fdf -rw-r--r-- 1 eh8122 g46 3126764256 May 5 08:08 zeroPG.fullBZ.WFSX -rw-r--r-- 1 eh8122 g46 1952965425 May 5 09:34 zeroPG.out -rw-r--r-- 1 eh8122 g46 5402 May 5 09:34 zeroPG.times -rw-r--r-- 1 eh8122 g46 25848 May 5 08:07 zeroPG.xtl -rw-r--r-- 1 eh8122 g46 15588 May 5 08:07 zeroPG.xyz

I used your latest command i got a major change

after making sure yellow is +0.03 e and blue
-0.03 e. It is getting better but for some reason both colors are mixed up. But that is better. May I ask if there are more tutorials for sisl other than: https://sisl.readthedocs.io/en/latest/scripts/sgrid.html

Otherwise NICK i really want to thank you :) I will ask VESTA and should get an answer soon on why the mixture of positive and negative charges on both graphene and phosphorene unless my calculations are wrong?

[zerothi]

It looks correct, on both graphene and phosphorene. The DeltaRho.grid.nc file contains the deflection charge. It basically gives the re-arrangement of charge due to the environment.
So what you see is how the electrons move due to the surrounding atoms. If there is nothing in the graphene layer it means that all atoms in the graphene sheet behaves as "lone C atoms" (i.e. no charge re-arrangement).

[zerothi]

Also, please update, 0.9.8 is out, and has some bug-fixes.

[el-abed]

Thank you very much Nick! Could you further explain the physics behind no charge rearrangement because I assumed because of the Voronoi charge analysis, i have got average of 0.03 e transfer from graphene to phosphorene. Yet sisl is giving a different explanation? Could you further the explanation please?

[zerothi]

I guess it isn't.

Your picture really depends a lot on your iso-surfaces, and the actual integration of the difference regions.

My guess is that if you do an integration of the plot around the graphene atoms, and around the P atoms, you would find the charge transfer. You have a VERY large unit-cell which means that the 0.3 electrons are distributed for the whole cell. So variations are very small. You should analyze much closer, (don't rely on your eyes ;) )

[el-abed]

Really appreciate your time! So doing just graphene on its own then phosphorene on its own should do the trick. So your argument about lone carbon is not the case it is just the factor that the charge average for either P or C is very small. Is that correct? Because I have done that and i got them! Here for example graphene and here for P:

As for VESTA, it seems they are convinced other wise and I quote: "isosurface is the surface around the atoms that means it can go between two atoms obviously, at which the given value of density is met" Again thank you and I understand if the issue is closed now :) Just wanted people in the future to know how to use it. PS Are there any further sisl tutorials?

[zerothi]

So doing just graphene on its own then phosphorene on its own should do the trick.

What trick?

So your argument about lone carbon is not the case it is just the factor that the charge average for either P or C is very small. Is that correct?

No, you misunderstood... If DeltaRho for pristine graphene gives you 0's at all grid points it means it is behaving as a lone C atom (but it doesn't, luckily ;)).

What you should analyze is how the actual charge distribution is happening in the individual layers. An iso-surface plot won't give you this information directly. Hence I claim that if you integrate the DeltaRho grid you would find exactly what the Voronoi gives.

Something like this will give you the average charge redistribution along the z axis.

import numpy as np
import sisl as si
grid = si.get_sile('DeltaRho.grid.nc').read_grid()
z_vals = grid.average(0).average(1).grid.ravel()
dz = np.linspace(0, grid.sc.length[2], grid.shape[2])
import matplotlib.pyplot as plt
plt.plot(dz, z_vals)
plt.show()

Because I have done that and i got them! Here for example graphene

If you integrate the DeltaRho grid you should get roughly 0. However, your iso-surfaces could have different surface areas as for graphene. You'll see exactly what you have, some blue blops and some very small red blops. The integrated value is ~ 0.

and here for P: ... As for VESTA, it seems they are convinced other wise and I quote:

Convinced "other wise"? I don't get what you mean?

"isosurface is the surface around the atoms that means it can go between two atoms obviously, at which the given value of density is met"

I completely agree, and I didn't say anything to the contrary :)

Again thank you and I understand if the issue is closed now :) Just wanted people in the future to know how to use it. PS Are there any further sisl tutorials?

There are more tutorials here, mainly targeted TS + TBtrans and tight-binding calculations. https://github.com/zerothi/ts-tbt-sisl-tutorial

I'll close. Thank you.

[el-abed]

Hello again Nick, Sorry to reopening the issue again but I have a follow up question: Since the best way to actually analyze how the actual charge distribution is happening in the individual layers i should be calculating the Delta RHO of the system minus each of the individual layers. My question is the following command correct:

sgrid DeltaRho.grid.nc --geometry zeroPG.XV --diff P.grid.nc --geometry zeroP.XV --diff G.grid.nc --geometry zeroG.XV diff.cube

where i intend to find the deltarho difference between the 2D phosphorene/Graphene unitcell and each sheet separately each after a siesta run having their own XV file. IS THAT CORRECT? THANK YOU AND LOOKING FORWARD TO YOUR REPLY!

[zerothi]

The sgrid command only takes a single geometry as input. This geometry will be the one that is associated with all structures. It has no influence of the grid values, but will be used when writing the cube file (i.e. the geometry in the resulting plot).

From what I can see what you do is:

1. Read in DeltaRho.grid.nc, we call this grid A
2. From A subtract the values in P.grid.nc.
3. From A subtract the values in G.grid.nc
4. Then finally write out the resulting values to diff.cube

I think you should do:

sgrid DeltaRho.grid.nc --geometry zeroPG.XV --diff P.grid.nc --diff G.grid.nc diff.cube

which results in the same.

Edit: Also, please do not use CAP letters, it is confusing and not necessary ;)

[el-abed]

Thank you a lot for the quick reply. The reason i am asking for multiply geometry is because each geometry is relaxed separately. I thought i should be using for each calculation of RHO its own XV file. But it does not seem to be the case, Why? Thank you again and sorry for opening multiple threads and capital letters. looking forward to your reply

[zerothi]

sisl only reads data in this case. The geometry is only used for writing out stuff. So it depends on which geometry you want to have in the final output file. Again, only one geometry is useful.