Differences

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

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision
Previous revision
forum:data:2020:stranger_things [2020/01/07 19:08] Riccardo Piomboforum:data:2020:stranger_things [2020/01/07 20:21] (current) Maurits W. Haverkort
Line 1: Line 1:
-====== STRANGER THINGS ======+====== Wrong (missing) ground state Part 2 ======
 ;;# ;;#
 asked by [[mailto:riccardo.piombo@gmail.com|Riccardo Piombo]] (2020/01/07 18:57) asked by [[mailto:riccardo.piombo@gmail.com|Riccardo Piombo]] (2020/01/07 18:57)
Line 6: Line 6:
 <WRAP center box 100%> <WRAP center box 100%>
 Hi every one,  Hi every one, 
-I wrote a simple code to calculate the chemical potential for the addition and for the removal of a particle. The system under consideration is a planar D4h structure composed by a transition metal surrounded by 4 non metal ligands. The problem is the following: when computing the ground state of the N-1 system Quanty seems to output ciclycally three different value for the GS energy which differ from one another by 1.6 eV while for the N and N+1 case it doesn't change (as it must be!). Here is the code could someone copy and paste it on his pc and just run it too check this strange thing!?+I wrote a simple code to calculate the chemical potential for the addition and for the removal of a particle. The system under consideration is a planar D4h structure composed of a transition metal surrounded by 4 non-metal ligands. The problem is the following: when computing the ground state of the N-1 system Quanty seems to output cyclically three different value for the GS energy which differs by almost 1.6 eV in the extreme case while for the N and N+1 case it doesn't change (as it must be!). Here is the code could someone copy and paste it on his pc and just run it too check this strange thing!?
  
 ---------------------------------------------------------------------------------------------- ----------------------------------------------------------------------------------------------
-"Verbosity(0)+<code>
  
-NF=20\\ +Verbosity(0)
-NB=0\\+
  
--- Silver 4d spin-orbitals\\ +NF=20 
-IndexDn_4d=0, 2, 4, 6, 8}\\ +NB=0
-IndexUp_4d={ 1, 3, 5, 7, 9}\\ +
--- Fluorine 2p spin-orbitals\\ +
-IndexDn_Ld={10,12,14,16,18}\\ +
-IndexUp_Ld={11,13,15,17,19}\\+
  
--- number of electrons in d-shell and Ligand-P shell\\ +-- Silver 4d spin-orbitals 
-n4d = 9\\ +IndexDn_4d={ 0, 2, 4, 6, 8} 
-nL = 10\\+IndexUp_4d={ 1, 3, 5, 7, 9} 
 +-- Fluorine 2p spin-orbitals 
 +IndexDn_Ld={10,12,14,16,18} 
 +IndexUp_Ld={11,13,15,17,19}
  
-Udd = 5.0\\ +-- number of electrons in d-shell and Ligand-P shell 
-Delta_pd =  1.5\\ +n4d 9 
-Tpp  0.4 \\     +nL 10
-TdLb1g 2.6  \\  +
  
--- Racah parameters \\ +Udd 5.0 
-= 0.09\\ +Delta_pd =  1.5 
-= 0.54\\+Tpp  0.4       
 +TdLb1g = 2.6    
  
--- Slater-Koster integrals for monopole and \\ +-- Racah parameters  
--- multipole part of Coulomb interaction\\ +B = 0.09 
-F2dd = 49.0*+ 7*C\\ += 0.54
-F4dd 441.0*C/35.0\\ +
-F0dd Udd + (F2dd + F4dd)*2.0/63.0  \\+
  
--- Onsite splitting of the 4d shell \\ +-- Slater-Koster integrals for monopole and  
-Eda1g  0.0\\ +-- multipole part of Coulomb interaction  
-Edb1g  0.0\\ +F2dd 49.0*B + 7*C 
-Edb2g  0.0\\ +F4dd 441.0*C/35.0 
-Edeg  =  0.0\\+F0dd Udd + (F2dd + F4dd)*2.0/63.0  
  
--- Onsite splitting of the Ligand P shell \\ +-- Onsite splitting of the 4d shell  
-ELa1g -Delta_pd - Tpp \\ +Eda1g  0.0 
-ELb1g -Delta_pd + Tpp\\ +Edb1g  0.0 
-ELb2g -Delta_pd - Tpp\\ +Edb2g  0.0 
-ELeg  -Delta_pd\\+Edeg   0.0
  
