Dear Shanyi, just a quick question: have you considered the Basis Set 
Superposition Error (BSSE) in your SIESTA calculations?

----- Original Message -----
From: "shangy" <[email protected]>
To: [email protected]
Sent: Wednesday, April 24, 2013 3:20:16 PM
Subject: [SIESTA-L] The noticeable discrepancies between SIESTA and VASP 
calculations



Dear Siesta developers and users : 



(I am sorry for the same mail again. The last one was rejected by the mail 
server. And the pseudopotentails were not attached in last mail) 




I found there are noticeable discrepancies between SIESTA and VASP calculations 
when I calculated the formation energies of binary compounds in Mg-Y-Zn system, 
such as MgZn 2 , Mg 2 Zn 11 , Mg 24 Y 5 , Mg 2 Y , MgY , Zn 12 Y , Zn 17 Y 2 , 
Zn 3 Y and ZnY. The VASP calculations are consistent with the published results 
obtained by VASP. However, the formation energies calculated by SIESTA for 
these compounds are not comparable to these available results, as showed in 
following table 1. The formation energies (per atom) of various compounds with 
composition of Mg x Zn y Y z are calculated by: 

Delta_E f = E(Mg x Zn y Y z ) - E(Mg)*x/(x+y+x) - E(Zn)*y/(x+y+x) - 
E(Y)*z/(x+y+x), 

where E(Mg x Zn y Y z ) is the total energy of per atom of Mg x Zn y Y z , and 
E( i , i =Mg, Zn, Y) is the energy of per i atom in bulk. 



Table 1. The formation energies calculated by SIESTA and VASP for various 
binary compounds 
        

Formular        

Delta_E f (meV/atom) 
        

VASP    

SIESTA  

Scripta Mater., 67 (2012) 798 (with VASP) 
        

MgZn 2  

-140.8  

-84.8   

-138 
        

Mg 2 Zn 11      

-71.7   

-42.3   

-68 
        

Mg 24 Y 5       

-59.3   

-90.0   

-59 
        

Mg 2 Y  

-93.5   

-129.4  

-93 
        

MgY     

-109.8  

-140.6  

-109 
        

Zn 12 Y         

-187.8  

-196.5  

-184 
        

Zn 17 Y 2       

-244.5  

-259.9  

-241 
        

Zn 3 Y  

-373.6  

-462.1  

-370 
        

ZnY     

-366.2  

-473.6  

-363 



The PAW-PBE pseudopotentials supplied by VASP and an energy cutoff of 300 eV 
were used in all VASP calculations. The attached Troullier-Martins PBE 
pseudopotentials of Mg, Zn and Y were used in SIESTA calculations. The 
pseudopotentials of Mg, Zn and Y were validated via the their lattices and 
cohesive energies, as showed in following table, which all are comparable to 
experimental data (the values within the bracket in table 2). 



Table 2. The lattices and cohesive energies calculated by SIESTA for Mg, Zn and 
Y bulk, the experimental data are given in the bracket. 
        

Pseudopotentials        

a (angstrom)    

c (angstrom)    

Cohesive (eV) 
        

Mg      

3.234 (3.209)   

5.158 (5.211)   

1.447 (1.51) 
        

Zn      

2.697 (2.665)   

5.184 (4.947)   

1.434 (1.35) 
        

Y       

3.661 (3.648)   

5.729 (5.732)   

4.543 (4.37) 



The k-point meshes and convergence criteria of SIESTA calculations are 
identical to that of VASP calculations. The SIESTA input files for all 
calculations are similar except the k-points meshes and their atom coordinates. 
The optimzed structures and energies of various compounds are obtained with 
MD.VariableCell method and by restarting from previous calculations at least 
two times. Taking input file of MgZn 2 as an example, its input file is 
following: 



------------------siesta-input-file-for-MgZn 2 --------------------- 

SystemName MgZn2-vc-opt-re2 

SystemLabel MgZn2 

NumberOfSpecies 2 

NumberOfAtoms 12 



XC.functional GGA 

XC.authors PBE 

MaxSCFIterations 150 



#DM.UseSaveDM T 

#MD.UseSaveXV T 

#MD.UseSaveCG T 

#DM.MixSCF1 T 



DM.MixingWeight 0.01 

DM.NumberPulay 5 

DM.Tolerance 1.0E-4 

DM.Require.Energy.Convergence T 

DM.Energy.Tolerance 1.0E-5 eV 



NumberOfEigenStates 400 

Use.New.Diagk T 



MeshCutoff 400 Ry 

MeshSubDivision 6 



SolutionMethod diagon 

OccupationFunction MP 

#OccupationMPOrder 3 

ElectronicTemperature 1000 K 



Diag.ParallelOverK F 



MD.TypeOfRun CG 

MD.VariableCell T 

MD.NumCGsteps 500 

MD.MaxForceTol 0.01 eV/Ang 

MD.MaxStressTol 0.05 Gpa 

MD.TargetPressure 0.0 Gpa 



WriteForces T 

WriteMDxmol T 

