Chengcheng-Xiao / RotSph

Rotation Matrix for Real Spherical Harmonics
MIT License
0 stars 2 forks source link

Plot PDOS after finding Euler's angles automatically #1

Closed CrankyBot2001 closed 4 months ago

CrankyBot2001 commented 4 months ago

Thank you for the great code, It's the what I was looking for. I did DOS calculation with LORBIT = 14 and then run find_angle.py for my structure ( contain different type transition metals and those octahedrals are tilted with different angle). So I am looking for how to plot PDOS for each $d$ orbitals after finding Euler's angle. I tried to modify your code rot_pband_plot.py for plotting $m_l$ projected PDOS, but unable to do that. If you provide the modified script for PDOS, it will be very helpful for me. Thank you in advance.

Chengcheng-Xiao commented 4 months ago

Hi @CrankyBot2001, the rot_pband_plot.py script is only for projected band structure plots. I've uploaded a script that can give you the pDOS plot (https://github.com/Chengcheng-Xiao/RotSph/tree/master/example/automatic_rotation).

Please feel free to post your results (DOS plots before and after rotation) here.

I would really appreciate it if this repo can be cited using the following temporary BibTeX entry in any publication that utilizes this code:

@software{rotsph,
   author = {Chengcheng Xiao},
   title = {rotsph: a code to calculate the rotation matrix of real spherical harmonics},
   url = {https://github.com/Chengcheng-Xiao/RotSph},
   year = {2024}
}

I do have a plan to merge this functionality to pymatgen in the future so this BibTeX entry will be updated in due time.

CrankyBot2001 commented 4 months ago

Thank you @Chengcheng-Xiao , I will let you know after calculating PDOS. I will definitely cite your code in my research. I want to confirm few things, that is the Fermi energy shifted to 0eV in this PDOS plot? As my system contains 66 atoms and so many octahedrals , I hope this code able to project $d$ orbitals along z axis for all transition metal. I want to select more than one atom and select all $d$ orbitals to plot $t_{2g}$-$e_g$ resolved DOS and also want to plot spin-polarized DOS. Again thank you for your great code.

Chengcheng-Xiao commented 4 months ago

I hope this code able to project orbitals along z axis for all transition metal.

By z axis I assume you mean in the rotated frame?

I want to select more than one atom and select all orbitals to plot t2g-eg resolved DOS and also want to plot spin-polarized DOS.

rot_pdos_plot.py output the results in pdos.dat, For your purpose, you can run it multiple times with different input parameters and manually add the results together and plot the sum yourself.

Spin polarized calculation require a small re-write of the script, which can be downloaded here.

The input should be:

0      |  6   |  Euler alpha | Euler beta | Euler gamma | dos E min | dos E max | dos points | smearing Sigma| spin |
Atom 0 |  dz2 |  -135        | 0          | -54.75      | -2        | -2        | 100        | 0.1           | 1 (up) or -1 (down) |
CrankyBot2001 commented 4 months ago

Yes I meant that. Yes I will sum those orbitals contribution to get t2g-2g pdos, and what about the Fermi energy, is it shifted to zero in your code?

I think it's not, so I subtracted the Fermi energy value from DOSCAR to shift at 0 eV. Thank you.

CrankyBot2001 commented 4 months ago

Look @Chengcheng-Xiao, after performing auto rotation of my 1st transition metal atom and plot pdos(E-F=0eV) for that particular atom, still I am getting plot like this where I couldn't distinguish which orbital is contributing where explicitly, getting mixed things at a particular energy point, high_quality_plot(TiCrO2)

Chengcheng-Xiao commented 4 months ago

To me this plot actually looks pretty good: $eg$ orbitals are located at ~3 eV and $t_2 g$ orbitals are at ~2 eV.

I guess you are confused by $eg$ orbitals being composed by $d{z^2}$ and $d{xy}$ instead of $d_{x^2-y^2}$. This is because another 45-degree rotation along the optimum z-axis is needed to get you exactly what you want.

While the optimum z-axis is easy to find using the auto rotation routine I've developed, the exact optimum sets of angles are harder to find. This is because $d{z^2}$ is not symmetrically related to any other d orbitals while $d{xy}$ and $d{x^2-y^2}$ orbitals and $d{xz}$ and $d_{yz}$ orbitals are related to each other by a rotation along the z-axis.

At the moment I'm only searching for a rotation frame that gives the largest TOTAL projection coefficient. Because of the symmetry relation I just mentioned, this might not give you the exact best rotation frame. I have some ideas to solve this but more playing around is needed.

Finally, in case you are also wondering, the small peak located at around -5 eV are due to hybridization with other orbitals.

Chengcheng-Xiao commented 4 months ago

Now sitting down and thinking about this, I'm not sure another 45-degree rotation along the optimum z-axis should be applied.

Because I'm using the square of the projection coefficients, any addition rotation along the optimum z-axis should show up as a sizable reduction to the total projection coefficient. As such, I would tend to believe that my algorithm gives the correct results.

What I think is happening here is that $d{z^2}$ and $d{xy}$ form a set of two degenerate orbitals while $d{xz}$ and $d{yz}$ form another set of two degenerate orbitals. That leaves $d_{x^2-y^2}$ to form the final set:


                                       - - dz2 and dxy

 - - - - -  d orbitals    =>            -  dx^2-y2

                                       - - dxz dyz

These sets can be clearly seen on your plot as orbitals within the same set show the identical pDOS. Similar claim has been reported in: https://journals.aps.org/prl/pdf/10.1103/PhysRevLett.112.117601

Note that this is not your typical $eg$ and $t_2 g$ orbital sets, and I think might be caused by the site symmetry of your metal atom.

CrankyBot2001 commented 4 months ago

Yes, I get your points. Yes for a single atom from my system(total 9 $Ti$ atom) I able to distinguish each orbitals. What I showed you it was for close to $Ti^{4+}$ single atom. Now I know in my system some $Ti$ contains $+3$ , and some $+4$ states. So I want to select all $+4$ active $Ti$ atoms to plot $t_{2g}-e_g$ resolved pdos. So when I am selecting 2 $+4$ active atoms and plotted pdos by adding up those two atoms contribution, then I am getting like this, high_quality_plot(TiCrO2) So, to report in paper should I plot single atom resolved pdos or all $+4$ active pdos and $+3$ active pdos separately?

Chengcheng-Xiao commented 4 months ago

Based on the limited information I have. I would suggest you find Ti atoms with the same site symmetries (not just valence) and only plot a single representative pDOS for each category. I don't think plotting combined pDOS is a good idea because a 90-degree rotation along the local z-axis might still lead to interchangeable $d{xz}$ and $d{yz}$ labeling.

In the meantime, I recommend against using $eg$ and $t_2 g$ tables because the previous pDOS plot indicates other symmetries are more suited.

I'm sorry but I don't think I can provide any more help analyzing specific chemical properties of your system, as it is out of the scope of this issue and repo. If you require more help from me with your system, you can contact me via email and we can discuss possible solutions. Thanks.

CrankyBot2001 commented 4 months ago

Yes, now I am doing each octahedral pdos except sum over pdos. But what I found for some Ti4+ atom $d{xy}$ orbital under $t{2g}$ and for another Ti4+ atom I saw $d{x^2-y^2}$ orbital under $t{2g}$, I couldn't understood what is the the physical significance of these two orbitals interchange. Furthermore, after finding Euler's angle for $Cr^{+3}$, I got this PDOS, Cr2 here, $d^{x^2-y^2}$ orbital and $d_z^2$ are degenrate under $eg$ that's ok. But why $d{xy}$ orbital not present in occupied section and why this and $e_g$ orbitals are in same energy point?

Chengcheng-Xiao commented 4 months ago

I have two suggestions:

  1. $e_g$ and $t2g$ categorization is NOT a law. It is based strictly on octahedral symmetry, as such, any distortion on the octahedral or second-near neighbor might break this symmetry. E.g., in your case, $d{xy}$ and $d_{x^2-y^2}$ exchange places probably means the atoms have different site symmetry.
  2. The rotation will only work when d orbital Hilbert space is complete. If you find some orbital contributions are missing, most of the time you just need to include more empty bands.

Again, I think this issue is veering towards a scientific collaboration that's not very related to how the code behaves. And I think it is better to (potentially) continue this conversation via emails so proper attribution can be acknowledged.

Chengcheng-Xiao commented 4 months ago

Closed due to inactivity.