--- Hopping between the Ligand and 4d shell\\  +-- Onsite splitting of the Ligand shell  
-TdLa1g TdLb1g/sqrt(3.0)\\ +ELa1g -Delta_pd - Tpp 
-TdLb2g TdLb1g/2.0\\ +ELb1g -Delta_pd + Tpp 
-TdLeg  TdLb1g/(2.0*sqrt(2.0))\\+ELb2g = -Delta_pd - Tpp 
 +ELeg  -Delta_pd
  
--- turning U and Delta to onsite energies\\  +-- Hopping between the Ligand and 4d shell  
-Delta_ZSA Delta_pd + Tpp/5.0\\ +TdLa1g TdLb1g/sqrt(3.0) 
-eL (Udd*n4d*(1-n4d)/2.0 - n4d*(Delta_ZSA -Udd*n4d))/(n4d+nL)\\ +TdLb2g TdLb1g/2.0 
-ed = eL + Delta_ZSA - Udd*n4d \\+TdLeg  = TdLb1g/(2.0*sqrt(2.0))
  
--- D4h Crystal field on 4d states\\+-- turning U and Delta to onsite energies  
 +Delta_ZSA = Delta_pd + Tpp/5.0 
 +eL = (Udd*n4d*(1-n4d)/2.0 - n4d*(Delta_ZSA -Udd*n4d))/(n4d+nL) 
 +ed = eL + Delta_ZSA - Udd*n4d  
 + 
 +-- D4h Crystal field on 4d states
 Akm_4d = {{0, 0, ed} , Akm_4d = {{0, 0, ed} ,
-          {2, 0, Eda1g + (-1)*(Edb1g) + (-1)*(Edb2g) + Edeg} , +{2, 0, Eda1g + (-1)*(Edb1g) + (-1)*(Edb2g) + Edeg} , 
-          {4, 0, (3/10)*((6)*(Eda1g) + Edb1g + Edb2g + (-8)*(Edeg))} , +{4, 0, (3/10)*((6)*(Eda1g) + Edb1g + Edb2g + (-8)*(Edeg))} , 
-          {4,-4, (3/2)*((sqrt(7/10))*(Edb1g + (-1)*(Edb2g)))} , +{4,-4, (3/2)*((sqrt(7/10))*(Edb1g + (-1)*(Edb2g)))} , 
-          {4, 4, (3/2)*((sqrt(7/10))*(Edb1g + (-1)*(Edb2g)))} }+{4, 4, (3/2)*((sqrt(7/10))*(Edb1g + (-1)*(Edb2g)))} }
  
-OppCF_4d = NewOperator("CF", NF, IndexUp_4d, IndexDn_4d, Akm_4d)\\ +OppCF_4d = NewOperator("CF", NF, IndexUp_4d, IndexDn_4d, Akm_4d) 
-OppCF_4d.Name = "Crystal Field on 4d-states"\\+OppCF_4d.Name = "Crystal Field on 4d-states"
  
--- D4h Crystal field on Ligand-P states\\+-- D4h Crystal field on Ligand-P states
 Akm_L = {{0, 0, eL}, Akm_L = {{0, 0, eL},
-         {2, 0, ELa1g + (-1)*(ELb1g) + (-1)*(ELb2g) + ELeg} , +{2, 0, ELa1g + (-1)*(ELb1g) + (-1)*(ELb2g) + ELeg} , 
-         {4, 0, (3/10)*((6)*(ELa1g) + ELb1g + ELb2g + (-8)*(ELeg))} , +{4, 0, (3/10)*((6)*(ELa1g) + ELb1g + ELb2g + (-8)*(ELeg))} , 
-         {4,-4, (3/2)*((sqrt(7/10))*(ELb1g + (-1)*(ELb2g)))} , +{4,-4, (3/2)*((sqrt(7/10))*(ELb1g + (-1)*(ELb2g)))} , 
-         {4, 4, (3/2)*((sqrt(7/10))*(ELb1g + (-1)*(ELb2g)))} }+{4, 4, (3/2)*((sqrt(7/10))*(ELb1g + (-1)*(ELb2g)))} }
  