WriteCoorXmol T 

WriteCoorStep T 



PAO.BasisSize DZP 

PAO.EnergyShift 80 meV 



%block PAO.Basis 

Zn 2 

n=4 0 2 P 

0.0 0.0 

n=3 2 2 

0.0 0.0 

%endblock PAO.Basis 



%block ChemicalSpeciesLabel 

1 12 Mg 

2 30 Zn 

%endblock ChemicalSpeciesLabel 



LatticeConstant 1 Ang 

%block LatticeVectors 

5.2512841225 0.0000000000 0.0000000000 

-2.6256420612 4.5477454525 0.0000000000 

0.0000000000 0.0000000000 8.4452314377 

%endblock LatticeVectors 



%block kgrid_Monkhorst_Pack 

8 0 0 0.5 

0 8 0 0.5 

0 0 5 0 

%endblock kgrid_Monkhorst_Pack 



AtomicCoordinatesFormat Fractional 

%block AtomicCoordinatesAndAtomicSpecies 

0.33333299 0.66666698 0.06241800 1 

0.66666698 0.33333299 0.56241798 1 

0.66666698 0.33333299 0.93758202 1 

0.33333299 0.66666698 0.43758199 1 

0.00000000 0.00000000 0.00000000 2 

0.00000000 0.00000000 0.50000000 2 

0.17036100 0.34072199 0.75000000 2 

0.82963902 0.17036100 0.25000000 2 

0.34072199 0.17036100 0.25000000 2 

0.65927798 0.82963902 0.75000000 2 

0.17036100 0.82963902 0.75000000 2 

0.82963902 0.65927798 0.25000000 2 


%endblock AtomicCoordinatesAndAtomicSpecies 





And for Mg 2 Y, its input file is following: 

------------------siesta-input-file-for-Mg 2 Y--------------------- 

SystemName Mg2Y-vc-opt-re2 

SystemLabel Mg2Y 

NumberOfSpecies 2 

NumberOfAtoms 12 



XC.functional GGA 

XC.authors PBE 

MaxSCFIterations 150 



#DM.UseSaveDM T 

#MD.UseSaveXV T 

#MD.UseSaveCG T 

#DM.MixSCF1 T 



DM.MixingWeight 0.01 

DM.NumberPulay 5 

DM.Tolerance 1.0E-4 

DM.Require.Energy.Convergence T 

DM.Energy.Tolerance 1.0E-5 eV 



NumberOfEigenStates 400 

Use.New.Diagk T 



MeshCutoff 400 Ry 

MeshSubDivision 6 



SolutionMethod diagon 

OccupationFunction MP 

#OccupationMPOrder 3 

ElectronicTemperature 1000 K 



Diag.ParallelOverK F 



MD.TypeOfRun CG 

MD.VariableCell T 

MD.NumCGsteps 500 

MD.MaxForceTol 0.01 eV/Ang 

MD.MaxStressTol 0.05 Gpa 

MD.TargetPressure 0.0 Gpa 



WriteForces T 

WriteMDxmol T 

WriteCoorXmol T 

WriteCoorStep T 



PAO.BasisSize DZP 

PAO.EnergyShift 80 meV 



%block ChemicalSpeciesLabel 

1 12 Mg 

2 39 Y 

%endblock ChemicalSpeciesLabel 



LatticeConstant 1 Ang 

%block LatticeVectors 

6.0594434738 0.0000000000 0.0000000000 

-3.0297217369 5.2476319811 0.0000000000 

0.0000000000 0.0000000000 9.8398962021 

%endblock LatticeVectors 



%block kgrid_Monkhorst_Pack 

7 0 0 0 

0 7 0 0 

0 0 5 0 

%endblock kgrid_Monkhorst_Pack 



AtomicCoordinatesFormat Fractional 

%block AtomicCoordinatesAndAtomicSpecies 

0.00000000 0.00000000 0.00000000 1 

0.00000000 0.00000000 0.50000000 1 

0.82894897 0.65789801 0.25000000 1 

0.17105100 0.82894897 0.75000000 1 

0.65789801 0.82894897 0.75000000 1 

0.34210199 0.17105100 0.25000000 1 

0.82894897 0.17105100 0.25000000 1 

0.17105100 0.34210199 0.75000000 1 

0.33333299 0.66666698 0.06474100 2 

0.66666698 0.33333299 0.56474101 2 

0.66666698 0.33333299 0.93525899 2 

0.33333299 0.66666698 0.43525901 2 

%endblock AtomicCoordinatesAndAtomicSpecies 

--------------------------------------------------------------------------------------------
 



I want to know why there are such noticeable discrepancies between SIESTA and 
VASP calculations on formation energies for these compounds. Did I miss some 
important settings or assign some unreasonable parameters in the input files? 
Or do these discrepancies arise from the description of interactions among Mg, 
Zn and Y atoms, which are different that of VASP. 

Any ideas and suggestions are highly appreciated! 

Thanks in advance! 



Shanyi Ma 

Institute of Metal Research, CAS 

Responder a