I wanted to know if it is possible to apply U in two orbitals of the same
atom (for example d and f orbitals of the same atom). From the user guide
and mailing list discussions it seems to be possible. But when I tried to
do such a calculation I encountered the following problem. I am using
Wien2k_08.2 version.

The first scf cycle was OK. In the second cycle it stops at 'orb -up'
after lapw0 saying 'error in Vorb' and in the dayfile it says 'can't open
file 'case.dmatup' whereas in the case.dmatup the density matrices for
both d and f orbitals are written. I don't really understand where is the
problem.

My input files are given below:

case.inorb

1  1  0                     nmod, natorb, ipr
PRATT  1.0                   BROYD/PRATT, mixing
  1 2 2 3                    iatom nlorb, lorb
  1                          nsic 0..AFM, 1..SIC, 2..HFM
   0.37 0.00        U J (Ry)   Note: we recommend to use U_eff =
   0.52 0.00        U J

case.indm
-9.                  Emin cutoff energy
1                   number of atoms for which density matrix is calculated
1  2  2  3          index of 1st atom, number of L's, L1
0 0           r-index, (l,s)index


Thanks for any help.

Is your case.dmatup file void ?
Did you run a complex calculation which needs
case.indmc file?

No the case.dmatup file is not empty. The density matrices are written in
it. I am not doing a complex calculation.

Thanks for any help.

-best regards

Did you execute LDA+U calculation with different
configuration before(e.g. single orbital u)?
Maybe, case.dmat* files are damaged. please try to update them.

For example,

$ runsp -dm -i 1
$ runsp -orb

or

$ rm *.dmat*
$ runsp -orb

My group is considering buying some hardware that would be dedicated
mostly to running Wien2k as efficiently as possible. We are working on
iron-based magnetic systems. Typically, we would like to look at say 128
atoms supercells. Our computer officer suggested buying 5 blades with each
2 quadcore 2.5GHz Intel E5540 CPUs and 24 GB of RAM (gigabit
interconnet...infiniband might be added later, but I think I will stick
to k-point parallelisation only). So that would give us 40 cores with 3GB
of memory per core. Does this sound like a sensible choice for the type
and size of system we want to look at ?

Comments:
1) It sounds reasonable. Get your vendor(s) to allow you to benchmark
on different systems. This will not just give you ideas about
comparable performance, it will tell you which ones you can trust.
2) You will get reasonably good performance for matrix sizes up to
about 18K if you run with 8 nodes per k-point. If you need to go much
larger (e.g. 24K up) you will need infiniband or equivalent.

The CPU is one of the new "Nehalem I7" cpus, so thats what I would also buy.
Compare prices of different vendors! For those new high end systems
there could be quite some difference.

Most likely you can stay without infiniband, since you deal with
metallic systems and will need k-points anyway. One can use
parallelization with mkl (OMP_NUM_THREAD) which works very well for 2
cores (less efficient for 4 cores) and mpi in lapw0 (just 4-8 cores)
in case memory is an issue.

I have a problem about spin-orbit coupling calculation. I followed the
procedure listed below:
1) Do a general scf calculation:  runsp_lapw -eece -it0 -p -i 100 -ec
0.000001 -NI
2) initso
3) runsp_lapw -so -eece -dm -it0 -p -i 100 -ec 0.000001 -NI
But the calculation stoped with following message in step 3:

>   stop error: the required input file hybrid.indmc for the next step could
not be found

Do I need .indmc file? I didn't find any description of this file in
user-guide.
Any advice is appreciated.

Best regards,

The LAPWSO will always generate complex Kohn-Sham vectors even for
systems with inverse symmetry, so the input files for lapw2 and lapwdm
will be always in2c and indmc, respectively.
"c"  is automatically added to  case.in2  by initso, but for "dmat" you
have to add it by yourself.

You have to provide case.indm, case.indmc and case.inorb and edit it
properly and then do initso_lapw and then run the command.
for just spinorbit coupling run runsp_lapw -so  -orb -i 100 -ec 0.0001 -cc
0.0001 -fc 1.0.hope this helps

Thanks to ***and **. Then how to prepare case.indmc? Simply copy
another case.indm and add "c" to its filename? I hope someone give an
description or reference about the content of this file.

