Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Next revision
Previous revision
forum:data:2023:directional_x-ray_spectroscopy_potential_bug [2023/04/07 23:56] – Created from the form at forum:start Charles Cardotforum:data:2023:directional_x-ray_spectroscopy_potential_bug [2023/04/08 01:00] (current) Charles Cardot
Line 7: Line 7:
 Hello, Hello,
  
-I have been experimenting with calculating direction dependent x-ray emission with Quanty with DFT+MLFT, but I seem to be running into an issue with the mean field operators calculated from the Tight Binding (TB) object.+I have been experimenting with calculating direction dependent x-ray emission using Quanty and the DFT+MLFT framework. However, I seem to be running into an issue with the mean field operators calculated from the Tight Binding (TB) object.
  
-The directional dependence is encoded in the DFT step, where I distort the crystal lattice slightly (~1%) along a particular axis. When I stretch along the a-axis, I get an HDFT that looks like+The directional dependence is encoded in the DFT step, where I distort the crystal lattice slightly (~2%) along a particular axis.
  
-a_stretched_HDFT: https://drive.google.com/file/d/1_pMPkm1Pxw_f7F8whgx8hkA1boeQX_5-/view?usp=share_link+   # lattice constants; 
 +   @l@ 4.27633759 4.18633759 4.18633759 
 +   # set axis angles 
 +   @a@ 90. 90. 90. 
 +   # setup Wyckoff positions 
 +   @n@2 
 +   # Now, give list of ALL !!! Wyckoff positions. 
 +   @1@ Ni @ 0.00000000  0.00000000  0.00000000 
 +   @2@ O  @ 0.00000000  0.00000000  0.50000000
  
-which makes sense in the context of the distorted octahedral crystal field. The *xy* and *xz* orbitals are still degenerate, while the *yz* is now lower energy. Also the z^2-r^2 and x^2-y^2 orbital become mixed. The MLFT calculation within Quanty proceeds much the same way as described in the 2022 and 2019 Heidelberg tutorials for DFT + MLFT, where the mean field components Coulomb interaction is subtracted from HDFT to avoid double counting. The full code can be found here (https://drive.google.com/drive/folders/1G0iSTfIMT8i39CWWg1StEs5Kj26ZYQ1W?usp=share_link), but an excerpt is shown below. 
  
-    -- DEFINE ALL HAMILTONIANS+When I stretch along the a-axis and go through the regular Wannierization process to create a tight binding Hamiltonian, I get an HDFT that looks like
  
-    print("--Define Intermediate State Hamiltonian--\n")+a_stretched_HDFT: https://drive.google.com/file/d/1sE72v2rbbPM9Bof4Hg_CgNgvqmtjKjYY/view?usp=share_link
  
 +which makes sense in the context of the distorted octahedral crystal field. The //xy// and //xz// orbitals are still degenerate, while the //yz// is now lower energy. Also the //z^2-r^2// and //x^2-y^2// orbitals become mixed. The MLFT calculation within Quanty proceeds much the same way as described in the 2022 and 2019 Heidelberg tutorials for DFT + MLFT, where the mean field components Coulomb interaction is subtracted from HDFT to avoid double counting. The full code that I am running can be found here (https://drive.google.com/drive/folders/1G0iSTfIMT8i39CWWg1StEs5Kj26ZYQ1W?usp=share_link), but an excerpt is shown below.
 +    
 +
 +    print("--Define Intermediate State Hamiltonian--\n")
     --Hamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4     --Hamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4
     Hamiltonian = HDFT - F0dd * OppF0MFDFT - F2dd * OppF2MFDFT - F4dd * OppF4MFDFT     Hamiltonian = HDFT - F0dd * OppF0MFDFT - F2dd * OppF2MFDFT - F4dd * OppF4MFDFT
Line 26: Line 37:
     Hamiltonian = Hamiltonian/2     Hamiltonian = Hamiltonian/2
     Hamiltonian = Hamiltonian + ConjugateTranspose(Hamiltonian)     Hamiltonian = Hamiltonian + ConjugateTranspose(Hamiltonian)
 +    
     esinter      = -nd*Usd     esinter      = -nd*Usd
     edinter      = -((-20*Delta + 19*nd*Udd+nd*nd*Udd+40*Usd)/(2*(10+nd)))     edinter      = -((-20*Delta + 19*nd*Udd+nd*nd*Udd+40*Usd)/(2*(10+nd)))
     eLinter      = nd*((1+nd)*Udd/2 - Delta + 2*Usd)/(10+nd)     eLinter      = nd*((1+nd)*Udd/2 - Delta + 2*Usd)/(10+nd)
 +    
     OperatorSetTrace(Hamiltonian,esinter,Index["Ni_1s"])     OperatorSetTrace(Hamiltonian,esinter,Index["Ni_1s"])
     OperatorSetTrace(Hamiltonian,edinter,Index["Ni_3d"])     OperatorSetTrace(Hamiltonian,edinter,Index["Ni_3d"])
     OperatorSetTrace(Hamiltonian,eLinter,Index["Ligand_d"])     OperatorSetTrace(Hamiltonian,eLinter,Index["Ligand_d"])
- +    
     print("--Define XES-Hamiltonian--\n")     print("--Define XES-Hamiltonian--\n")
     --XESHamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4     --XESHamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4
Line 45: Line 55:
     XESHamiltonian = XESHamiltonian/2     XESHamiltonian = XESHamiltonian/2
     XESHamiltonian = XESHamiltonian + ConjugateTranspose(XESHamiltonian)     XESHamiltonian = XESHamiltonian + ConjugateTranspose(XESHamiltonian)
 +    
     epfinal      = -nd*Upd     epfinal      = -nd*Upd
     edfinal      = -((-20*Delta + 19*nd*Udd + nd*nd*Udd + 120*Upd)/(2*(10 + nd)))     edfinal      = -((-20*Delta + 19*nd*Udd + nd*nd*Udd + 120*Upd)/(2*(10 + nd)))
     eLfinal      = nd*(-2*Delta + Udd + nd*Udd + 12*Upd)/(2*(10 + nd))     eLfinal      = nd*(-2*Delta + Udd + nd*Udd + 12*Upd)/(2*(10 + nd))
 +    
     OperatorSetTrace(XESHamiltonian,epfinal,Index["Ni_2p"])     OperatorSetTrace(XESHamiltonian,epfinal,Index["Ni_2p"])
     OperatorSetTrace(XESHamiltonian,edfinal,Index["Ni_3d"])     OperatorSetTrace(XESHamiltonian,edfinal,Index["Ni_3d"])
     OperatorSetTrace(XESHamiltonian,eLfinal,Index["Ligand_d"])     OperatorSetTrace(XESHamiltonian,eLfinal,Index["Ligand_d"])
 +    
     print("--Compute eigenstates--\n")     print("--Compute eigenstates--\n")
     -- we now can create the lowest Npsi eigenstates:     -- we now can create the lowest Npsi eigenstates:
Line 60: Line 70:
     {DeterminantString(NFermi,Index["Ni_3d"],Index["Ligand"]),16+nd,16+nd},     {DeterminantString(NFermi,Index["Ni_3d"],Index["Ligand"]),16+nd,16+nd},
     {DeterminantString(NFermi,Index["Ni_2p"]),6,6}}     {DeterminantString(NFermi,Index["Ni_2p"]),6,6}}
