TURBOMOLE Modules > Molecular Properties, Wavefunction Analysis, and Interfaces to Visualization Tools
Smaller grid for calculation of electron density
prasanta13:
Hi there, I am trying to calculate the electron density of H2 molecule. I have ran a simple scf calculation. Control file is provided below. I got the output as td.xyz.
--- Code: ---$title
$symmetry c1
$user-defined bonds file=coord
$coord file=coord
$optimize
internal off
redundant off
cartesian on
global off
basis off
$atoms
h 1-2 \
basis =h def2-QZVPPD
$basis file=basis
$scfmo file=mos
$closed shells
a 1 ( 2 )
$scfiterlimit 30
$scfconv 7
$thize 0.10000000E-04
$thime 5
$scfdamp start=0.300 step=0.050 min=0.100
$scfdump
$scfintunit
unit=30 size=0 file=twoint
$scfdiis
$maxcor 500 MiB per_core
$scforbitalshift automatic=.1
$drvopt
cartesian on
basis off
global off
hessian on
dipole on
nuclear polarizability
$interconversion off
qconv=1.d-7
maxiter=25
$coordinateupdate
dqmax=0.3
interpolate on
statistics 5
$forceupdate
ahlrichs numgeo=0 mingeo=3 maxgeo=4 modus=<g|dq> dynamic fail=0.3
threig=0.005 reseig=0.005 thrbig=3.0 scale=1.00 damping=0.0
$forceinit on
diag=default
$energy file=energy
$grad file=gradient
$forceapprox file=forceapprox
$pointval fmt=xyz
$rundimensions
natoms=2
nbf(CAO)=76
nbf(AO)=66
$last step dscf
$orbital_max_rnorm 0.44543416007246E-05
$last SCF energy change = -1.1335538
$charge from dscf
0.000 (not to be modified here)
$dipole from dscf
x 0.00000000000002 y -0.00000000000017 z 0.00000000000088 a.u.
| dipole | = 0.0000000000 debye
$end
--- End code ---
The td.xyz starts as.
--- Code: ---#origin 0.000000 0.000000 0.000000
#vector1 1.000000 0.000000 0.000000
#vector2 0.000000 1.000000 0.000000
#vector3 0.000000 0.000000 1.000000
#grid1 start -5.313977 delta 0.205059 points 44
#grid2 start -7.744278 delta 0.203075 points 46
#grid3 start 4.056177 delta 0.205458 points 40
#title for this grid 101
#density
#plotdata
# cartesian coordinates x,y,z and f(x,y,z)
-5.31397670 -7.74427833 4.05617653 0.00000006
-5.10891784 -7.74427833 4.05617653 0.00000008
-4.90385899 -7.74427833 4.05617653 0.00000010
-4.69880013 -7.74427833 4.05617653 0.00000013
-4.49374127 -7.74427833 4.05617653 0.00000017
-4.28868241 -7.74427833 4.05617653 0.00000021
-4.08362355 -7.74427833 4.05617653 0.00000026
-3.87856469 -7.74427833 4.05617653 0.00000032
-3.67350583 -7.74427833 4.05617653 0.00000039
--- End code ---
What I get is that delta is the grid spacing, for each direction, x, y and z. Can I obtain a smaller grid spacing? Rather than 0.20 I want to get 0.02 or even smaller?
antti_karttunen:
Hi,
yes, you can set the density of the grid yourself. This is done with sub-keywords grid1, grid2, and grid3 under the $pointval keyword. The option "points" determines the number of grid points in each direction:
--- Code: ---$pointval
grid1 points 200
grid2 points 200
grid3 points 200
--- End code ---
You can also adjust the grid vectors and ranges this way. More information in the manual section 22.2.26 Keywords for wave function analysis and generation of plotting data.
Best,
Antti
uwe:
Hi,
for a smaller grid size and if you do not need a specific number, just add 'quick' to the $pointval line.
The number of points is set to half of the default in each direction. Might be too sparse, but TmoleX uses it when you generate 3D grid data for visualization by clicking on the 'quick' button instead of the 'high-res' one. Usually sufficient for visualization purposes.
Best, Uwe
prasanta13:
Actually, not for plotting purposes but for numerical reasons, I am trying to integrate the total electron density so that it matches the total electron number of the molecule.
For H2 molecule, the total density value should be 2.
Without changing the grid start, and origin I am just increasing the grid points (thus decreasing the delta or grid size) such that I reach the value closest to 2. I am showing the results,
delta density
0.03 1.976
0.025 1.980
0.015 1.9869
Are there any quick hacks to get the electron density very close to the actual one? Nay special keyword to write in control?
uwe:
Hi,
two things you could/might have to do:
1. grid
To catch all of the density you probably need to use a bigger box which reaches further out in space - depending on the basis set.
(a bit off topic):
To find out which atomic orbitals contribute to the density (by checking the occupied molecular orbitals and printing the contribution of each atomic orbital) try:
tm2molden mostat <list-of-MOs> above <threshold>
So if you have 40 occupied orbitals (eiger quickly tells you how many you have) and you want to see all AOs which contribute to each MO more than, say, 0.01% just call on the command line:
tm2molden mostat 1-40 above 0.01
Enter a few times <enter> and then you will get a file named mostat.out which looks like that:
MO Index (from eiger): 3
Sym= 3a
Energy [Hartree]= -.99177169241775E+01
Energy [eV]= -.26987482366088E+03
Spin= Alpha
Occup= 2.000000
contribution of AO 2 c 1 s is: 24.939 %
contribution of AO 6 c 1 s is: 24.939 %
contribution of AO 3 c 1 s is: 24.939 %
contribution of AO 5 c 1 s is: 24.939 %
contribution of AO 2 c 2 s is: 0.052 %
and so on. That is useful to see how many diffuse AOs contribute to the occupied orbitals and thus to the density.
2. numerical accuracy
Another source of numerical inaccuracy is the format of either plt or cube files. They are by default single-precision data. As the density is written to a grid and then summed up, just from pure numerics the sum will never exactly be the number of electrons (well, never is not correct, it could coincidentally be just an integer number...).
Try to add fmt=dtx to the $pointval input line instead. This will print double precision data as ASCII text (dtx = double precision text file). Looks like that:
#origin 0.000000 0.000000 0.000000
#vector1 1.000000 0.000000 0.000000
#vector2 0.000000 1.000000 0.000000
#vector3 0.000000 0.000000 1.000000
#grid1 start -5.767470 delta 0.202135 points 66
#grid2 start -5.767470 delta 0.202829 points 59
#grid3 start -5.020450 delta 0.204532 points 67
#title for this grid 101
#density
#plotdata
# cartesian coordinates x,y,z and f(x,y,z)
-5.7674700000000E+00 -5.7674700000000E+00 -5.0204500000000E+00 3.9793580979130E-10
-5.5653350769231E+00 -5.7674700000000E+00 -5.0204500000000E+00 5.8733532389287E-10
-5.3632001538462E+00 -5.7674700000000E+00 -5.0204500000000E+00 8.5021125663768E-10
-5.1610652307692E+00 -5.7674700000000E+00 -5.0204500000000E+00 1.2071193363237E-09
...
The last column then needs to be summed up. Whether or not that will increase accuracy - it's at least worth a try.
Navigation
[0] Message Index
[#] Next page
Go to full version