Yes. case.indmc is same as case.indm. Just make a copy of it.

case.indm can be edited like these.........
1st line ...no chanfe
2nd line no of atoms for which spin orbit coupling to be added
For third and subsequent lines, index of atom (see your struct file e.g. 1,
2 etc)
last line ,no change (see user guide for specific case to edit this line)

here i have give n an example for my case which contains 5 inequivalent
atoms, and i have added SO for one orbital and that is f(l=3).

-9                                     emin cutoff
5                                       no of
atms........................(see case.indm in src_templates)
1   1   3                            index of 1st atom, ..............
2   1   3
3   1   3
4   1   3
5    1  3
0 0    rindex

just give the same copy in case.indmc (as spin orbit is complex) and edit
your case.inorb accordingly.hope this will help.

when performing LDA+U calculations I encountered some errors.
To do the LDA+U calculation, I did as follows according the UG:
1 set up a new case ,with the exchange-correlation energy selected as LSDA and the separation energy -6 Ry.
2 In the initialization process, select to perform spin-polarized calculation
3 copy the case.indm and case.inorb from $WIEN2K/SRC_templates to the current case directory
4 edit the case.indm as follows:
  -9.0                     Emin cutoff energy
1                       number of atoms for which density matrix is calculated
1  1  3      index of 1st atom, number of L's, L1
0 0           r-index, (l,s)index
(I want to add orbital dependent potential only on l=3 orbit of the first atom)
5 edit the case.inorb as follows:
   1  1  0                     nmod, natorb, ipr
PRATT  1.0                    BROYD/PRATT, mixing
  1 1 3                          iatom nlorb, lorb
  0                              nsic 0..AFM, 1..SIC, 2..HFM
   0.52 0.00        U J (Ry)   Note: we recommend to use U_eff = U-J and J=0
After the above 5 steps being finished, I issued the command in the terminal:
runsp_lapw -orb -cc 0.0001 -NI
No errors occurred in the first 2 cycles of SCF
However, when the 3rd SCF, the following errors occurred:

cycle 3         (2009 04 14 11:11:46 CST)         (38/97 to go)

>   lapw0         (11:11:46) 6.4u 0.2s 0:07.03 95.7% 0+0k 0+2864io 0pf+0w
>   orb -up         (11:11:53) 0.0u 0.0s 0:00.03 100.0% 0+0k 0+32io 0pf+0w
>   orb -dn         (11:11:53) 0.0u 0.0s 0:00.00 0.0% 0+0k 0+32io 0pf+0w
>   lapw1  -up   -orb        (11:11:53) 334.0u 3.7s 5:43.18 98.4% 0+0k 0+2728io 0pf+0w
>   lapw1  -dn   -orb        (11:17:37) 333.5u 4.0s 5:41.38 98.8% 0+0k 0+2712io 0pf+0w
>   lapw2 -up          (11:23:18) 0.0u 0.0s 0:00.11 100.0% 0+0k 0+280io 0pf+0w
error: command   /home/wings/wien2k/lapw2 uplapw2.def   failed

>   stop error
I do not know where is the wrong?
any help will  be appreciated!

check if there was strong oscilation of the 4f element charge, see NTO*** and
QTL***. For the 4f narrow band, you may need to use a small mixing parameter,
e.g 0.02 or even less, or Fermi broadening. In addition, the 4f orbital
symmetry may also be of a concern: sometimes the symmetry needs to be
lowered when working with LDA+U.

I am trying to compute srtuctural properties of SmTe.  For this I have defined the structure using the experimental lattice constant (a=6.6008 Ang) with RMT value 2.96 and 2.8 for Sm and Te.  When I viewed the structure using Xcrysden, all the atoms of Sm and Te are at the correct postion but it doesn't show any bonding.  At the value smaller than the expt. lattice constant (a=6 Ang), I viewed the structure then it show bonding.

So, May I request you to kindly help me to resolve the problem.

You probably mean: Xcrysden does not draw a bond by default.

Xcrysden has some default distances and draws bonds only if they are fulfilled.
You can change the defaults using   Modify: Atomic radii: Chemical connectivity or covalent radius

I am using the latest WIEN2k_08.2 release on a unix platform.