- 
     psiList = Eigensystem(Hamiltonian, StartRestrictions, Npsi, {{'Zero',1e-12},{'Epsilon',1e-12}})     psiList = Eigensystem(Hamiltonian, StartRestrictions, Npsi, {{'Zero',1e-12},{'Epsilon',1e-12}})
     psiList = Chop(psiList)     psiList = Chop(psiList)
     print(StartRestrictions)     print(StartRestrictions)
 +    
     Npsiline = {}     Npsiline = {}
     for i=1,Npsi,1     for i=1,Npsi,1
Line 101: Line 110:
     print("Finished")     print("Finished")
  
-When the output is plotted it looks like+ 
 +The plotted output looks like
  
 XES_including_MFDFT: https://drive.google.com/file/d/1hYutZdeujtK8lohFPPD7sZ-7z5EN4asl/view?usp=share_link XES_including_MFDFT: https://drive.google.com/file/d/1hYutZdeujtK8lohFPPD7sZ-7z5EN4asl/view?usp=share_link
  
-where clearly the x, y, and z polarized emissions are different. I would expect the y and z spectra to be identical, with the x spectra being slightly different, as the distortion is along the a-axis. In investigating this, I noticed that omitting the subtraction of the Mean Field DFT contribution, +where clearly the x, y, and z polarized emissions are different. Given that the distortion is along the a-axis, I expect the y and z spectra to be identical, with the x spectra being slightly different. In investigating this, I noticed that omitting the subtraction of the Mean Field DFT contribution,
- +
-print("--Define Intermediate State Hamiltonian--\n")+
  
 +    print("--Define Intermediate State Hamiltonian--\n")
     Hamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4     Hamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4
     --Hamiltonian = HDFT - F0dd * OppF0MFDFT - F2dd * OppF2MFDFT - F4dd * OppF4MFDFT     --Hamiltonian = HDFT - F0dd * OppF0MFDFT - F2dd * OppF2MFDFT - F4dd * OppF4MFDFT