-OppCF_Ld = NewOperator("CF", NF, IndexUp_Ld, IndexDn_Ld, Akm_L)\\ +OppCF_Ld = NewOperator("CF", NF, IndexUp_Ld, IndexDn_Ld, Akm_L) 
-OppCF_Ld.Name = "Crystal Field on Ligand-P states"\\+OppCF_Ld.Name = "Crystal Field on Ligand-P states"
  
--- Hopping t_pd\\+-- Hopping t_pd
 Akm_dL = {{0, 0, (1/5)*(TdLa1g + TdLb1g + TdLb2g + (2)*(TdLeg))} , Akm_dL = {{0, 0, (1/5)*(TdLa1g + TdLb1g + TdLb2g + (2)*(TdLeg))} ,
          {2, 0, TdLa1g + (-1)*(TdLb1g) + (-1)*(TdLb2g) + TdLeg} ,          {2, 0, TdLa1g + (-1)*(TdLb1g) + (-1)*(TdLb2g) + TdLeg} ,
Line 90: Line 92:
  
  
-OppHopLd4d  = NewOperator("CF", NF, IndexUp_4d,IndexDn_4d, IndexUp_Ld,IndexDn_Ld,Akm_dL) + NewOperator("CF", NF, IndexUp_Ld,IndexDn_Ld, IndexUp_4d,IndexDn_4d,Akm_dL)\\ +OppHopLd4d  = NewOperator("CF", NF, IndexUp_4d,IndexDn_4d, IndexUp_Ld,IndexDn_Ld,Akm_dL) +   
-OppHopLd4d.Name = "P-d hybridization"\\+              NewOperator("CF", NF, IndexUp_Ld,IndexDn_Ld, IndexUp_4d,IndexDn_4d,Akm_dL) 
 +OppHopLd4d.Name = "P-d hybridization"
  
--- Udd repulsion \\ +-- Udd repulsion  
-OppF0 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {1,0,0})\\ +OppF0 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {1,0,0}) 
-OppF2 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {0,1,0})\\ +OppF2 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {0,1,0}) 
-OppF4 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {0,0,1})\\ +OppF4 = NewOperator("U", NF, IndexUp_4d, IndexDn_4d, {0,0,1}) 
-OppUdd = F0dd*OppF0 + F2dd*OppF2 + F4dd*OppF4\\+OppUdd = F0dd*OppF0 + F2dd*OppF2 + F4dd*OppF4
  
--- Silver hamiltonian\\ +-- Silver hamiltonian 
-H_Ag = OppUdd + OppCF_4d \\+H_Ag = OppUdd + OppCF_4d 
  
--- Fluorine hamiltonian\\ +-- Fluorine hamiltonian 
-H_F = OppCF_Ld\\+H_F = OppCF_Ld
  
--- Ag + F Hamiltonian\\ +-- Ag + F Hamiltonian 
-H_Cluster = H_Ag + H_F\\ +H_Cluster = H_Ag + H_F 
-H_tot = H_Cluster + (-1)*OppHopLd4d\\+H_tot = H_Cluster + (-1)*OppHopLd4d
  