I am attempting to calculate the eigenvalues for my system in the full
Brillouin zone. I wish to do so using LDA+U with S-O coupling.

I am running the execution commands in the following sequence as
recommended so that I may ramp up the U slowly: runsp_lapw save_lapw case
initso_lapw runsp_lapw -so save_lapw case runsp_lapw -so -orb

The question I have relates to a phrase in the UG on p.83 under the ORB
function: "The presence of such an orbital field may lower the symmetry. In
particular the complex version of LAPW1 must be used."

My question is, if the symmetry is lowered by the orbital field will the
code automatically detect this and run lapw1c? In my case.dayfile this is
not occuring for the -orb run:
    cycle 1         (Mon Aug 18 18:15:55 BST 2008)         (40/99 to go)

>   lapw0         (18:15:55) 4.580u 0.420s 0:06.08 82.2%        0+0k 0+0io 182pf+0w
>   orb -up         (18:16:01) 0.000u 0.000s 0:01.14 0.0%        0+0k 0+0io 98pf+0w
>   orb -dn         (18:16:03) 0.010u 0.010s 0:00.08 25.0%        0+0k 0+0io 0pf+0w
>   lapw1 -up (18:16:03) 352.730u 7.760s 6:02.01 99.5% 0+0k 0+0io 302pf+0w
>   lapw1  -dn          (18:22:05) 352.650u 7.660s 6:00.82 99.8%        0+0k 0+0io 0pf+0w
>   lapwso -up -orb (18:28:06) 177.240u 5.500s 3:05.50 98.5% 0+0k 0+0io
> 266pf+0w
>   lapw2 -c -up -so (18:31:11) 115.440u 3.340s 2:01.17 98.0% 0+0k 0+0io
> 239pf+0w
>   lapw2 -c -dn -so (18:33:13) 115.430u 3.500s 1:59.11 99.8% 0+0k 0+0io
> 0pf+0w
>   lapwdm -up -so -c (18:35:12) 5.250u 0.620s 0:07.51 78.1% 0+0k 0+0io
> 158pf+0w
>   lcore -up        (18:35:19) 0.050u 0.090s 0:01.60 8.7%        0+0k 0+0io 95pf+0w
>   lcore -dn        (18:35:21) 0.060u 0.080s 0:00.45 31.1%        0+0k 0+0io 0pf+0w
>   mixer         (18:35:22) 0.360u 0.410s 0:02.45 31.4%        0+0k 0+0io 139pf+0w
:ENERGY convergence:  0 0.0001 0
:CHARGE convergence:  0 0.0000 0

Otherwise how do I know if I need to direct lapw1c to run and how would I
do it?

When adding SO AND orb; the latter is added in lapwso, which is
"complex" anyway.

Still you have the symmetry which is compatible with a) your crystal
symmetry, and b) with the SO interaction (direction of magnetization).

There are cases (mainly rare-earth compounds), where you may want to
reduce the symmetry further (to get a "localized, atomic-like state",
which does not obey crystal symmetry. But this requires that you
"understand" the symmetries and their consequences (degeneracy of states).

thanks for the help. I am considering a system with 4f moments
and so further reductions of symmetry to allow for the formation of
atomic-like orbitals may be valid.

I am a bit unsure of the technical procedure to achieve this in Wien2k. If
I would like to reduce the system symmetry further after the steps
runsp_lapw save_lapw case initso_lapw runsp_lapw -so ie. after the
spin-orbit part has been added and before I start ramping up the U, how
would I go about it?

If I manually lower the symmetry group chosen in the .struct file then I
could just add the atoms into the unit cell manually that are not produced
by the remaining symmetry operations. Is this the logical way to approach
this?

Alternatively should I start the entire series of calculations with the
lower symmetry?

Usually one should keep as much symmetry as possible!

In general it can lead to all kind of unwanted results when one breaks
symmetry more than necessary.

Thus first check if you really have some symmetry, where some f-orbitals
would still be degenerate (like octahedral or tetrahedral or axial symmetry).

The most simple way to break symmetry is usually to introduce a slicght change in
oneof the positions, i.e. if one atom sits in a high symmetry position (like (0,0,0)),
move it to (dx,dy,dz), dx can be small (0.0001); run the initializations (init and initso)
and before the scf, you can change back the coordinates.

I am attempting to calculate intraband contributions to optical
properties of bcc sodium metal (using the OPTIC/JOINT/KRAM programs for
the first time - under WIEN2k_08.2 (Release 21/4/2008)) and ran into a
problem. Everything runs fine when I calculate the interband
contributions (SWITCH 4 in line 4 of case.injoint) but when I use SWITCH
6 for intraband contributions the KRAM program crashes immediately. I
have searched the archives, read the userguide and FAQs, found at least
one person with what appears to be the identical problem
(http://zeus.theochem.tuwien.ac.a ... -April/009089.html)
which was unresolved in the list. I am sure it must be something simple
I am missing but I'm a bit of a beginner and can’t see what.


The process I followed is as follows:

1) executed run_lapw

2) set 200000 kpts (4720 irreducible)