Line 116: Line 125:
     Hamiltonian = Hamiltonian/2     Hamiltonian = Hamiltonian/2
     Hamiltonian = Hamiltonian + ConjugateTranspose(Hamiltonian)     Hamiltonian = Hamiltonian + ConjugateTranspose(Hamiltonian)
 +    
     esinter      = -nd*Usd     esinter      = -nd*Usd
     edinter      = -((-20*Delta + 19*nd*Udd+nd*nd*Udd+40*Usd)/(2*(10+nd)))     edinter      = -((-20*Delta + 19*nd*Udd+nd*nd*Udd+40*Usd)/(2*(10+nd)))
     eLinter      = nd*((1+nd)*Udd/2 - Delta + 2*Usd)/(10+nd)     eLinter      = nd*((1+nd)*Udd/2 - Delta + 2*Usd)/(10+nd)
 +    
     OperatorSetTrace(Hamiltonian,esinter,Index["Ni_1s"])     OperatorSetTrace(Hamiltonian,esinter,Index["Ni_1s"])
     OperatorSetTrace(Hamiltonian,edinter,Index["Ni_3d"])     OperatorSetTrace(Hamiltonian,edinter,Index["Ni_3d"])
     OperatorSetTrace(Hamiltonian,eLinter,Index["Ligand_d"])     OperatorSetTrace(Hamiltonian,eLinter,Index["Ligand_d"])
- +    
     print("--Define XES-Hamiltonian--\n")     print("--Define XES-Hamiltonian--\n")
     XESHamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4     XESHamiltonian = HDFT + F0dd * OppF0 + F2dd * OppF2 + F4dd * OppF4
Line 135: Line 143:
     XESHamiltonian = XESHamiltonian/2     XESHamiltonian = XESHamiltonian/2
     XESHamiltonian = XESHamiltonian + ConjugateTranspose(XESHamiltonian)     XESHamiltonian = XESHamiltonian + ConjugateTranspose(XESHamiltonian)
 +    
     epfinal      = -nd*Upd     epfinal      = -nd*Upd
     edfinal      = -((-20*Delta + 19*nd*Udd + nd*nd*Udd + 120*Upd)/(2*(10 + nd)))     edfinal      = -((-20*Delta + 19*nd*Udd + nd*nd*Udd + 120*Upd)/(2*(10 + nd)))
     eLfinal      = nd*(-2*Delta + Udd + nd*Udd + 12*Upd)/(2*(10 + nd))     eLfinal      = nd*(-2*Delta + Udd + nd*Udd + 12*Upd)/(2*(10 + nd))
 +    
     OperatorSetTrace(XESHamiltonian,epfinal,Index["Ni_2p"])     OperatorSetTrace(XESHamiltonian,epfinal,Index["Ni_2p"])
     OperatorSetTrace(XESHamiltonian,edfinal,Index["Ni_3d"])     OperatorSetTrace(XESHamiltonian,edfinal,Index["Ni_3d"])
     OperatorSetTrace(XESHamiltonian,eLfinal,Index["Ligand_d"])     OperatorSetTrace(XESHamiltonian,eLfinal,Index["Ligand_d"])
  
-the issue went away and the spectra exhibited the expected behavior, with the y and z polarized emission lying exactly on top of each otherand the x emission being distinct from the other two directions.+causes the issue to go away and the spectra exhibits the expected behavior, 
 + 
 +XES_excluding_MFDFT: https://drive.google.com/file/d/157Vyqxn2Ld20tajT0w3anxgFYQB9p8xu/view?usp=share_link,
  
-XES_excluding_MFDFT: https://drive.google.com/file/d/157Vyqxn2Ld20tajT0w3anxgFYQB9p8xu/view?usp=share_link+with the y and z polarized emission lying exactly on top of each other, and the x emission being distinct from the other two directions.
  
-What this tells me is that the crystal structure that is reflected in the Tight Binding Hamiltonain (HDFTis indeed correct, and produces the expected directional dependence, but something about the Coulomb correction is not behaving as expected. I'm not sure if this is an issue of numerical accuracy (the distortion is small) or a potential bug in the implementation of the MeanFieldOperator. Any help that can be provided would be much appreciated.+What this tells me is that the crystal structure that is reflected in the tight binding HamiltonainHDFTis indeed correct, and produces the expected directional dependence. However, something about the Coulomb correction is not behaving as expected. I'm not sure if this is an issue of numerical accuracy (the distortion is small) or a potential bug in the implementation of the MeanFieldOperator. Any help that can be provided would be much appreciated.
 </WRAP> </WRAP>
  
 ~~DISCUSSION|Answers~~ ~~DISCUSSION|Answers~~
  
Print/export