visualizing dftb orbital density

[Muhammad Radifar]

Dear Macicris,

So far I've never done any DFTB calculation but I think the way to
produce the file that can show the electron density or molecular
orbital should be the same with any other calculation. Here is the
excerpt of input file to produce the electron density file:
&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../BASIS_SET
    POTENTIAL_FILE_NAME ../POTENTIAL
    &MGRID
      CUTOFF 150
    &END MGRID
    &QS
      EPS_DEFAULT 1.0E-8
    &END QS
    &SCF
      EPS_SCF 1.0E-4
      SCF_GUESS ATOMIC
    &END SCF
    &XC
      &XC_FUNCTIONAL Pade
      &END XC_FUNCTIONAL
    &END XC
    &PRINT
      &E_DENSITY_CUBE
#      STRIDE 1 1 1
      &END E_DENSITY_CUBE
      &TOT_DENSITY_CUBE
#      STRIDE 1 1 1
      &END TOT_DENSITY_CUBE
    &END PRINT
  &END DFT
......

Notice the sections and keywords inside the print section, the
e_density_cube is used to produce the electron density cube, you can
view the cube file using any molecule viewer that support cube file,
my favorite for this is VMD. While the tot_density_cube is used to
produce the total density cube file, the total density is the density
of core + electron. The stride keyword is to adjust the smoothness of
your cube file, the default is 2 2 2, 1 1 1 will produce the smoother
one but it uses more harddisk space.

To produce the molecular orbital here is the excerpt of input file
that I use:
&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../BASIS_SET
    POTENTIAL_FILE_NAME ../POTENTIAL
    &MGRID
      CUTOFF 100
    &END MGRID
    &QS
      EPS_DEFAULT 1.0E-8
    &END QS
    &SCF
      EPS_SCF 1.0E-4
      SCF_GUESS ATOMIC
    &END SCF
    &XC
      &XC_FUNCTIONAL Pade
      &END XC_FUNCTIONAL
    &END XC
    &PRINT
      &MO_CUBES
#      STRIDE 1 1 1
        NHOMO 1
        NLUMO 1
      &END MO_CUBES
    &END PRINT
.....

With the mo_cubes I found that the result is the HOMO and LUMO for
overall system, and I found that using the localized one can show the
HOMO and LUMO for individual molecule:
&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../BASIS_SET
    POTENTIAL_FILE_NAME ../POTENTIAL
    &MGRID
      CUTOFF 100
    &END MGRID
    &QS
      EPS_DEFAULT 1.0E-8
    &END QS
    &SCF
      EPS_SCF 1.0E-4
      SCF_GUESS ATOMIC
    &END SCF
    &XC
      &XC_FUNCTIONAL Pade
      &END XC_FUNCTIONAL
    &END XC
    &LOCALIZE
    &END LOCALIZE
    &PRINT
      &LOCALIZATION
        &WANNIER_CUBES
#         STRIDE 1 1 1
        &END WANNIER_CUBES
      &END LOCALIZATION
    &END PRINT

Watch out that there will be so many orbital CUBE file produced,
usually the HOMO(or LUMO) is the one that we're interested with. So I
think you have to find out which one is the HOMO by try to run it with
run_type energy then open the orbital file one by one. And then you
specify which is the HOMO you're looking for, like this:
    &LOCALIZE
    &END LOCALIZE
    &PRINT
      &LOCALIZATION
        &WANNIER_CUBES
#        STRIDE 1 1 1
          CUBES_LIST 13 14 15 16
        &END WANNIER_CUBES
      &END LOCALIZATION
    &END PRINT

The cubes_list specify which HOMO that you want. If you don't get it,
you might want to run this input file and view the cube file it
produce:
&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../BASIS_SET
    POTENTIAL_FILE_NAME ../POTENTIAL
    CHARGE +1
    &MGRID
      CUTOFF 100
    &END MGRID
    &QS
      EPS_DEFAULT 1.0E-8
    &END QS
    &SCF
      EPS_SCF 1.0E-4
      SCF_GUESS ATOMIC
    &END SCF
    &XC
      &XC_FUNCTIONAL Pade
      &END XC_FUNCTIONAL
    &END XC
    &LOCALIZE
    &END LOCALIZE
    &PRINT
      &LOCALIZATION
        &WANNIER_CUBES
#        STRIDE 1 1 1
          CUBES_LIST 13 14 15 16
        &END WANNIER_CUBES
      &END LOCALIZATION
    &END PRINT
  &END DFT
  &SUBSYS
    &CELL
      ABC 8.0 8.0 8.0
    &END CELL
    &COORD
 O         4.972295        2.976983        4.667921
 H         5.092295        3.036983        3.677921
 H         4.012295        2.726983        4.827922
 O         5.112295        5.365985        5.587922
 H         5.452294        5.585985        6.457922
 H         4.932295        6.245984        5.207923
 O         2.522294        2.386984        5.077923
 H         2.098294        2.915983        5.602922
 H         2.212294        1.496982        5.287923
 O         4.832295        3.230984        2.240921
 H         3.929295        3.330984        2.052921
 H         5.224295        3.257984        1.286921
 H         5.072295        3.906982        5.017922
    &END COORD
    &KIND H
      BASIS_SET DZVP-GTH-PADE
      POTENTIAL GTH-PADE-q1
    &END KIND
    &KIND O
      BASIS_SET DZVP-GTH-PADE
      POTENTIAL GTH-PADE-q6
    &END KIND
  &END SUBSYS
&END FORCE_EVAL
&GLOBAL
  PROJECT H2O
  RUN_TYPE MD
  PRINT_LEVEL LOW
&END GLOBAL
&MOTION
  &MD
    ENSEMBLE NVE
    STEPS 2
    TIMESTEP 0.5
    TEMPERATURE 300.0
  &END MD
&END MOTION

One more thing, there is one caveat in cube file printing, read this:
http://groups.google.com/group/cp2k/browse_thread/thread/f94b4924e091443b

I hope that can help :)

Best regards,

Radif


On Dec 21, 4:10 am, Maricris <mlodr... at gmail.com> wrote:
> Hi,
>
> I would like to know how to view the electron density from DFTB
> results. How can I generate the file that contains the charge
> distribution, molecular orbitals, etc from the DFTB calculation? Is
> this possible?
>
> Any input would be highly appreciated.
>
> Best,
> Maricris