--- Restrictions = {filling of nd electrons in 4d states, filling of nL electrons in Ligand-P states}\\ +-- Restrictions = {filling of nd electrons in 4d states, filling of nL electrons in Ligand-P states} 
-print("")\\ +print(""
-print("N-body states energies computation")\\ +print("N-body states energies computation"
-Npsi_i = 1\\ +Npsi_i = 1 
-StartRestrictions1 = {NF, NB, {"1111111111 0000000000",n4d,n4d}, {"0000000000 1111111111",nL,nL}}\\ +StartRestrictions1 = {NF, NB, {"1111111111 0000000000",n4d,n4d}, {"0000000000 1111111111",nL,nL}} 
-psiList_N = Eigensystem(H_tot, StartRestrictions1, Npsi_i)\\ +psiList_N = Eigensystem(H_tot, StartRestrictions1, Npsi_i) 
---psiList_N = {psiList_N}\\ +psiList_N = {psiList_N} 
-print("Done")\\+print("Done")
  
-E_gs_N = psiList_N[1]*(H_tot)*psiList_N[1]\\ +E_gs_N = psiList_N[1]*(H_tot)*psiList_N[1] 
-print("")\\ +print(""
-print("Egs N particle = ", E_gs_N)\\+print("Egs N particle = ", E_gs_N)
  
  
--- initialization of the variables\\ +-- initialization of the variables 
-remove_mu = 0\\ +remove_mu = 0 
-add_mu  = 0\\ +add_mu  = 0 
-mean_mu = 0\\ +mean_mu = 0 
-E_gs_Nminus1 = 0\\ +E_gs_Nminus1 = 0 
-E_gs_Nplus1 =  0\\+E_gs_Nplus1 =  0
  
---################################################################################## +-- chemical potential μ- to remove a particle 
--- HERE SEEMS TO BE THE PROBLEM +-- ########################################################### 
--- chemical potential μ- to remove a particle\\ +-- THE PROBLEM SEEMS TO BE HERE 
-print("")\\ +print(""
-print("(N-1)-body states energies computation")\\ +print("(N-1)-body states energies computation"
-StartRestrictions_rem = {NF, NB, {"1111111111 0000000000",n4d-1,n4d-1}, {"0000000000 1111111111",nL,nL}}\\ +StartRestrictions_rem = {NF, NB, {"1111111111 0000000000",8,8}, {"0000000000 1111111111",10,10}} 
-psiList_Nminus1 = Eigensystem(H_tot, StartRestrictions_rem, 1)\\ +psiList_Nminus1 = Eigensystem(H_tot, StartRestrictions_rem, 1) 
-psiList_Nminus1 = {psiList_Nminus1}\\ +psiList_Nminus1 = {psiList_Nminus1} 
-print("Done")\\+print("Done")
  
-E_gs_Nminus1 = psiList_Nminus1[1]*(H_tot)*psiList_Nminus1[1]\\ +E_gs_Nminus1 = psiList_Nminus1[1]*(H_tot)*psiList_Nminus1[1] 
-remove_mu = E_gs_N - E_gs_Nminus1  \\ +remove_mu = E_gs_N - E_gs_Nminus1   
-print("N minus 1", E_gs_Nminus1)\\ +print("N minus 1", E_gs_Nminus1) 
-print("μ- = ", remove_mu)\\ +print("μ- = ", remove_mu) 
-print("")\\ +print(""
---##################################################################################+-- ###########################################################
  
--- chemical potential μ+ to add a particle\\ +-- chemical potential μ+ to add a particle 
-print("(N+1)-body states energies computation")\\ +print("(N+1)-body states energies computation"
-StartRestrictions_add = {NF, NB, {"1111111111 0000000000",n4d+1,n4d+1}, {"0000000000 1111111111",nL,nL}}\\ +StartRestrictions_add = {NF, NB, {"1111111111 0000000000",10,10}, {"0000000000 1111111111",10,10}} 
-psiList_Nplus1 = Eigensystem(H_tot, StartRestrictions_add, 1)\\ +psiList_Nplus1 = Eigensystem(H_tot, StartRestrictions_add, 1) 
-psiList_Nplus1 = {psiList_Nplus1}\\+psiList_Nplus1 = {psiList_Nplus1}
 print("Done") print("Done")
-\\ 
-E_gs_Nplus1 = psiList_Nplus1[1]*(H_tot)*psiList_Nplus1[1]\\ 
-add_mu = E_gs_Nplus1 - E_gs_N\\ 
-print("N plus 1", E_gs_Nplus1)\\ 
-print("μ+ = ", add_mu)\\ 
-print("")\\ 
  
-print("mean μ = [(μ+) + (μ-)]/2 = ", mean_mu)\\ +E_gs_Nplus1 = psiList_Nplus1[1]*(H_tot)*psiList_Nplus1[1] 
-"+add_mu = E_gs_Nplus1 - E_gs_N 
 +print("N plus 1", E_gs_Nplus1) 
 +print("μ+ = ", add_mu) 
 +print(""
 + 
 +-- mean μ  --> Fermi energy is not defined in the experimental values 
 +-- so we have to fix it and then freely translate the experimental spectrum along w axis 
 +mean_mu = (add_mu + remove_mu)/2.0 
 +print("mean μ = [(μ+) + (μ-)]/2 = ", mean_mu) 
 + 
 +</code> 
 + 
 +P.S: after some simulations, the problem seems to appear only when n4d-1 is equal to 8 
 +the GS energy becomes -15.628822671742, -14.433720341628 or -14.042673360744  
 </WRAP> </WRAP>
  
 ~~DISCUSSION|Answers~~ ~~DISCUSSION|Answers~~
  
Print/export