3) edited case.in2 to read:

FERMI (TOT,FOR,QTL,EFG,FERMI)
-9.0 9.0 0.50 0.05 EMIN, NE, ESEPERMIN, ESEPER0
TETRA 101.000 (GAUSS,ROOT,TEMP,TETRA,ALL eval)
0 0 4 0 4 4 6 0 6 4
12.00 GMAX
NOFILE FILE/NOFILE write recprlist

4) executed run_lapw -s lapw1 -e lcore

5) edited case.inop to read:

200000 1 number of k-points, first k-point
-9.0 3.0 Emin, Emax for matrix elements
2 number of choices (columns in *outmat) - 0: MME into case.mme
1 Re xx
3 Re zz
OFF ON/OFF writes MME to unit 4

[Emin, Emax from case.in1]

6) executed x optic

7) edited case.injoint to read:


1 41 41 : LOWER,UPPER and (optional) UPPER-VAL BANDINDEX
0.0000 0.00100 3.0000 : EMIN DE EMAX FOR ENERGYGRID IN ryd
eV : output units eV / ryd / cm-1
6 : SWITCH
2 : NUMBER OF COLUMNS
0.1 0.1 0.3 : BROADENING (FOR DRUDE MODEL - switch 6,7 -
ONLY)

[upper band index from bottom of case.output2]

8) executed x joint

9) edited case.inkram to read:

0.1 Gamma: broadening of interband spectrum
0.0 energy shift (scissors operator)
1 add intraband contributions? yes/no: 1/0
5.8146 plasma frequencies (from joint, opt 6)
0.20 Gammas for Drude terms

[plasma frequency from case.outputjoint]

10) executed x kram, it immediately crashes with this error:

forrtl: severe (64): input conversion error, unit 10, file
/home/alazicki/Nabcc_0GPa/Nabcc_0GPa.joint
Image PC Routine Line Source
kram 0000000000449FBF Unknown Unknown Unknown
kram 000000000044849E Unknown Unknown Unknown
kram 000000000042D0FC Unknown Unknown Unknown
kram 0000000000408C12 Unknown Unknown Unknown
kram 0000000000408873 Unknown Unknown Unknown
kram 0000000000417457 Unknown Unknown Unknown
kram 0000000000402BDD MAIN__ 82 kram.f
kram 00000000004026AA Unknown Unknown Unknown
libc.so.6 000000359151C4BB Unknown Unknown Unknown
kram 00000000004025EA Unknown Unknown Unknown
0.000u 0.001s 0:00.00 0.0% 0+0k 0+0io 0pf+0w
error: command /home/alazicki/WIEN2k/kram kram.def failed


The output files from optic and joint appear fine except case.intra is
empty (?).


Other problems:

The case.sumrules file is always empty for me (when using switch 4 in
case.injoint and successfully running KRAM). Under what conditions
should we get something written there?

Is it just me or does it take about 5 minutes to perform a search of the
archives? It would be nice if we could get a better search engine somehow.

Thank you very much for any advice!

For crystals having cubic symmetry one only needs to calculate the
optical properties along one direction. So, I see no reason to
have two components in case.inop file. If you do so, you should have
2 plasma frequencies in the case.inkram file. This is the reason the program
crashes

Thanks very much for your help!

I have tried changing the number of columns in case.inop and
case.injoint to 1:

___case.inop___
200000 1       number of k-points, first k-point
-9.0 3.0      Emin, Emax for matrix elements
1             number of choices (columns in *outmat) - 0: MME into case.mme
1             Re xx
3             Re zz
OFF           ON/OFF   writes MME to unit 4

