[Pw_forum] Application of LDA+U to Ce psp with zero 4f occupation

Thu Dec 15 04:22:34 CET 2011

I am trying to apply Hubbard U (DFT+U) to Ce in cerium based oxides
including CeO2 and Ce2O3. However, as you see from header of Ce psp (see [I]
below), our Ce psp (RRKJ, generated by OPIUM) is designed to have zero 4f
occupation. When I applied +U to Ce using QE V4.3.1, I got following error
message:

from offset_atom_wfc : error #     9

wrong offset

So, I searched previous threads and I figured out that I needed to modify QE
code (http://www.democritos.it/pipermail/pw_forum/2010-July/017622.html). 

I did modification and recompiled QE code. It looked fine but I found weird
energetics for reduction of CeO2 to Ce2O3 in comparison with data from PRB
75, 035109 (2007). So, I looked into output files and found out that
occupation matrix for both CeO2 and Ce2O3 are somewhat unphysical (at least
for me):

*U_eff=0.01 eV 

Tr[ns(na)] 

CeO2  = 1.1690424

Ce2O3 = 1.3908759

*U_eff=3.0 eV 

Tr[ns(na)] 

CeO2 =  0.9963736

Ce2O3 = 1.2650403 

*U_eff=9.0 eV 

Tr[ns(na)] 

CeO2 =  0.5972449 

Ce2O3 = 1.1056369

I thought that even with non-orthogonal terms, at least Tr[ns(na)] for CeO2
should be very close to 0.0 as Ce exists as Ce4+ in CeO2 (with ionic
picture) and should have empty 4f state. With analyzing Ce psp itself, I
wanted to know whether the modification of QE discussed above is safe or
not. 

Answer of Dr. Sclauzero to my question was as follows:

"Briefly, I think that your problem stems exactly from that fix which I
proposed and is now included in the latest PWscf version. That fix will
actually introduce a problem if the atomic wave function which you want to
use for making projections have zero occupations. I thought that this would
never be the case, but probably one should not exclude it a priori. A
further problem arises if you generate such a pseudopotential (Ce with zero
4f occupation) and want to use it for making the LDA+U projections on the
empty states, since the LD1 code usually saves in the PP file only the
wavefunctions with nonzero occupations. First, you should check if you PP
has the 4f wavefunctions. If it does, you can remove the check upf(nt)%oc(n)
> 0.D0 as a quick fix. If not, then you need to regenerate the PP".

Based on his comment, I did some test calculations after modifying our Ce
psp--simply made zero 4f occupation as nonzero, 0.01 (In (I), 430 0.00 -->
430 0.01). However, there is no total E change between the calculation with
Ce psp of zero 4f occupation and that with nonzero 4f occupation (E
difference <0.001/CeO2). Also,  change in Tr[ns(na)] is minute (e.g.
0.9938046 --> 0.9963736)

My questions to QE users are:

- Maybe I need to fill 4f more? 

- Do you have any explanations about why is occupation matrix sum close to
1.0 not 0.0? Again, while there are some non-orthogonal terms, I think that
for CeO2 Tr[ns(na)] for Ce should be close to zero not one. Did I miss
something?

- How about using starting_ns_eigenvalue? By the way, while I tried to use
this, I couldn't understand instruction. For instance, in example below
([II]), if I want to make initial occupations as zero, how can I do? In
starting_ns_eigenvalue(m,ispin,I), I can see that I should have
starting_ns_eigenvalue(m,1,1) and starting_ns_eigenvalue(m,2,1) if U is
applied to spin-up and spin-down for atomic species 1. Then how about m?
m=1,2,3,4,5,6, and 7? If not, what will be desirable values for m?

Thank you for your advice and comments in advance.

P.S. I post this while there is existing post regarding the questions above
because:

1) The title of the existing post is not suitable and will never be easy for
QE users to search. 

2) I am not comfortable with adding threads under the post which contains my
personal emails (not for the forum) which violate my privacy. 

++++++++++++++++++++++++++++++

[I]

[Atom]

Ce

15              # norb: number of orbitals

100 2.00 -643.0 # nlm occ eigen(- means auto-generate)

200 2.00 -77.4

210 6.00 -68.4

300 2.00 -10.0

310 6.00  -6.0

320 10.0  -6.4

400 2.00  -6.4

410 6.00  -6.0

420 10.0  -6.0

500 2.00  -2.4

510 6.00  -2.0

600 0.00  -1.0

610 0.00  -1.0

520 0.00  -1.0

430 0.00  -1.9

+++++++++++++++++++++++++++++++

[II]

U( 1) =  3.0000      U( 2) =  0.0000

alpha( 1) =  0.0000      alpha( 2) =  0.0000

atom    1   Tr[ns(na)]=   2.0000000

atom    1   spin  1

eigenvalues:  0.1428571 0.1428571 0.1428571 0.1428571 0.1428571 0.1428571
0.1428571

eigenvectors

1   1.0000000  0.0000000  0.0000000  0.0000000  0.0000000  0.0000000
0.0000000

2   0.0000000  1.0000000  0.0000000  0.0000000  0.0000000  0.0000000
0.0000000

3   0.0000000  0.0000000  1.0000000  0.0000000  0.0000000  0.0000000
0.0000000

4   0.0000000  0.0000000  0.0000000  1.0000000  0.0000000  0.0000000
0.0000000

5   0.0000000  0.0000000  0.0000000  0.0000000  1.0000000  0.0000000
0.0000000

6   0.0000000  0.0000000  0.0000000  0.0000000  0.0000000  1.0000000
0.0000000

7   0.0000000  0.0000000  0.0000000  0.0000000  0.0000000  0.0000000
1.0000000

occupations

0.143  0.000  0.000  0.000  0.000  0.000  0.000

0.000  0.143  0.000  0.000  0.000  0.000  0.000

0.000  0.000  0.143  0.000  0.000  0.000  0.000

0.000  0.000  0.000  0.143  0.000  0.000  0.000

0.000  0.000  0.000  0.000  0.143  0.000  0.000

0.000  0.000  0.000  0.000  0.000  0.143  0.000

0.000  0.000  0.000  0.000  0.000  0.000  0.143

-------
Chan-Woo Lee, Ph.D.

Department of Chemistry
University of Pennsylvania
231 South 34th Street
Philadelphia, PA 19104-6323 
Phone: 1-215-898-3564 (Office)
Email:  <mailto:leechanw at sas.upenn.edu> leechanw at sas.upenn.edu /
<mailto:cwandtj at gmail.com> cwandtj at gmail.com

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://www.democritos.it/pipermail/pw_forum/attachments/20111214/f3815e8a/attachment-0001.htm 