[I have also tried deleting line 5]

___case.injoint___
    1  41     41                : LOWER,UPPER and (optional) UPPER-VAL
BANDINDEX
   0.0000    0.00100   3.0000 : EMIN DE EMAX FOR ENERGYGRID IN ryd
eV                            : output units  eV / ryd  / cm-1
     6                        : SWITCH
     1                        : NUMBER OF COLUMNS
   0.1  0.1  0.1              : BROADENING (FOR DRUDE MODEL - switch 6,7 -
ONLY)

[I have also tried leaving the number of columns here equal to 2 on line 5]

___case.inkram___
  0.1    Gamma: broadening of interband spectrum
  0.0    energy shift (scissors operator)
  1      add intraband contributions? yes/no: 1/0
5.8146  plasma frequencies  (from joint, opt 6)
  0.20   Gammas for Drude terms

and I have the same problem when executing x kram.

I have tried 2 columns in case.inop and case.injoint, and listing 2
plasma frequencies in case.inkram this way:

___case.inkram___
  0.1    Gamma: broadening of interband spectrum
  0.0    energy shift (scissors operator)
  1      add intraband contributions? yes/no: 1/0
5.8146 5.8146  plasma frequencies  (from joint, opt 6)
  0.20   Gammas for Drude terms

and I have the same problem when executing x kram.

....you must run x joint with option 6 and get the plasma frequency then
you run x joint again with option 4. After that, you run x kram with
the plasma frequency and intraband contribution to get the optical properties.
best,

I have finished spin polarised calculation for Gd(2*2*2 supercell) with Cd
impurity.But When I started with same structure for unpolarised calculation
with RKMAX=7.5:GMAX=14 and 300kpoints(to begin with),the computation
doesn;t seems to converge.After 200 cycle we noticed that both charge and
energy are oscillating.So please help me in this regard.I have also
calculated
unpolarised calculation for pure Gd.So what may be the problem?

I'm not surprised that a Gd-supercell does NOT converge in an
unpolarized calculation.

You have to consider the underlying physics: Gd metal has 7 4f
electrons, thus in a polarized calculation this gives a rather stable
situation, i.e. "all" 7 spin-up states are occupied, but all 7 spin-dn
states empty.
In an unpolarized case you must put 7 electrons into "3.5 states" (each
state takes now 2 electrons!) and it is rather naturally that this will
cause problems.

The VERY efficient mixer seems to manage this for pure Gd, but in the
supercell there are many more degrees of freedom and thus convergence is
much more difficult. Eventually, switching to TEMP and a rather large
broadening will help to converge.

In any case, you have to ask yourself: what should this calculation be
good for ?? Even above the Neel/Curie-point one has local magnetic
moments (only the long range order has gone).

I have calculated spin polarised calculation in Gd(2*2*2)for Cd
impurity and the results,hyperfine field, magnetic moment,electric
field gradient reasonably matches with available literature.actually
I am interested to find out the electric field gradient of Cd in Gd at room
temperature which is paramagnetic, and thats why I have started unpolarised
calculation,and faced the problem.So how to proceed in this regard?

A paramagnetic state of a 4f compound does NOT mean that one should do a
unpolarized calculation.
Even in the paramagnetic case, one has local moments on each Gd site, but
they are no longer long range ordered, but randomly oriented, thus there
is no net magnetization.

Eventually you can try the effects of magnetic disorder and flip the spins
on two Gd atoms in your supercell (so that the net moment of the cell is again
zero, but you do not have the strict AF state anymore.

PS: I'd expect a rather small influence on the EFG on Cd by that. Other effects
(lattice expansion and phonons) may play a role too.

I am working on the magnetic properties of rare earth–iron compounds. I have a question in structure optimization.
How to set the U and J parameters when using LDA+U method, since U and J can’t be calculated by ab initio way? So I just set U=J=0. If it is reasonable, then is it necessary to runsp_lapw with ‘–orb’?

Best wishes,

> I am working on the magnetic properties of rare earth–iron   
> compounds. I have a question in structure optimization.
> How to set the U and J parameters when using LDA+U method, since U   
> and J can’t be calculated by ab initio way? So I just set U=J=0. If   
> it is reasonable, then is it necessary to runsp_lapw with ‘–orb’?

With U=J=0, you are doing a LDA calculation... That will not produce  
good results for rare earth-iron compounds.

Check the literature to see what other people found to be meaningful  
U-values for rare earth-iron systems, for instance

Torumba et al., PRB 74, 014409 (2006),  
http://link.aps.org/abstract/PRB/v74/e014409
Torumba et al., PRB 77, 155101 (2008),  
http://link.aps.org/abstract/PRB/v77/e085123
Tran et al., PRB 77, 85123 (2008),  
http://link.aps.org/abstract/PRB/v77/e155101

and references therein. If these U-values do not produce good results  
for your compounds, modify them (slightly) and see whether that  
improves the situation for you.

there is a paper about U determination by Madsen and Novak on the wien2k page:
http://www.wien2k.at/reg_user/textbooks/

They suggest that you put J=0 and choose some positive value of U.
Also there was a discussion on the list about U determination some
months ago, search the list.

I am using Wien2k_08.1. I'd like to calculate the optical properties of
GaAs using the OPTIC function in Wien. I have run the calculations
without including the interband transitions and it is fine but when I
include the interband contributions, the programme crashes at the Kram
stage. I think the reason for it is that when the interband
contributions are included the eps1 and eps2 outputs of joint are no
longer valid and hence Kram crashes. This is the procedure that I have
followed:

1) run_lapw
   
  2) x kgen (with 20000 k points)
   
  3) edit case.in2c, changing TOT to FERMI and writing 101.0 instead of
0.000 in TETRA
   
  4) x lapw1
   
  5) x lapw2
   
  6) Modify case.inop:
   
  99999 1       number of k-points, first k-point
   
  0.0 2.2       Emin, Emax for matrix elements
   
  2             number of choices (columns in *outmat) - 0: MME into
case.mme
   
  1             Re xx
   
  3             Re zz
   
  OFF           ON/OFF   writes MME to unit 4
   
  7)  x opticc
   
  8)  Modify case.injoint
   
  1  9999 9999           : LOWER,UPPER and (optional) UPPER-VAL
BANDINDEX
   
   0.0000    0.00100   2.0000 : EMIN DE EMAX FOR ENERGYGRID IN ryd
   
  eV                            : output units  eV / ryd  / cm-1
   
  6                        : SWITCH
   
  2                        : NUMBER OF COLUMNS
   
  0.1  0.1  0.3             : BROADENING (FOR DRUDE MODEL - switch 6,7
-ONLY
   
  9) x joint
   
  10) Modify case.inkram
   
  11) It crashes immediately giving the error message:

    forrtl: severe (64): input conversion error, unit 10,      file
/home/sk701/Wien2k/GaAs_12-02-08/GaAs_12-02-08.joint
Image              PC                Routine            Line
Source            
kram               000000000044BE97  Unknown               Unknown
Unknown
kram               000000000044A25E  Unknown               Unknown
Unknown
kram               000000000042E274  Unknown               Unknown
Unknown
kram               00000000004094E3  Unknown               Unknown
Unknown
kram               000000000040914C  Unknown               Unknown
Unknown
kram               0000000000418218  Unknown               Unknown
Unknown
kram               0000000000402CD0  Unknown               Unknown
Unknown
kram               000000000040272A  Unknown               Unknown
Unknown
libc.so.6          00002AAAAAE4B5AA  Unknown               Unknown
Unknown
kram               000000000040266A  Unknown               Unknown
Unknown
0.000u 0.000s 0:00.00 0.0% 0+0k 0+0io 0pf+0w
error: command   /usr/local/wien2k_new/kram kram.def   failed

I have searched through the mailing list and other people seemed to have
had the same problem but I couldn't find an answer to anyone's problem
posted on the list. I have run the calculations for TiC as well but the
same problem occurs. I would appreciate your suggestions on how to solve
this problem.

You have made already the proper analysis:
You need to run the intraband contributions first (and get the Plasma
frequency out of it, which you should place into case.inkram).
Then calculate the interband part and kram.

PS: Intraband contributions are for metals ! Thus you should not even
attempt them for GaAs.

So don't use switch 6 incase.injoint