The ADF package contains a series of sample runs. Provided are UNIX scripts to run the calculations and the resulting output files. In most directories, there are also files for ADFinput present.
The examples serve:
Where references are made to the operating system (OS) and to the file system on your computer,
the terminology of a UNIX type OS is used and a hierarchical structure of directories is assumed.
All sample files are stored in subdirectories under $ADFHOME/examples/, where $ADFHOME is the main directory of the ADF package. There are two main subdirectories in examples/: adf/ for calculations with the molecular code ADF (and related utility programs) and band/ for calculations with the periodic structures code BAND. Each sample run has its own directory (under adf/ or band/ respectively). For instance, $ADFHOME/examples/adf/HCN/ contains an ADF calculation on the HCN molecule. Each sample subdirectory contains:
Notes:
Many of the provided samples have been devised to be short and simple, at the expense of physical or chemical relevance and precision or general quality of results. They serve primarily to illustrate the use of input, necessary files, and type of results. The descriptions have been kept brief. Extensive information about using keywords in input and their implications is given in the User's Manuals (ADF and BAND).
When you compare your own results with the sample outputs, you should check in particular (as far as applicable):
General remarks about comparisons:
Default settings of print options result in a considerable amount of output. This is also the case in some of the sample runs, although in many of them quite a bit of 'standard' output is suppressed by inserting applicable print control keys in the input file. Consult the User's Guide about how to regulate input with keys in the input file.
The Survey of Applications follows a survey of the main application topics with references to related sample runs is given. A sample run usually involves several calculations, for instance a few CREATE runs (with ADF), then a molecular calculation (also ADF), and finally a NMR calculation (with the NMR program) to compute chemical shifts. The samples are identified in this documentation by the name of the directory they reside in. The samples are indicated by these directory names. For instance, GO_H2O refers to the directory GO_H2O/ (in $ADFHOME/adf/), where in this case GO stands for Geometry Optimization.
For property calculations, xc potentials with asymptotically correct (-1/r) behavior outside the molecule, the results tend to be superior to regular LDA or GGA calculations. This is especially true for small molecules and for properties that depend heavily on the proper description of the outer region of the molecule. In the example, all-electron basis sets are used. This is mandatory for the SAOP potential.
In the first example, excitation energies are calculated with the GRACLB potential. This potential requires one number as argument: the experimental ionization potential in atomic units. This number can be either based on an experimental value, or on previous GGA total energy calculations.
$ADFBIN/adf <<EOR title CO excitations grac potential INTEGRATION 6.0 XC Model GRACLB 0.515 End Atoms O 0 0 0 C 1.128205364 0 0 end Excitation Lowest 10 Onlysing End Basis Type TZ2P Core None End end input EOR rm TAPE21 logfile
The same calculation with the SAOP xc potential would differ in the XC block only:
XC Model SAOP End
SAOP depends on the orbitals which makes it more expensive to evaluate than GRAC for large molecules, but is much easier to use, since it does not require an ionization potentiali parameter as input.
Sample directory adf/OH_MetaGGA
First two calculations on OH are performed which use, respectiveley, the hybrid meta-GGA TPSSh and the meta-GGA TPPS during the SCF. They require, respectively, the following XC input:
XC MetaHybrid TPSSh END
XC MetaGGA TPSS END
Next large even-tempered basis sets are used in the calculation of the atomization energy of OH using various modern GGA, meta-GGA and hybrid post-SCF energy expressions.
In the Create runs, a large even-tempered basis set is selected for O and H, which should give results closer to the basis set limit than the regular ADF basis sets. For both atoms, a second atomic calculation follows the Create run, in order to enable a comparison to the true atoms, rather than the artificial spherically symmetric atom from the Create run. This is achieved by specifying the keywords
unrestricted charge 0 2 symmetry C(lin) occupations sigma 3 // 3 pi 2 // 0 end
in the case of oxygen. This fixes the proper occupations. The result files of both the Create runs and the atomic correction runs are stored.
In the molecular calculation, the symmetry of the molecule is explicitly broken and the occupations are specified in order to avoid the fractional occupations that ADF would otherwise choose. Although it is not said that such a solution would be inferior, the integer occupation solution is the one which allows direct comparison to literature results obtained with other programs.
One of the new GGA potentials has been specified for the xc potential and the keyword METAGGA implies that a series of GGA and meta-GGA xc energies is to be calculated and compared to those energies from the atomic calculations. Specifying HARTREEFOCK also enables calculation of PostSCF energies using hybrid functionals.
METAGGA symmetry C(lin) xc GGA PBE end HARTREEFOCK
A fairly high numerical integration has been specified. For meta-GGA calculations we do recommend this, at least 6 for the time being, as the numerical stability of the results tends to be somewhat lower than for regular GGA calculations.
The block key ENERGYFRAG
ENERGYFRAG O t21.unr.O H t21.unr.H END
implies that the meta-GGA result must not only be compared to the spherically symmetric results from the Create runs, but also to the non-spherical atoms.
The molecular output file prints the PBE Total Bonding energy as usual (in various energy units).
Then a prints a list of 'Total Bonding Energies' for many different Exc functionals, including PBE. Because the numerical approach to obtain the two PBE results is somewhat different, small differences may occur between the two numbers. You now have an overview of the bonding energies of all (meta)GGA functionals currently implemented in ADF. This should give a good indication of the theoretical error bar or the uncertainty in the xc approximation.
Total Bonding Energy: -0.286127457276205 -7.7859 -179.55 -751.23
TOTAL BONDING ENERGIES FROM VARIOUS XC FUNCTIONALS
with respect to fragments in FRAGMENTS input block
hartree eV kcal/mol kJ/mol
Total Bonding Energy with respect to FRAGMENTS
XC Energy Functional
====================
FR: KCIS-modified [1] = -0.2755742057 -7.4987587362 -172.9254430523 -723.5200549587
FR: KCIS-original [2] = -0.2777894828 -7.5590395194 -174.3155506035 -729.3362649626
FR: PKZB [3] = -0.2815570432 -7.6615600946 -176.6797306630 -739.2279943483
FR: VS98 [4] = -0.3017049511 -8.2098127875 -189.3227350810 -792.1263249228
FR: LDA(VWN) [5] = -0.2887564297 -7.8574654492 -181.1974143810 -758.1299830563
FR: PW91 [6] = -0.2876922977 -7.8285089331 -180.5296614163 -755.3361046473
FR: BLYP [7] = -0.2770745036 -7.5395839361 -173.8668943006 -727.4590869882
FR: BP [8] = -0.2855241909 -7.7695117221 -179.1691537365 -749.6437405057
FR: PBE [9] = -0.2858734106 -7.7790144775 -179.3882924288 -750.5606167957
.....
The same energy comparison is done with respect to the fragments (which most currently be atomic) in the ENERGYFRAG block. These are the numbers which should be comparable to experimental numbers.
Finally, the references for the various Exc functionals are printed in the output file.
XC Energy Functional ==================== EF: KCIS-modified [1] = -0.1713622482 -4.6630059333 -107.5314455812 -449.9115690750 EF: KCIS-original [2] = -0.1701706820 -4.6305817515 -106.7837263654 -446.7831118709 EF: PKZB [3] = -0.1716508948 -4.6708604097 -107.7125740668 -450.6694106602 EF: VS98 [4] = -0.1712676117 -4.6604307410 -107.4720602503 -449.6631008503 EF: LDA(VWN) [5] = -0.1980694328 -5.3897456994 -124.2904587006 -520.0312800855 EF: PW91 [6] = -0.1759694023 -4.7883730257 -110.4224787188 -462.0076517434 EF: BLYP [7] = -0.1748768123 -4.7586421272 -109.7368680765 -459.1390568111 EF: BP [8] = -0.1785853781 -4.8595573769 -112.0640284617 -468.8758958794 EF: PBE [9] = -0.1751227104 -4.7653333576 -109.8911714787 -459.7846622469 ....
Similar calculations can be done to obtain energy differences between different molecules. In that case the ENERGYFRAG keyword is not operational though. No detailed breakdown of the bonding energy is currently available for these new energy functionals. Experience shows that the energy values depend only mildly on the chosen xc functional for the xc potential.
Example shows a Hartree-Fock calculation with a non-relativistic, scalar relativistic ZORA, and a spin-orbit coupled ZORA Hamiltonian. In this case ADF also calculates the electric field gradient (EFG) at the H and I nuclei (keyword QTENS).
First the non-relativistic calculation. Note that in this case the all-electron basis sets are obtained from the $ADFRESOURCES/ZORA directory.
$ADFBIN/adf << eor Atoms H 0 0 0 I 0 0 1.609 End qtens xc hartreefock end integration 5 Basis Type ZORA/TZ2P Core None End End input eor
Next the scalar relativistic ZORA calculation. Note that in this case the all-electron basis sets are also obtained from the $ADFRESOURCES/ZORA directory, but this is default place where the key BASIS will search for basis sets in case of ZORA. ADF will also calculate the EFG including the small component density, also called SR ZORA-4.
$ADFBIN/adf << eor Atoms H 0 0 0 I 0 0 1.609 End qtens xc hartreefock end Relativistic Scalar ZORA integration 5 Basis Type TZ2P Core None End End input eor
Next the spin-orbit coupled relativistic ZORA calculation. Note that in this case the all-electron basis sets are also obtained from the $ADFRESOURCES/ZORA directory, but again this is default place where the key BASIS will search for basis sets in case of ZORA. If one calculates this molecule with symmetry nosym, ADF will also calculate the EFG including the small component density, also called ZORA-4.
$ADFBIN/adf << eor Atoms H 0 0 0 I 0 0 1.609 End qtens xc hartreefock end Relativistic Spinorbit ZORA symmetry nosym integration 5 Basis Type TZ2P Core None End End input eor
Sample directory: adf/H2PO_B3LYP/
Example shows an unrestricted B3LYP calculation. In this case ADF also calculates the hyperfine interactions at H, P, and O nuclei (keyword ESR).
The 'DEPENDENCY' key is set to 1e-4. Note that for hybrids and Hartree-Fock the dependency key is always set. The default value in that case is 4e-3. By explicitely setting the 'DEPENDENCY' key we can use a lower value, which is possible in this case. One should check that the results remain reliable if one uses a smaller value for the 'DEPENDENCY' key.
$ADFBIN/adf << eor
Title hfs H2PO B3LYP TZ2P
Atoms
O 1.492 0.000 0.000
P 0.000 0.000 0.000
H -0.600 -0.650 1.100
H -0.600 -0.650 -1.100
End
xc
hybrid B3LYP
end
Basis
Type TZ2P
Core None
End
dependency bas=1e-4
integration 5
esr
end
unrestricted
charge 0 1
end input
eor
For the hyperfine interactions it is important to use all-electron basis sets on the interesting nuclei. One can get more accurate results if one uses a larger basis set, like the QZ4P basis set, which is present in the $ADFRESOURCES/ZORA directory. The Basis key should then be:
Basis Type ZORA/QZ4P Core None End
The QZ4P results for the isotropic value of the A-tensor are approximately: -24.77 MHz for 17O, 962.02 MHz for 31P, and 110.72 MHz for 1H.
You may want to compare the results with previous B3LYP results by N. R. Brinkmann and I. Carmichael, J. Phys. Chem. A (2004), 108, 9390-9399, which give for the Isotropic Fermi Contact Couplings (MHz) for the 2A' State of H2PO using B3LYP, with an aug-cc-pCVQZ basis set: -24.24 MHz for 17O, 963.33 MHz for 31P, and 111.51 MHz for 1H.
Sample directory: adf/MM_Dispersion/
Summary:
First example shows a geometry optimization of a van der Waals complex of two benzene molecules, connected to each other with a hydrogen molecule. With the MMDISPERSION keyword an extra empirical force (of similar form as in molecular mechanics) is added to the interaction between the three fragments, where one benzen molecule is fragment 1 (FD=1), the other benzene molecule is fragment 2 (FD=2), and the hydrogen molecule is fragment 3 (FD=3).
The atomic parameters are read from the file $ADFRESOURCES/MMDispersion/disp-param. The PBE functional and the TZP basis set are used, which is necessary if one wants to use the TZ parameters for the damping function, which are optimized for this combination of functional and basis set.
$ADFBIN/adf << eor basis type TZP core small End XC GGA PBE End geometry converge grad=0.001 iterations 5 end Integration 4.5 SCF Iterations 60 Converge 1.0E-06 1.0E-6 End mmdispersion damping sigm damp_param tz combi s-k file_name $ADFRESOURCES/MMDispersion/disp-param nodefault end noprint sfo Atoms cartesians C.ctr 0.000000000000 3.050000000000 1.391500000000 FD=1 H.h 0.000000000000 3.050000000000 2.471500000000 FD=1 C.ctr 1.205074349366 3.050000000000 0.695750000000 FD=1 H.h 2.140381785453 3.050000000000 1.235750000000 FD=1 C.ctr 1.205074349366 3.050000000000 -0.695750000000 FD=1 H.h 2.140381785453 3.050000000000 -1.235750000000 FD=1 C.ctr -0.000000000000 3.050000000000 -1.391500000000 FD=1 H.h -0.000000000000 3.050000000000 -2.471500000000 FD=1 C.ctr -1.205074349366 3.050000000000 -0.695750000000 FD=1 H.h -2.140381785453 3.050000000000 -1.235750000000 FD=1 C.ctr -1.205074349366 3.050000000000 0.695750000000 FD=1 H.h -2.140381785453 3.050000000000 1.235750000000 FD=1 C.ctr -1.205074349366 -3.050000000000 -0.695750000000 FD=2 H.h -2.140381785453 -3.050000000000 -1.235750000000 FD=2 C.ctr -0.000000000000 -3.050000000000 -1.391500000000 FD=2 H.h -0.000000000000 -3.050000000000 -2.471500000000 FD=2 C.ctr 1.205074349366 -3.050000000000 -0.695750000000 FD=2 H.h 2.140381785453 -3.050000000000 -1.235750000000 FD=2 C.ctr 1.205074349366 -3.050000000000 0.695750000000 FD=2 H.h 2.140381785453 -3.050000000000 1.235750000000 FD=2 C.ctr -0.000000000000 -3.050000000000 1.391500000000 FD=2 H.h -0.000000000000 -3.050000000000 2.471500000000 FD=2 C.ctr -1.205074349366 -3.050000000000 0.695750000000 FD=2 H.h -2.140381785453 -3.050000000000 1.235750000000 FD=2 H.h 0.0 0.35 0.0 FD=3 H.h 0.0 -0.35 0.0 FD=3 End End Input
The part of the bond energy that is due to the Grimme dispersion corrected functional is only inter-molecular (atom-atom contributions for which the fragment numbers FD are different).
In the second example a structure with 2 benzene molecules and a hydrogen molecule is optimized with the Grimme dispersion corrected PBE. Needed is the subkey DISPERSION in the key XC. If one starts with atomic fragments the part of the bond energy that is due to the Grimme dispersion corrected functional is both inter-molecular as well as intra-molecular. In this case the subargument FD= in the ATOMS block key word is not used, which was only used in the old MM dispersion calculation.
$ADFBIN/adf << eor Title Geometry optimization with Grimme dispersion correction for GGA basis type TZP core small End XC GGA PBE DISPERSION End geometry converge grad=0.001 Branch OLD iterations 50 end Integration 4.5 Atoms cartesians C 0.000000000000 3.050000000000 1.391500000000 H 0.000000000000 3.050000000000 2.471500000000 C 1.205074349366 3.050000000000 0.695750000000 H 2.140381785453 3.050000000000 1.235750000000 C 1.205074349366 3.050000000000 -0.695750000000 H 2.140381785453 3.050000000000 -1.235750000000 C -0.000000000000 3.050000000000 -1.391500000000 H -0.000000000000 3.050000000000 -2.471500000000 C -1.205074349366 3.050000000000 -0.695750000000 H -2.140381785453 3.050000000000 -1.235750000000 C -1.205074349366 3.050000000000 0.695750000000 H -2.140381785453 3.050000000000 1.235750000000 C -1.205074349366 -3.050000000000 -0.695750000000 H -2.140381785453 -3.050000000000 -1.235750000000 C -0.000000000000 -3.050000000000 -1.391500000000 H -0.000000000000 -3.050000000000 -2.471500000000 C 1.205074349366 -3.050000000000 -0.695750000000 H 2.140381785453 -3.050000000000 -1.235750000000 C 1.205074349366 -3.050000000000 0.695750000000 H 2.140381785453 -3.050000000000 1.235750000000 C -0.000000000000 -3.050000000000 1.391500000000 H -0.000000000000 -3.050000000000 2.471500000000 C -1.205074349366 -3.050000000000 0.695750000000 H -2.140381785453 -3.050000000000 1.235750000000 H 0.0 0.35 0.0 H 0.0 -0.35 0.0 End End Input
In the last example first three molecules (2 benzene molecules and a hydrogen molecule) are calculated with the Grimme dispersion corrected PBE. Needed again is the subkey DISPERSION in the key XC. The one for H2 is given below:
$ADFBIN/adf << eor Title Grimme dispersion-corrected GGA basis type TZP core small End XC GGA PBE DISPERSION End SCF Iterations 60 Converge 1.0E-06 1.0E-6 End Atoms H 0.000000 0.000000 -0.377906 H 0.000000 0.000000 0.377906 End End Input eor mv TAPE21 h2.t21
Note that even for such a molecule there is a contribution from the so called Dispersion energy in the bonding energy (although it will be very small in this case).
Next a structure is calculated in which the three calculated molecules in it. If one starts with molecular fragments the part of the bond energy that is due to the Grimme dispersion corrected functional is only inter-molecular.
$ADFBIN/adf << eor Title Grimme dispersion-corrected GGA Fragments b1 benzene1.t21 b2 benzene2.t21 h2 h2.t21 End XC GGA PBE DISPERSION End Atoms C 0.000000 1.398973 -3.054539 f=b1 H 0.000000 2.490908 -3.049828 f=b1 C 1.211546 0.699486 -3.054539 f=b1 H 2.157190 1.245454 -3.049828 f=b1 C 1.211546 -0.699486 -3.054539 f=b1 H 2.157190 -1.245454 -3.049828 f=b1 C 0.000000 -1.398973 -3.054539 f=b1 H 0.000000 -2.490908 -3.049828 f=b1 C -1.211546 -0.699486 -3.054539 f=b1 H -2.157190 -1.245454 -3.049828 f=b1 C -1.211546 0.699486 -3.054539 f=b1 H -2.157190 1.245454 -3.049828 f=b1 C -1.211546 -0.699486 3.054539 f=b2 H -2.157190 -1.245454 3.049828 f=b2 C 0.000000 -1.398973 3.054539 f=b2 H 0.000000 -2.490908 3.049828 f=b2 C 1.211546 -0.699486 3.054539 f=b2 H 2.157190 -1.245454 3.049828 f=b2 C 1.211546 0.699486 3.054539 f=b2 H 2.157190 1.245454 3.049828 f=b2 C 0.000000 1.398973 3.054539 f=b2 H 0.000000 2.490908 3.049828 f=b2 C -1.211546 0.699486 3.054539 f=b2 H -2.157190 1.245454 3.049828 f=b2 H 0.000000 0.000000 -0.377906 f=h2 H 0.000000 0.000000 0.377906 f=h2 End End Input eor
Sample directory: adf/methane_dimer_dispersion/
The density-dependent dispersion energy correction, dDsC, by S.N. Steinmann and C. Corminboeuf, is used to calculate interactions between nonoverlapping densities, which standard density functional approximations cannot accurately describe. The example is for the methane dimer.
The 'debug dispersion' is included such taht in the output one can find more details on the exact parameters that are used.
$ADFBIN/adf << eor TITLE methane-dimer non-relativistic debug dispersion ATOMS C -0.000959 0.000775 1.853082 H -0.747186 0.712608 1.489389 H 0.987865 0.294742 1.490258 H -0.241511 -0.998876 1.480724 H -0.002970 -0.005330 2.946903 C 0.000962 -0.000776 -1.853082 H 0.004264 0.004684 -2.946903 H -0.989749 -0.289559 -1.491241 H 0.743184 -0.716138 -1.488107 H 0.246099 0.997870 -1.481022 End Basis Type TZP End addDiffusefit INTEGRATION 6.0 XC GGA Becke LYP Dispersion DDsC End Geometry End eprint geostep Energy GradientTerms Gradients Upd end symmetry nosym End input eor
Sample directory: adf/Au2_ZORA/
A relativistic geometry optimization with the ZORA formalism. Both a scalar relativistic geometry optimization, as well as a spin-orbit coupled relativistic geometry optimization are performed. Spin-orbit effects on the geometry are often not so large for closed shell molecules, but takes typically a factor 4 times as much as a scalar relativistic calculation.
The input for the spin-orbit coupled relativistic geometry optimization is given below. The only diffeerence bewteen this inpu and the input for the scalar relativistic geometry optimization is the key RELATIVISTIC. Use 'Relativistic scalar ZORA' in case of sclar relativistic ZORA calculations.
$ADFBIN/adf << eor Title Au2 relativistic optimization: spinorbit ZORA Integration 5 5 Atoms Zmat Au 0 0 0 Au 1 0 0 2.5 End Basis Type TZ2P Core Small End Relativistic SpinOrbit ZORA Geometry convergence grad=1e-4 End End Input eor
Sample directory: adf/SO_ZORA_Bi2/
Application of the Spin-Orbit relativistic option (using double-group symmetry) to Bismuth (atom and dimer).
For comparision with the full double-group calculation, the 'standard' unrestricted calculation on Bismuth is carried out, using the scalar relativistic option.
A net spin polarization of 3 electrons is applied (key charge).
$ADFBIN/adf <<eor title Bi unrestricted integration 6.0 relativistic scalar ZORA ATOMS Bi 0.000000 0.000000 0.00000000 end Basis Type TZ2P Core None end unrestricted charge 0 3 xc GGA becke perdew end end input eor
The CHARGE key, in conjunction with the UNRESTRICTED key is used to specify that 3 electrons must be unpaired (second value of the CHARGE key), while the system is neutral (first value of the CHARGE key).
Next we do a Spin-Orbit calculation on the Bismuth atom.
Note that it is a 'unrestricted' run using the noncillinear approximation, and SYMMETRY NSOYM. The electronic charge density is spin-polarized.
$ADFBIN/adf <<eor title Bi spinorbit unrestricted noncollinear integration 6.0 relativistic spinorbit ZORA ATOMS Bi 0.000000 0.000000 0.00000000 end symmetry nosym unrestricted noncollinear Basis Type TZ2P Core None end xc GGA becke perdew end end input eor
Because an all electron basis set is used, the bond energy is huge, due to the very large higher order spin-orbit effect on the 2p orbitals.
Now we turn to the dimer Bi2: a series of Single Point calculations, all with the same inter atomic distance.
First the scalar relativistic run.
$ADFBIN/adf <<eor title Bi2, scalar relativistic integration 6.0 relativistic scalar ZORA ATOMS Bi 0.0 0.0 1.33 Bi 0.0 0.0 -1.33 end Basis Type TZ2P Core None end xc GGA becke perdew end end input eor mv tape21 t21Bi2
The calculated scalar relativistic atomization energy will be close to 2.74 eV. This is the bond energy of the dimer minus 2 times the bond energy of the unrestricted scalar relativistic atom.
The result file tape21 is used as reference in subsequent calculations: run the spin-orbit case starting from the just completed dimer calculation as a fragment. The resulting 'bonding energy', ie the energy w.r.t. the scalar relativistic ZORA dimer, gives directly the effect of the full-relativistic ZORA versus the scalar relativistic ZORA option: the energy is lowered by huge amount, again mainly due to the large spin-orbit effect on the 2p orbitals.
$ADFBIN/adf <<eor title Bi2 from fragment Bi2, with SpinOrbit coupling PRINT SpinOrbit integration 6.0 relativistic spinorbit ZORA ATOMS Bi 0.0 0.0 1.33 f=Bi2 Bi 0.0 0.0 -1.33 f=Bi2 end fragments Bi2 t21Bi2 end xc GGA becke perdew end end input eor rm TAPE21 logfile
If one looks at the SFO analysis in the output of this calculation, one can observe that a first-order spin-orbit splitting of the scalar relativistic fragment orbitals is a good approximation to many of the calculated valence spinors.
A final consistency check: run the spin-orbit dimer from single-atom fragments. The bonding energy should equal the sum of the bonding energies of the previous two runs: scalar relativistic dimer w.r.t. single atom fragments plus spin-orbit dimer w.r.t. the scalar relativistic dimer.
$ADFBIN/adf <<eor title Bi2 from atomic fragments, SpinOrbit coupling PRINT SpinOrbit integration 6.0 relativistic spinorbit ZORA ATOMS Bi 0.0 0.0 1.33 Bi 0.0 0.0 -1.33 end Basis Type TZ2P Core None end xc GGA becke perdew end end input eor
The calculated spin-orbit coupled relativistic atomization energy will be close to 2.18 eV. This is the bond energy of the dimer minus 2 times the bond energy of the unrestricted non-collinear spin-orbit coupled relativistic atom. Note that one has to subtract huge numbers.
Sample directory: adf/Tl_noncollinear/
Application of the Spin-Orbit relativistic option (using double-group symmetry, in this case NOSYM) to Tl using the collinear and non-collinear approximation for unrestricted Spin-Orbit calculations
Note: For the collinear and the non-collinear approximation one should use symmetry NOSYM and use the key UNRESTRICTED.
The non-collinear example:
$ADFBIN/adf << eor Title Tl spinorbit noncollinear Atoms Tl 0 0 0 End Relativistic Spinorbit ZORA XC gradients becke perdew end symmetry nosym unrestricted noncollinear Fragments Tl t21.Tl End End input eor
If one replaces the key NONCOLLINEAR with COLLINEAR the collinear approximation will be used instead of the non-collinear approximation. In the case of the collinear approximation default the direction of the magnetization is in the direction of the z-axis. In the non-collinear approximation the magnetization can differ in each point in space.
Sample directory: adf/AuH_analyse_exciso/
Calculation of the excitation energies of AuH including spin-orbit coupling.
$ADFBIN/adf << eor Title [AuH] Atoms Au .0000 .0000 1.5238 H .0000 .0000 0.0000 End relativistic scalar zora Basis Type TZ2P Core None End symmetry C(7v) EPRINT SFO eig ovl END integration 6.0 Excitations lowest 40 End End input eor mv TAPE21 t21.frag rm logfile $ADFBIN/adf << eor Title [AuH] Atoms Au .0000 .0000 1.5238 f=Frag H .0000 .0000 0.0000 f=Frag End relativistic spinorbit zora symmetry C(7v) EPRINT SFO eig ovl END integration 6.0 Excitations lowest 40 End Fragments Frag t21.frag End STCONTRIB End input eor
ADF can not handle ATOM and linear symmetries in excitation calculations. Therefore a subsymmetry is used, in this case symmetry C(7v).
A relatively small TZ2P basis set is used, which is not sufficient for excitations to Rydberg-like orbitals, one needs more diffuse functions.
The key STCONTRIB is used, which will give a composition of the spin-orbit coupled excitation in terms of singlet-singlet and singlet-triplet scalar relativistic excitations. In order to use the key STCONTRIB the scalar relativistic fragment should be the complete molecule.
Starting from ADF2008.01 one needs to include the subkey SFO of the key EPRINT with arguments eig and ovl in order to get the SFO MO coefficients and SFO overlap matrix printed on standard output.
Sample directory: adf/Solv_HCl/
Computing solvent effects, with the COSMO model, is illustrated in the HCl example.
After a non-solvent (reference) calculation, which is omitted here, two solvent runs are presented, with somewhat different settings for a few input parameters. The block key Solvation controls all solvent-related input.
All subkeys in the SOLVATION block are discussed in the User's Guide. Most of them are rather technical and should not severely affect the outcome. Physically relevant is the specification of the solute properties, by the SOLVENT subkey: the dielectric constant and the effective radius of the solvent molecule.
Note that a non-electrostatic terms as a function of surface area is included in the COSMO calculation, by setting the values for CAV0 and CAV1 in the subkey SOLVENT of the key SOLVATION. In ADF2010 one should explicitely include such values for CAV0 and CAV1, otherwise this non-electrostatic term will be taken to be zero, since the defaults have changed in ADF2010.
A rather strong impact on the computation times has the method of treating the 'C-matrix'. There are 3 options (see the User's Guide): EXACT is the most expensive, but presumably most accurate. POTENTIAL is the cheapest alternative and is usually quite adequate. EXACT uses the exact charge density for the Coulomb interaction between the molecular charge distribution and the point charges (on the Van der Waals type molecular surface) which model the effects of the solvent. The alternatives, notably 'POTENTIAL', use the fitted charge density instead. Assuming that the fit is a fairly accurate approximation to the exact charge density, the difference in outcome should be marginal.
$ADFBIN/adf << eor TITLE HCl(1) Solv-excl surfac; Gauss-Seidel (old std options) SYMMETRY NOSYM ATOMS Cartesian H 0.000000 0.000000 0.000000 R=1.18 Cl 1.304188 0.000000 0.000000 R=1.75 END Fragments H t21.H Cl t21.Cl End SOLVATION Solvent epsilon=78.8 radius=1.4 cav0=1.321 cav1=0.0067639 SurfaceType esurf DivisionLevel ND=4 min=0.5 Ofac=0.8 ChargeUpdate Method=Gauss-Seidel DiscAttributes SCale=0.01 LEGendre=10 TOLerance=1.0d-2 SCF Variational C-Matrix Exact END NOPRINT Bas EigSFO EKin SFO, frag, functions EPRINT SCF NoEigvec END END INPUT eor rm TAPE21 logfile
In the second solvent run, another (technical) method is used for determining the charge distribution on the cavity surface (conjugate-gradient versus Gauss-Seidel in the previous calculation), and the POTENTIAL variety is used for the C-matrix handling. The results show that it makes little difference in outcome, but quite a bit in computation times.
$ADFBIN/adf << eor TITLE HCl(9) NoDisk and Cmatrix potential FRAGMENTS H t21.H Cl t21.Cl END ATOMS Cartesian H 0.000000 0.000000 0.000000 R=1.18 Cl 1.304188 0.000000 0.000000 R=1.75 END SOLVATION Solvent epsilon=78.8 radius=1.4 cav0=1.321 cav1=0.0067639 SurfaceType esurf DivisionLevel ND=4 min=0.5 Ofac=0.8 ChargeUpdate Method=conjugate-gradient SCF Variational C-Matrix POTENTIAL END NOPRINT Bas EigSFO EKin SFO, frag, functions EPRINT SCF NoEigvec END END INPUT eor
Sample directory: adf/3DRISM-Glycine/
Computing solvent effects with the 3D-RISM model is illustrated on the glycine example.
All subkeys in the RISM block are discussed in the User's Guide. The things to pay attention to here are SigU and EpsU parameters for each atom in the ATOMS block, the solvent parameters in the SOLVENT sub-block and the FFT box parameters in the SOLUTE sub-block. Both SigU and EpsU values as well as the solvent parameters may be obtained from force field parameter lists. Parameters for some common solvents are available in the ADF User's Guide.
One should take into account the following when choosing FFT box parameters in the SOLUTE block:
$ADFBIN/adf << eor
Title 3D-RISM test
SYMMETRY C(s)
Geometry
Branch Old
End
Define
rco=1.208031
rcoh=1.341959
rcc=1.495685
rcn=1.427005
roh=0.992780
rch1=1.107716
rnh1=1.028574
aoco=123.553475
acco=124.769221
ancc=115.495309
ahoc=105.645766
ach1=107.591718
ah1=109.800726
dch1=123.973836
dh1=57.697485
dc=180.
dn=0.0
doh=0.0
End
ATOMS internal
C 0 0 0 0.00 0.00 0.00 SigU=3.50 EpsU=0.066
O 1 0 0 rco 0.00 0.00 SigU=2.96 EpsU=0.200
O 1 2 0 rcoh aoco 0.00 SigU=2.96 EpsU=0.200
C 1 2 3 rcc acco dc SigU=3.50 EpsU=0.066
N 4 1 2 rcn ancc dn SigU=3.25 EpsU=0.170
H 3 1 2 roh ahoc doh SigU=1.00 EpsU=0.046
H 4 1 2 rch1 ach1 dch1 SigU=1.00 EpsU=0.046
H 4 1 2 rch1 ach1 -dch1 SigU=1.00 EpsU=0.046
H 5 4 1 rnh1 ah1 dh1 SigU=1.00 EpsU=0.046
H 5 4 1 rnh1 ah1 -dh1 SigU=1.00 EpsU=0.046
End
Basis
Type DZP
Core small
End
XC
LDA
End
RISM glycine 1N
RISM1D
subend
SOLVENT water
UNITS uWeight=g/mol ULJsize=A ULJenergy=kcal/mol Ucoord=A Udens=1/A3
Parameters Weight=18.015 nAtoms=2
1 -0.8476 3.166 0.1554 0.000000 0.00000 0.000000
2 0.4238 1.000 0.0460 -0.816490 0.00000 0.577359
0.816490 0.00000 0.577359
DenSpe=0.03333
SUBEND
SOLUTE CO
BOXSIZE 32.0 32.0 32.0
BOXGRID 64 64 64
SUBEND
END
End input
eor
rm TAPE21 logfile
Sample directories:adf/Efield.PntQ_N2/ and adf/Field_PtCO
Two illustrations of applying the very useful BASIS keyword and of application of an Electric Field.
For N2, three calculations are provided: 1) a normal N2 run as a reference with the BASIS keyword, 2) with a homogeneous electric field, 3) with a point charge.
In this example, no Create run is needed in the input file, because the first molecular calculation uses the BASIS keyword. If the $ADFBIN/adf script finds this keyword, it will first generate a new input file which will then be executed. The new input file will contain the required Create run for the N atom in this case. The proper xc functional and relativistic options will automatically be selected by the BASIS keyword. This includes Dirac calculations in case of relativistic runs. The output files is identical to what would have appeared if one would provide the Create runs explicitly in the input file. It also copies the atomic input, so that everything can be checked.
$ADFBIN/adf -n1 << eor title N2 reference for comparison with E-Field runs atoms N 0 0 -.55 N 0 0 +.55 end Basis Type DZP Core Small End end input eor rm TAPE21 logfile $ADFBIN/adf << eor scf conv 1e-8 end title N2 in a homogeneous electric field atoms N 0 0 -.55 N 0 0 +.55 end fragments N t21.N end EField 0 0 0.01 end input eor rm TAPE21 logfile $ADFBIN/adf << eor title N2 polarized by a point charge on the axis EField 0 0 3.0 1.0 end atoms N 0 0 -.55 N 0 0 .55 end Fragments N t21.N end endinput eor
In the second n2 run the homogeneous field is supplied with the key efield, used as simple key: one record, data on the same line as the keyword. The field strength is specified in atomic units.
Homogeneous electric fields can be used to study the polarizability: for sufficiently small fields the dipole moment should respond linearly.
For point charges, the third calculation, the block form of the key efield must be used. The program first tries to find data on the same line as the keyword (defining a homogeneous field). If this is absent, a data block is expected with point-charge specifications: x, y, z and q.
The coordinates are in the same units as in the atoms block (angstrom by default) (but always Cartesian). Q is the charge in elementary units (+1 for a proton).
Point charges can be used for instance to simulate crystal fields (Madelung potential).
Note: the symmetry will be determined automatically by the program as C(lin), rather than D(lin), in the two runs that involve an electric field: the fields break the symmetry.
For PtCO, a fairly large electric field is applied in combination with a tight SCF convergence criterion.
The BASIS keyword in this example illustrates how different choices can be made for different atoms (in this case a frozen core for Pt).
Basis Type DZ Core None Pt Pt.4d END
Sample directory: adf/FDE_H2O_128/
This example demonstrates how to use FDE in combination with a large environment, that is modeled as a superposition of the densities of isolated molecules. Here, the excitation energies of a water molecule surrounded by an environment of 127 water molecules. For details, see C.R. Jacob, J. Neugebauer, L. Jensen, L. Visscher, Phys. Chem. Chem. Phys., 2006 8: 2349.
This calculation consists of two steps:
To reduce the amount of output the next lines are included in the adf calculations:
EPRINT SFO NOEIG NOOVL NOORBPOP SCF NOPOP END NOPRINT BAS FUNCTIONS
First, a prototype water molecule is calculated. The density of this isolated water molecules will afterwards be used to model the environment. Since this molecule will be used as a frozen fragment that is rotated and translated, the option NOSYMFIT has to be included.
$ADFBIN/adf << eor Title Input generated by modco UNITS length bohr angle degree END XC LDA END SYMMETRY NOSYM GEOMETRY sp END SCF iterations 50 converge 1.0e-6 1.0e-6 mixing 0.2 lshift 0.0 diis n=10 ok=0.5 cyc=5 cx=5.0 cxx=10.0 END INTEGRATION 5.0 5.0 FRAGMENTS O t21.DZP.O H t21.DZP.H END ATOMS O -11.38048700000000 -11.81055300000000 -4.51522600000000 H -13.10476265095705 -11.83766918322447 -3.96954531282721 H -10.51089289290947 -12.85330720999229 -3.32020577897331 END ENDINPUT eor mv TAPE21 t21.mol_1
Afterwards, the FDE calculation is performed. In this FDE calculation, there is one nonfrozen water molecule and the previously prepared water molecule is included as a frozen fragment that is duplicated 127 times. For this frozen fragment, the more efficient fitted density is used.
$ADFBIN/adf << eor
Title Input generated by modco
UNITS
length bohr
angle degree
END
XC
MODEL SAOP
END
SYMMETRY NOSYM
SCF
iterations 50
converge 1.0e-6 1.0e-6
mixing 0.2
lshift 0.0
diis n=10 ok=0.5 cyc=5 cx=5.0 cxx=10.0
END
EXCITATION
ONLYSING
LOWEST 5
END
INTEGRATION 4.0 4.0
FRAGMENTS
O t21.DZP.O
H t21.DZP.H
frag1 t21.mol_1 type=fde &
fdedenstype SCFfitted
SubEnd
END
ATOMS
O 0.00000000000000 0.00000000000000 0.00000000000000
H -1.43014300000000 0.00000000000000 1.10739300000000
H 1.43014300000000 0.00000000000000 1.10739300000000
O -11.38048700000000 -11.81055300000000 -4.51522600000000 f=frag1/1
H -13.10476265095705 -11.83766918322447 -3.96954531282721 f=frag1/1
H -10.51089289290947 -12.85330720999229 -3.32020577897331 f=frag1/1
O -1.11635000000000 9.11918600000000 -3.23094800000000 f=frag1/2
H -2.82271357869859 9.71703285239153 -3.18063201242303 f=frag1/2
H -0.12378551814273 10.53819303003839 -2.70860866559857 f=frag1/2
...
O 5.96480100000000 4.51370300000000 3.70332800000000 f=frag1/127
H 5.24291272273548 3.06620845434369 2.89384293177905 f=frag1/127
H 4.73614594944492 5.00201400735317 4.93765482424434 f=frag1/127
END
FDE
PW91K
END
ENDINPUT
eor
Sample directory: adf/FDE_HeCO2_freezeandthaw/
This example demonstrates how a freeze-and-thaw FDE calculation can be performed. As test system, a He-CO2van der Waals complex is used. It will further be shown how different exchange-correlation potential can be used for different subsystems, and how different basis set expansions can be employed. For details, see C.R. Jacob, T.A. Wesolowski, L. Visscher, J. Chem. Phys. 123 (2005), 174104. It should be stressed that the basis set and integration grid used in this example are too small to obtain good results.
Summary:
In the first part, the PW91 functional will be used for both the He and the CO2 subsystems. In this part, the FDE(m) basis set expansion is used, i.e., basis functions of the frozen subsystem are not included in the calculation of the nonfrozen subsystem.
First, the CO2 molecule is prepared. In this calculation, the C2v symmetry of the final complex is used, and the NOSYMFIT option has to be included because this molecule will be rotated as a frozen fragment.
$ADFBIN/adf << eor Title TEST 1 -- Preparation of frozen CO2 Units Length Bohr end Atoms C 0.000000 0.000000 0.000000 O -2.192000 0.000000 0.000000 O 2.192000 0.000000 0.000000 end Symmetry C(2V) NOSYMFIT Fragments C t21.C O t21.O End integration 5.0 xc GGA pw91 end End Input eor mv TAPE21 t21.co2.0
Afterwards, the FDE calculation is performed. In this calculation, the He atom is the nonfrozen system, and the previously prepared CO2 molecule is used as frozen fragment. For this frozen fragment the RELAX option is specified, so that the density of this fragment is updated in freeze-and-thaw iteration (a maximum number of three iteration is specified).
$ADFBIN/adf << eor Title TEST 1 -- Embedding calulation: He + frozen CO2 density -- freeze-and-thaw Units Length Bohr end Atoms He 0.000000 0.000000 6.019000 f=He C 0.000000 0.000000 0.000000 f=co2 O -2.192000 0.000000 0.000000 f=co2 O 2.192000 0.000000 0.000000 f=co2 end Fragments He t21.He co2 t21.co2.0 type=fde & fdeoptions RELAX SubEnd End NOSYMFIT integration 5.0 xc GGA pw91 end FDE PW91K FULLGRID RELAXCYCLES 3 end End Input eor
In this second part, the above example is modified such that PW91 is employed for the CO2 subsystem, while the SAOP potential is used for He. This can be achieved by choosing SAOP in the XC key (this sets the functional that will be used for the nonfrozen subsystem). Additionally, for the frozen fragment the XC option is used to chose the PW91 functional for relaxing this fragment. Furthermore, the PW91 functional is chosen for the nonadditive exchange-correlation functional that is used in the embedding potential with the GGAPOTXFD and GGAPOTCFD options in the FDE key.
$ADFBIN/adf << eor Title TEST 2 -- Embedding calulation: He + frozen CO2 density -- freeze-and-thaw Units Length Bohr end Atoms He 0.000000 0.000000 6.019000 f=He C 0.000000 0.000000 0.000000 f=co2 O -2.192000 0.000000 0.000000 f=co2 O 2.192000 0.000000 0.000000 f=co2 end Fragments He t21.He co2 t21.co2.0 type=fde & fdeoptions RELAX XC GGA PW91 SubEnd End NOSYMFIT integration 5.0 xc MODEL SAOP end FDE PW91K FULLGRID GGAPOTXFD PW91x GGAPOTCFD PW91c RELAXCYCLES 3 end End Input eor
In this third part, the PW91 functional is applied for both subsystems again, but in contrast to part 1, now the FDE(s) basis set expansion is used, i.e., the basis functions of the frozen subsystem are included in the calculation of the nonfrozen subsystem. This can be achieved by employing the USEBASIS option. This option can be combined with the RELAX option.
$ADFBIN/adf << eor Title TEST 3 -- Embedding calulation: He + frozen CO2 density -- freeze-and-thaw Units Length Bohr end Atoms He 0.000000 0.000000 6.019000 f=He C 0.000000 0.000000 0.000000 f=co2 O -2.192000 0.000000 0.000000 f=co2 O 2.192000 0.000000 0.000000 f=co2 end Fragments He t21.He co2 t21.co2.0 type=fde & fdeoptions RELAX USEBASIS SubEnd End NOSYMFIT integration 5.0 xc GGA pw91 end FDE PW91K FULLGRID RELAXCYCLES 3 end End Input eor eor
The example continues with the same calculation where partly the SAOP potential is used.
Sample directory: adf/FDE_Energy_NH3-H2O/
This is example for a calculation of FDE interaction energies in ADF in case of closed shell fragments.
It performs single point runs for H2O and NH3 with LDA/DZ (all-electron) and uses these fragments in:
Integration accuracy is 6.0 which should give total energies for the fragments accurate at least up to 10**(-4) atomic units.
$ADFBIN/adf << EOF
Title H2O LDA/DZ single point
ATOMS
O 1.45838 0.10183 0.00276
H 0.48989 -0.04206 0.00012
H 1.84938 -0.78409 -0.00279
END
SYMMETRY tol=1e-2
BASIS
Type DZ
Core None
END
XC
LDA
END
INTEGRATION
accint 6.0
END
NOSYMFIT
EOF
rm logfile
mv TAPE21 t21.water
EOF
In a similar way the N3 fragment is calculated. Next the FDE calculation is performed. The subkey ENERGY of the key FDE is used, such that the total FDE energy and FDE interaction energy is calculated. First an FDE energy embedding calculation calculation in which the energy of water in presence of a frozen ammonia is computed. This requires a supermolecular integration grid.
$ADFBIN/adf << EOF
Title NH3-H2O LDA/Thomas-Fermi/DZ FDE single point with interaction energy
ATOMS
O 1.45838 0.10183 0.00276 f=frag1
H 0.48989 -0.04206 0.00012 f=frag1
H 1.84938 -0.78409 -0.00279 f=frag1
N -1.51248 -0.03714 -0.00081 f=frag2
H -1.71021 0.95994 -0.11003 f=frag2
H -1.96356 -0.53831 -0.76844 f=frag2
H -1.92899 -0.35123 0.87792 f=frag2
END
SYMMETRY tol=1e-2
FRAGMENTS
frag1 t21.water
frag2 t21.ammonia type=FDE
END
XC
LDA
END
INTEGRATION
accint 6.0
END
EXACTDENSITY
FDE
THOMASFERMI
FULLGRID
ENERGY
END
EOF
Next a fully variational FDE energy calculation (with freeze-and-thaw) is performed.
$ADFBIN/adf << EOF
Title NH3-H2O LDA/Thomas-Fermi/DZ FDE single point with interaction energy
ATOMS
O 1.45838 0.10183 0.00276 f=frag1
H 0.48989 -0.04206 0.00012 f=frag1
H 1.84938 -0.78409 -0.00279 f=frag1
N -1.51248 -0.03714 -0.00081 f=frag2
H -1.71021 0.95994 -0.11003 f=frag2
H -1.96356 -0.53831 -0.76844 f=frag2
H -1.92899 -0.35123 0.87792 f=frag2
END
SYMMETRY tol=1e-2
FRAGMENTS
frag1 t21.water
frag2 t21.ammonia type=FDE &
fdeoptions RELAX
SubEnd
END
XC
LDA
END
INTEGRATION
accint 6.0
END
EXACTDENSITY
SAVE TAPE21
FDE
THOMASFERMI
RELAXCYCLES 3
ENERGY
END
EOF
Sample directory: adf/FDE_Energy_H2O-Ne_unrestricted/
This is example for a calculation of FDE interaction energies in ADF for an open-shell frozen fragment.
It performs single point runs for H2O and Ne, the latter unrestricted with LDA/DZ (all-electron) and uses these fragments in an FDE energy embedding calculation in which the energy of water in presence of a frozen (open-shell) neon atom is computed. This is a bit of an artificial example but it serves its purpose.
No freeze-thaw is done, this is at present not possible with unrestricted (open shell) fragments, but has to be done manually.
Integration accuracy is 6.0 which should give total energies for the fragments accurate at least up to 10**(-4) atomic units.
This test has been checked to yield the same energy as a run with a closed- shell (restricted) Ne atom (just comment UNRESTRICTED in the input below). First the Ne and H2O fragments are calculated.
$ADFBIN/adf << EOF
Title Ne LDA/DZ single point, unrestricted
ATOMS
Ne -1.51248 -0.03714 -0.00081
END
UNRESTRICTED
BASIS
Type DZ
Core None
END
INTEGRATION
accint 6.0
END
SCF
iterations 100
converge 1.0e-06 1.0e-06
END
EXACTDENSITY
NOSYMFIT
EOF
rm logfile
mv TAPE21 t21.ne
EOF
In a similar way the H2O fragment is calculated. Next the FDE calculation is performed. The subkey ENERGY of the key FDE is used, such that the total FDE energy and FDE interaction energy is calculated.
$ADFBIN/adf << EOF
Title Ne-H2O LDA/Thomas-Fermi/DZ FDE single point with interaction energy
ATOMS
O 1.45838 0.10183 0.00276 f=frag1
H 0.48989 -0.04206 0.00012 f=frag1
H 1.84938 -0.78409 -0.00279 f=frag1
Ne -1.51248 -0.03714 -0.00081 f=frag2
END
SYMMETRY tol=1e-2
FRAGMENTS
frag1 t21.water
frag2 t21.ne type=FDE
END
INTEGRATION
accint 6.0
END
SCF
iterations 100
converge 1.0e-06 1.0e-06
END
EXACTDENSITY
FDE
THOMASFERMI
FULLGRID
ENERGY
END
EOF
Sample directory: adf/GO_FDE_H2O-Li/
This examples checks the gradient implementation for FDE. It performs a structure optimization H2O-Li(+) with LDA/DZP.
First, the fragments are made, Li+, and water. Next the FDE geometry optimization is performed with:
$ADFBIN/adf << eor TITLE H2O-Li(+) FDE/LDA/DZP GO New Optimizer starting at too short Li-O distance ATOMS Li 0.000000000000 0.000000000000 -0.054032208082 O 0.000000000000 0.000000000000 -1.534032208080 f=water H -0.778216093965 0.000000000000 -2.135966332900 f=water H 0.778216093965 0.000000000000 -2.135966332900 f=water END CHARGE 1.0 FRAGMENTS Li t21.Li.LDA.DZP water t21.water.LDA.DZP type=fde END XC LDA VWN END FDE ThomasFermi END GEOMETRY Optim Delocalized iterations 15 Converge e=1.0e-3 grad=1.0e-3 END GEOSTEP GradientTerms INTEGRATION 5.0 5.0 5.0 eor
Sample directory: adf/GO_FDE_NH3-H2O/
This examples performs a structure optimization of H2O in presence of frozen NH3 (via optimization of selected coordinates) with LDA and DZ basis. We need a high accint of 6.0 here because the potential energy surface is rather flat and small errors might lead to discrepancies in final structures. It uses (at present) the old branch optimizer for this purpose.
First, the NH3 fragment is made. Next the FDE geometry optimization is performed with:
$ADFBIN/adf << eor
TITLE NH3-H2O dimer FDE LDA DZ structure optimization of H2O
ATOMS
N -1.01393958 -0.15260815 0.00000000 f=nh3
H -1.16290010 -1.15738765 0.00000000 f=nh3
H -1.49925696 0.21074929 0.81414267 f=nh3
H -1.49925696 0.21074929 -0.81414267 f=nh3
O Ox Oy Oz
H H1x H1y H1z
H H2x H2y H2z
END
GEOVAR
Ox 2.25288687
Oy -0.00423586
Oz 0.00000000
H1x 1.28270504
H1y 0.05211069
H1z 0.00000000
H2x 2.54788803
H2y 0.90516678
H2z 0.00000000
END
FRAGMENTS
O t21.O_LDA_DZ
H t21.H_LDA_DZ
nh3 t21.NH3_LDA_DZ type=fde
END
FDE
THOMASFERMI
END
XC
LDA VWN
END
GEOMETRY
Branch Old
Optim Selected
Iterations 100 ! (default is 30)
END
INTEGRATION 6.0 6.0 6.0
eor
Sample directory: adf/FDE_NMR_relax/
This examples demonstrates both the calculation of NMR shieldings using FDE, and how the approximate environment density can be improved by partial relaxation of individual solvent molecules.
The test system is a cluster of acetonitrile and 12 solvent water molecules, of which for two the densities are relaxed, while for the remaining 10 the frozen density of the isolated water is used.
For details, see Refs.
C. R. Jacob, J. Neugebauer, and L. Visscher, A flexible implementation of frozendensity embedding for use in multilevel simulation, submitted, 2007.
R. E. Bulo, Ch. R. Jacob, and L. Visscher, NMR Solvent Shifts of Acetonitrile from Frozen-Density Embedding Calculation, to be submitted, 2007
First, the isolated solvent water molecule is prepared. Again, because this will be rotated and translated afterwards, the option NOSYMFIT has to be included.
$ADFBIN/adf << eor
UNITS
Length Angstrom
Angle Degree
END
ATOMS
O -1.46800 2.60500 1.37700
H -0.95200 3.29800 0.96500
H -1.16100 1.79900 0.96100
END
FRAGMENTS
H t21.H.DZP
O t21.O.DZP
END
XC
LDA
END
INTEGRATION
accint 4.0
END
end input
eor
mv TAPE21 t21.h2o
Afterwards, the FDE calculation is performed. In addition to the nonfrozen acetonitrile molecule, three different fragments are used for the solvent water molecules. The first two fragments frag1 and frag2 are relaxed (in up to two freeze-and-thaw cycles), while the third fragment is used for the remaining 10 solvent molecules. Since a calculation of the shielding is performed afterwards, the option has to be included.
$ADFBIN/adf << eor
Title Input generated by PyADF
UNITS
Length Angstrom
Angle Degree
END
ATOMS
C 0.83000 0.66100 -0.44400
N 0.00000 0.00000 0.00000
C 1.87800 1.55900 -0.81900
H 1.78500 2.40300 -0.13500
H 1.76200 1.94900 -1.83000
H 2.82900 1.12200 -0.51300
O -1.46800 2.60500 1.37700 f=frag1/1
H -0.95200 3.29800 0.96500 f=frag1/1
H -1.16100 1.79900 0.96100 f=frag1/1
O 2.40400 -2.51000 -0.36200 f=frag2/1
H 2.70000 -3.41900 -0.40900 f=frag2/1
H 1.77500 -2.50000 0.35900 f=frag2/1
...
O -3.44400 2.36700 3.13700 f=frag3/10
H -2.70200 2.29200 2.53700 f=frag3/10
H -3.47300 3.29500 3.36800 f=frag3/10
END
FRAGMENTS
H t21.H.DZP
C t21.C.DZP
N t21.N.DZP
frag1 t21.h2o type=FDE &
fdeoptions RELAX
RELAXCYCLES 2
SubEnd
frag2 t21.h2o type=FDE &
fdeoptions RELAX
RELAXCYCLES 2
SubEnd
frag3 t21.h2o type=FDE &
FDEDENSTYPE SCFexact
SubEnd
END
XC
GGA BP86
END
INTEGRATION
accint 4.0
END
SAVE TAPE10
FDE
PW91k
END
End Input
eor
Finally, the calculation of the NMR shielding of the nitrogen atom is performed using the NMR program.
$ADFBIN/adf << eor NMR out tens iso nuc 3 END eor
Sample directory: adf/SUBEXCI_dimer/
This is example for coupled FDE calculation of excitation energies. The subsystem TDDFT code couples the monomer excitations to obtain the excited states of the total system (often denoted as coupled frozen density embedding, FDEc).
First the isolated fragments are calculated, the TAPe21's of these fragments are t21.iso.rho1 and t21.iso.rho2. Next uncoupled FDE excitation energies are calculated in which one fragment is frozen and the other active. The key ALLOW PARTIALSUPERFRAGS is necessary to be able to use subsystem information for only one subsystem from a TAPE21 file of a previous FDE calculation.
$ADFBIN/adf << eor Title Fragment no. 1; relaxed SYMMETRY NOSYM XC GGA BECKE88 PERDEW86 END EXCITATION ONLYSING LOWEST 20 CDSPECTRUM ANALYTIC VELOCITY END ALLOW PARTIALSUPERFRAGS FRAGMENTS rho1 t21.iso.rho1 subfrag=active rho2 t21.iso.rho2 subfrag=active type=fde END ATOMS 1 C 1.05754858422573 -1.70701086799077 -3.50 f=rho1 2 O 2.28164544472573 -1.70701086799077 -3.50 f=rho1 3 C 0.20221626882573 -0.49944933059077 -3.50 f=rho1 4 H 0.49106545372573 -2.67285236319077 -3.50 f=rho1 5 C -1.19491351307427 -0.64892031589077 -3.50 f=rho1 6 C 0.76942690052573 0.78743686120923 -3.50 f=rho1 7 C -2.02186702237427 0.47538393990923 -3.50 f=rho1 8 H -1.62606655117427 -1.65281003349077 -3.50 f=rho1 9 C -0.05719256647427 1.90851291410923 -3.50 f=rho1 10 H 1.85606600152573 0.87752625020923 -3.50 f=rho1 11 C -1.45186603427427 1.75322106580923 -3.50 f=rho1 12 H -3.10608163477427 0.35931991730923 -3.50 f=rho1 13 H 0.37787441672573 2.90894982040923 -3.50 f=rho1 14 H -2.09570824397427 2.63406412680923 -3.50 f=rho1 15 C 2.00708906832899 0.06235850568037 3.50 f=rho2 16 O 2.61913749857899 1.12245748356614 3.50 f=rho2 17 C 0.53364394260760 -0.07460023943380 3.50 f=rho2 18 H 2.56029077395134 -0.91115102374797 3.50 f=rho2 19 C -0.03547527794391 -1.35928561559301 3.50 f=rho2 20 C -0.29722687542061 1.06006167281502 3.50 f=rho2 21 C -1.42262807969967 -1.51329623449550 3.50 f=rho2 22 H 0.61834220104568 -2.23461995830645 3.50 f=rho2 23 C -1.68141695030640 0.90472624158027 3.50 f=rho2 24 H 0.16807297559397 2.04616343352651 3.50 f=rho2 25 C -2.24426699857796 -0.38074233566867 3.50 f=rho2 26 H -1.86422099386266 -2.51028564328820 3.50 f=rho2 27 H -2.33028723444571 1.78172375452933 3.50 f=rho2 28 H -3.32902057100121 -0.49790451479758 3.50 f=rho2 END FDE PW91K END END INPUT eor mv TAPE21 t21.emb.rho1 $ADFBIN/adf << eor Title Fragment no. 2; relaxed SYMMETRY NOSYM XC GGA BECKE88 PERDEW86 END EXCITATION ONLYSING LOWEST 20 CDSPECTRUM ANALYTIC VELOCITY END ALLOW PARTIALSUPERFRAGS FRAGMENTS rho1 t21.emb.rho1 subfrag=active type=fde rho2 t21.iso.rho2 subfrag=active END ATOMS 1 C 1.05754858422573 -1.70701086799077 -3.50 f=rho1 2 O 2.28164544472573 -1.70701086799077 -3.50 f=rho1 3 C 0.20221626882573 -0.49944933059077 -3.50 f=rho1 4 H 0.49106545372573 -2.67285236319077 -3.50 f=rho1 5 C -1.19491351307427 -0.64892031589077 -3.50 f=rho1 6 C 0.76942690052573 0.78743686120923 -3.50 f=rho1 7 C -2.02186702237427 0.47538393990923 -3.50 f=rho1 8 H -1.62606655117427 -1.65281003349077 -3.50 f=rho1 9 C -0.05719256647427 1.90851291410923 -3.50 f=rho1 10 H 1.85606600152573 0.87752625020923 -3.50 f=rho1 11 C -1.45186603427427 1.75322106580923 -3.50 f=rho1 12 H -3.10608163477427 0.35931991730923 -3.50 f=rho1 13 H 0.37787441672573 2.90894982040923 -3.50 f=rho1 14 H -2.09570824397427 2.63406412680923 -3.50 f=rho1 15 C 2.00708906832899 0.06235850568037 3.50 f=rho2 16 O 2.61913749857899 1.12245748356614 3.50 f=rho2 17 C 0.53364394260760 -0.07460023943380 3.50 f=rho2 18 H 2.56029077395134 -0.91115102374797 3.50 f=rho2 19 C -0.03547527794391 -1.35928561559301 3.50 f=rho2 20 C -0.29722687542061 1.06006167281502 3.50 f=rho2 21 C -1.42262807969967 -1.51329623449550 3.50 f=rho2 22 H 0.61834220104568 -2.23461995830645 3.50 f=rho2 23 C -1.68141695030640 0.90472624158027 3.50 f=rho2 24 H 0.16807297559397 2.04616343352651 3.50 f=rho2 25 C -2.24426699857796 -0.38074233566867 3.50 f=rho2 26 H -1.86422099386266 -2.51028564328820 3.50 f=rho2 27 H -2.33028723444571 1.78172375452933 3.50 f=rho2 28 H -3.32902057100121 -0.49790451479758 3.50 f=rho2 END FDE PW91K END END INPUT eor mv TAPE21 t21.emb.rho2
Finally a calculation in which the excitation energies may couple, using the SUBEXCI block keyword.
$ADFBIN/adf << eor Title COUPLED SUBSYSTEM EXCITATIONS SYMMETRY NOSYM XC GGA BECKE88 PERDEW86 END DIFFUSE ALLOW PARTIALSUPERFRAGS FRAGMENTS rho1 t21.emb.rho1 subfrag=active rho2 t21.emb.rho2 subfrag=active type=fde END ATOMS 1 C 1.05754858422573 -1.70701086799077 -3.50 f=rho1 2 O 2.28164544472573 -1.70701086799077 -3.50 f=rho1 3 C 0.20221626882573 -0.49944933059077 -3.50 f=rho1 4 H 0.49106545372573 -2.67285236319077 -3.50 f=rho1 5 C -1.19491351307427 -0.64892031589077 -3.50 f=rho1 6 C 0.76942690052573 0.78743686120923 -3.50 f=rho1 7 C -2.02186702237427 0.47538393990923 -3.50 f=rho1 8 H -1.62606655117427 -1.65281003349077 -3.50 f=rho1 9 C -0.05719256647427 1.90851291410923 -3.50 f=rho1 10 H 1.85606600152573 0.87752625020923 -3.50 f=rho1 11 C -1.45186603427427 1.75322106580923 -3.50 f=rho1 12 H -3.10608163477427 0.35931991730923 -3.50 f=rho1 13 H 0.37787441672573 2.90894982040923 -3.50 f=rho1 14 H -2.09570824397427 2.63406412680923 -3.50 f=rho1 15 C 2.00708906832899 0.06235850568037 3.50 f=rho2 16 O 2.61913749857899 1.12245748356614 3.50 f=rho2 17 C 0.53364394260760 -0.07460023943380 3.50 f=rho2 18 H 2.56029077395134 -0.91115102374797 3.50 f=rho2 19 C -0.03547527794391 -1.35928561559301 3.50 f=rho2 20 C -0.29722687542061 1.06006167281502 3.50 f=rho2 21 C -1.42262807969967 -1.51329623449550 3.50 f=rho2 22 H 0.61834220104568 -2.23461995830645 3.50 f=rho2 23 C -1.68141695030640 0.90472624158027 3.50 f=rho2 24 H 0.16807297559397 2.04616343352651 3.50 f=rho2 25 C -2.24426699857796 -0.38074233566867 3.50 f=rho2 26 H -1.86422099386266 -2.51028564328820 3.50 f=rho2 27 H -2.33028723444571 1.78172375452933 3.50 f=rho2 28 H -3.32902057100121 -0.49790451479758 3.50 f=rho2 END FDE PW91K END SUBEXCI CTHRES 10000.00 SFTHRES 0.00010000 COUPLBLOCK END END INPUT eor
pdb2adf: transforms a PDB file in a QM/MM adf-input file
See the pdb2adf example in the QM/MM manual.
QMMM_Butane: Basic QMMM Illustration
See the QMMM_Butane example in the QM/MM manual.
QMMM_CYT
See the QMMM_CYT example in the QM/MM manual.
QMMM_Surface: Ziegler-Natta catalysis
See the QMMM_Surface example in the QM/MM manual.
Summary:
Geometry optimization of the water molecule, using the (default) local density functional approximation (LDA)
Fair quality basis set: triple zeta with polarization. Four equivalent computations are carried out. The first optimization is done in delocalized coordinates, which requires that atomic coordinates in the input are Cartesian. In the three other optimizations the atomic coordinates are input in Z-matrix format. The optimization is carried out by optimization of the internal coordinates in the second and third calculations, and by optimizing the Cartesian coordinates in the fourth one. In calculation #3 the start-up Hessian is defined in the input file; in #1,2, and 4 the default start-up Hessian (from a force-field approximation) is applied.
As expected all final results are the same, within the range that might be expected from the convergence thresholds (here: the default values).
The '-n' flag, with value one (1) in the commands $ADFBIN/adf is used to control parallelization. 'adf' and other program names that you may find in $ADFBIN are not the executables themselves but (UNIX) scripts to control running the corresponding programs. If a parallel version has been installed and the machine configuration is right, you can carry calculations out in parallel by supplying a suitable value to the -n flag. Omitting the -n flag invokes the default value, which is given by the environment variable $NSCM. Finally, depending on the type of parallel platform, a file $ADFBIN/nodeinfo may define an upper bound on the parallel execution. See the Installation manual for details on the parallel installation.
In all subsequent examples, the set-shell and remove-file commands will be omitted, as well as any -n flags. Also any inputs for create runs will not be shown in other examples except when a special feature is involved or when it may help to clarify the example at hand.
$ADFBIN/adf <<eor
Title WATER Geometry Optimization with Delocalized Coordinates
Atoms
O 0.000000 0.000000 0.000000
H 0.000000 -0.689440 -0.578509
H 0.000000 0.689440 -0.578509
End
Basis
Type TZP
Core Small
End
Geometry
Optim Deloc
End
End Input
eor
A title is supplied. This title is printed in the output header. It is also written to any result files from the calculation and will be printed out when such a file is attached to another calculation, for instance as a fragment file. In addition, adf constructs a 'jobidentification' string that contains the adf release number and the date and time. The jobidentification is also printed in the output header and dumped on any result files.
The atomic positions are given with the key atoms. The Cartesian atomic coordinates are in Angstrom. The structure used here does not necessary imply that the two HO bonds must remain equal in the optimization. The symmetry will keep them equal.
The key geometry must be supplied to let the program do an optimization: otherwise a single point calculation would be carried out. The geometry data block is empty here, meaning that no default values are reset for the options that are controlled with this key.
No symmetry is specified by a Schönfliess type symbol (key symmetry). The program will use the true symmetry of the nuclear frame (accounting for any fields, if present). In this case that is C(2v). If such symmetry would not be acceptable for adf (not all pointgroups are supported!) or when you want to run in a lower symmetry, the symmetry to be used must be specified.
The fragment files are defined implicitly with the Basis keyword. In this case (as well as in most other samples) the fragment files reside in the local directory since they were created there in the same job. If they would have been located elsewhere you could specify a full path for each of the files, or alternatively (if all fragmentfiles are in one single directory) write the directory after the keyword fragments (on the same line).
The precision of numerical integration, to evaluate Hamiltonian matrix elements etc., is not specified and attains therefore the default value (4.0 in an optimization run).
$ADFBIN/adf <<eor Title WATER Geometry Optimization with Internal Coordinates Atoms Z-Matrix 1. O 0 0 0 2. H 1 0 0 rOH 3. H 1 2 0 rOH theta End Basis Type TZP Core Small End GeoVar rOH=0.9 theta=100 End Geometry End End Input eor
The atomic positions are given with the key atoms. Bond lengths are in Angstrom and angles in degrees. The key geometry assigns numerical starting values to the two variables. We could also have written numerical values directly in the atoms block. The structure used here implies that the two HO bonds are equal and must remain equal: they are associated with the same variable; this constraint would not have applied if numerical data had been put in the atoms section, although the symmetry would have kept them equal anyway.
$ADFBIN/adf <<eor Title WATER optimization with (partial) specification of Hessian Atoms Z-Matrix 1. O 0 0 0 2. H 1 0 0 rOH 3. H 1 2 0 rOH theta End GeoVar rOH=0.9 theta=100 End HessDiag rad=1.0 ang=0.1 Fragments H t21.H O t21.O End Geometry End End Input eor
All input is identical to the previous case, except for the key HessDiag. This defines here the start-up Hessian to be diagonal with values 1.0 and 0.1 for the entries related to bondlengths and angles respectively.
$ADFBIN/adf <<eor Title WATER Geometry Optimization in Cartesians Geometry Optim Cartesian End Define rOH=0.9 theta=100 End Atoms Z-Matrix 1. O 0 0 0 2. H 1 0 0 rOH 3. H 1 2 0 rOH theta End Fragments H t21.H O t21.O End End Input eor
In the last calculation the atomic coordinates are input in the same way as before, but the geometry block now specifies, with the subkey optim, that the cartesian coordinates are to be varied and monitored for convergence.
If different coordinates are specified in the optim instruction than were used for the input in the atoms block, no constraints can be used. The variable 'rOH' cannot be placed in the geovar block therefore, since that would imply a constraint: keep the two OH distances equal.
The placement of rOH (and theta) in the define block has a completely different meaning. define merely associates a numerical value with an identifier. Wherever the identifier occurs in input (not only in the gatoms block) it will be replaced by the numerical value. This means that there are now nine (9) variables: the x,y,z coordinates of the three atoms.
Pure translations and rotations will be filtered out by the program and the symmetry (explicitly specified or internally computed), C(2v) here, will be enfored on all developments so that the situation is equivalent to the previous calculation as regards the degrees of freedom of the system.
Remark: the define block must occur in input before the variables defined in it are used. This is one of the few cases where the relative position of keys in the input stream is relevant. The same does not hold for the geovar key used in the earlier example: geovar may be placed anywhere in the input, irrespective of the locations of atoms.
Sample directory: adf/GO_Formaldehyde/
In the input for the optimization run the atomic coordinates are in Z-matrix format while the optimization variables are the Cartesian coordinates. This is achieved with the optim subkey in the geometry block.
A single geovar variable is used for different coordinates. However, since the type of optimization variables (Cartesian) is not the same as the type of input coordinates (Z-matrix), no constraints are implied by this. In fact, the related coordinates do remain equal, but this is because they are symmetry related and the program preserves symmetry anyway.
NonLocal gradient corrections (gga: Generalized Gradient Approximation) according to the approach known as 'Becke' (for exchange) and 'Perdew' (correlation) are included self-consistently with the key xc.
$ADFBIN/adf << eor Title formaldehyde Geometry Optim cartes End XC GGA Becke Perdew END Symmetry C(2v) Atoms Z-matrix 1 O 0 0 0 0.0 0.0 0.0 2 C 1 0 0 r2 0.0 0.0 3 H 2 1 0 r3 a3 0.0 4 H 2 1 3 r3 a3 t4 End Basis Type DZP End Geovar r2 1.94 r3 0.95 a3 120 t4 -180 End integration 4.5 End Input eor
Sample directory: adf/DelocalGO_aspirin/
Geometry optimization of the aspirin molecule, using the delocalized coordinates.
$ADFBIN/adf <<eor
Title geometry optimization of aspirin in delocalized coordinates
Basis
Type DZ
End
GEOMETRY
OPTIM DELOCAL
END
ATOMS
C 0.000000 0.000000 0.000000
C 1.402231 0.000000 0.000000
C 2.091015 1.220378 0.000000
C 1.373539 2.425321 0.004387
C -0.034554 2.451759 0.016301
C -0.711248 1.213529 0.005497
O -0.709522 3.637718 0.019949
C -2.141910 1.166077 -0.004384
O -2.727881 2.161939 -0.690916
C -0.730162 4.530447 1.037168
C -0.066705 4.031914 2.307663
H -0.531323 -0.967191 -0.007490
H 1.959047 -0.952181 -0.004252
H 3.194073 1.231720 -0.005862
H 1.933090 3.376356 -0.002746
O -2.795018 0.309504 0.548870
H -2.174822 2.832497 -1.125018
O -1.263773 5.613383 0.944221
H -0.337334 4.693941 3.161150
H 1.041646 4.053111 2.214199
H -0.405932 3.005321 2.572927
END
END INPUT
Sample directory: adf/ZORA_GO_AuH/
A simple geometry optimization using the scalar relativistic ZORA option.
$ADFBIN/adf << eor title AuH relativistic ZORA optimization integration 5.5 5.5 5.5 atoms Zmat Au 0 0 0 H 1 0 0 1.5 end Basis Type TZ2P Core Small end xc GGA Becke Perdew end relativistic ZORA geometry convergence grad=1e-4 end end input eor
Sample directory: adf/GO_restraint/
The restraint does not have to be satisfied at the start of the geometry optimization. An extra force is added to restrain the bond length, angle, or dihedral angle to a certain value.
Example for angle restraint
$ADFBIN/adf << eor
title WATER geometry optimization with angle restraint
ATOMS
1.O 0.001356 0.000999 0.000000
2.H 0.994442 -0.037855 0.000000
3.H -0.298554 0.948531 0.000000
END
BASIS
Type DZP
END
INTEGRATION 4 4
RESTRAINT
ANGLE 3 1 2 125.0
END
GEOMETRY
END
endinput
eor
Example for bond length restraint
$ADFBIN/adf << eor
title WATER Geometry Optimization with bond length restraint
ATOMS
1.O 0.001356 0.000999 0.000000
2.H 0.994442 -0.037855 0.000000
3.H -0.298554 0.948531 0.000000
END
BASIS
Type DZP
END
INTEGRATION 4 4
RESTRAINT
DIST 1 2 1.03
DIST 1 3 1.03
END
GEOMETRY
END
endinput
eor
Example for dihedral angle restraint
$ADFBIN/adf << eor Title Restraining dihedral of ethane SYMMETRY NOSYM ATOMS 1.C -0.004115 -0.000021 0.000023 2.C 1.535711 0.000022 0.000008 3.H -0.399693 1.027812 -0.000082 4.H -0.399745 -0.513934 0.890139 5.H -0.399612 -0.513952 -0.890156 6.H 1.931188 0.514066 0.890140 7.H 1.931432 0.513819 -0.890121 8.H 1.931281 -1.027824 0.000244 END INTEGRATION 4 4 RESTRAINT DIHED 6 2 1 3 20.00 END BASIS type DZP END GEOMETRY END endinput eor
Sample directory: adf/GO_constraints/
The key CONSTRAINTS can only be used in case of the New branch for optimization of coordinates. The input for this key is very similar to that of the RESTRAINT keyword. The key CONSTRAINTS can, however, also be used to constrain Cartesian coordinates. Note that the key RESTRAINT and freezing of coordinates with the GEOVAR key can also be used in the New branch for optimization of coordinates. In ADF2007 the New branch for optimization can only be used in geometry optimizations and transition state searches. Note that before ADF2008.01 the key CONSTRAINTS was called NEWCONSTRAINTS.
The constraints do not have to be satisfied at the start of the geometry optimization.
Example for angle restraint
$ADFBIN/adf << eor
title WATER geometry optimization with angle constraint
ATOMS
1.O 0.001356 0.000999 0.000000
2.H 0.994442 -0.037855 0.000000
3.H -0.298554 0.948531 0.000000
END
BASIS
Type DZP
END
INTEGRATION 4 4
CONSTRAINTS
ANGLE 3 1 2 125.0
END
GEOMETRY
OPTIM DELOCAL
END
endinput
eor
Example for fixed-atom constraint. Note that the optimization should be done in Cartesian.
$ADFBIN/adf << eor
title WATER geometry optimization with fixed-atom constraint
ATOMS
1.O 0.001356 0.000999 0.000000
2.H 0.994442 -0.037855 0.000000
3.H -0.298554 0.948531 0.000000
END
BASIS
Type DZP
END
INTEGRATION 4 4
SYMMETRY NOSYM
CONSTRAINTS
ATOM 1 0.0 0.0 0.0
ATOM 2 1.0 0.0 0.0
END
GEOMETRY
OPTIM CARTESIAN
BRANCH NEW
END
endinput
eor
Example for bond length restraint.
$ADFBIN/adf << eor
title WATER Geometry Optimization with bond length constraint
ATOMS
1.O 0.001356 0.000999 0.000000
2.H 0.994442 -0.037855 0.000000
3.H -0.298554 0.948531 0.000000
END
BASIS
Type DZP
END
INTEGRATION 4 4
CONSTRAINTS
DIST 1 2 1.03
DIST 1 3 1.03
END
GEOMETRY
OPTIM CARTESIAN
BRANCH NEW
END
endinput
eor
Example for dihedral angle restraint
$ADFBIN/adf << eor Title Constraining dihedral of ethane SYMMETRY NOSYM ATOMS 1.C -0.004115 -0.000021 0.000023 2.C 1.535711 0.000022 0.000008 3.H -0.399693 1.027812 -0.000082 4.H -0.399745 -0.513934 0.890139 5.H -0.399612 -0.513952 -0.890156 6.H 1.931188 0.514066 0.890140 7.H 1.931432 0.513819 -0.890121 8.H 1.931281 -1.027824 0.000244 END INTEGRATION 4 4 CONSTRAINTS DIHED 6 2 1 3 20.00 END BASIS type DZP END GEOMETRY OPTIM DELOCAL END endinput eor
Example for Block constraint (with a dihedral constraint).
$ADFBIN/adf << eor Title Block constraints (with a dihedral constraint) SYMMETRY NOSYM ATOMS 1.C -0.004115 -0.000021 0.000023 b=b1 2.C 1.535711 0.000022 0.000008 b=b2 3.H -0.399693 1.027812 -0.000082 b=b1 4.H -0.399745 -0.513934 0.890139 b=b1 5.H -0.399612 -0.513952 -0.890156 b=b1 6.H 1.931188 0.514066 0.890140 b=b2 7.H 1.931432 0.513819 -0.890121 b=b2 8.H 1.931281 -1.027824 0.000244 b=b2 END INTEGRATION 4 4 CONSTRAINTS DIHED 6 2 1 3 20.00 BLOCK b1 BLOCK b2 END BASIS type DZP END GEOMETRY OPTIM DELOCAL END endinput eor
Sample directory: adf/GO_LiF_Efield/
In the first example a geometry optimization is performed with an external homogeneous electric field. In the second example a geometry optimization is performed with an external point charges
Note that SYMMETRY NOSYM should be used. In case of point charges it is important to use the QPNEAR subkeyword of the INTEGRATION key with a large enough value that would include some of the point charges.
$ADFBIN/adf << eor
Title LiF Cartesian Geometry Optimization in the presence of electric field
Symmetry NOSYM
Atoms
F 0.000000 0.800000 0.000000
Li 0.000000 -0.800000 0.000000
End
Basis
Type TZP
Core Small
End
Geometry
Optim Cartesian
Branch New
Converge 0.0000001
Iterations 100
End
Efield 0.0 0.0 0.01
End Input
eor
$ADFBIN/adf << eor
Title LiF Cartesian Geometry Optimization in the presence of point charges
Symmetry NOSYM
Atoms
F 0.000000 0.800000 0.000000
Li 0.000000 -0.800000 0.000000
End
Basis
Type TZP
Core Small
End
Geometry
Optim Cartesian
Branch New
Converge 0.001
Iterations 100
End
Efield &
0.0 0.0 5.3 0.5
0.0 0.0 -5.3 -0.5
End
integration
qpnear 20
end
End Input
eor
Sample directory: adf/EGO_CH2O_trip_constr/
Example for an excited triplet state geometry optimization with a constraint included.
Needed for such excited state optimizations are the key EXCITATIONS (to calculate excitation energies), the key GEOMETRY (to do a geometry optimization) and the key EXCITEDGO (to select for which excitation a geometryy optimization should be performed). In this case a Z-matrix input for the coordinates is used, and a constraint is used using the GEOVAR keyword.
$ADFBIN/adf <<eor
atoms zmatrix
C 0 0 0 0.0 0.0 0.0
O 1 0 0 R 0.0 0.0
H 1 2 0 1.0 110.0 0.0
H 1 2 3 1.0 110.0 170.0
end
XC
END
GEOVAR
R=1.2 F
END
GEOMETRY
END
BASIS
TYPE DZ
CORE NONE
end
excitations
DAVIDSON
LOWEST 5
onlytriplet
end
EXCITEDGO
STATE A'' 1
TRIPLET
OUTPUT=2
end
SYMMETRY C(S)
eor
The first calculation is a Linear Transit where the Hydrogen atom moves from one side of CN to the other by a parameterized step-by-step change of the angle H-C-N. The other coordinates of the system are optimized along the path.
In the atoms block, one coordinate value is represented by an identifier (th). In the geovar block this is asssigned two values, implying that it is a Linear Transit parameter. The initial and final values for the parameter are given.
Since the geometry block does not have OPTIM SELECTED, all other coordinates are optimized for each of the 10 Linear Transit points.
The subkey iterations in the geometry block carries two arguments: the first is the maximum number of optimization steps (per LT point). The second is the number of LT points to compute in this run: 4. This implies that only a part of the 10-point path defined by the LT parameter(s) will be scanned. The remainder will be done in a follow-up run to illustrate usage of the restart facility.
$ADFBIN/adf <<eor Title HCN Linear Transit, first part NoPrint SFO, Frag, Functions, Computation Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.3 0 0 3 H 1 2 0 1.0 th 0 End Basis Type DZP End Symmetry NOSYM Integration 6.0 6.0 Geometry Branch Old LinearTransit 10 Iterations 30 4 Converge Grad=3e-2, Rad=3e-2, Angle=2 END Geovar th 180 0 End End Input eor mv TAPE21 t21.LT rm logfile
The NoPRINT key turns off a lot of default output. There are several PRINT and NOPRINT options; see the User's Guides for details.
Since the geometry changes from linear to planar (and finally back to linear again), the symmetry must be given explicitly in the input file. Otherwise the program would find a C(lin) symmetry for the initial geometry and assume that this symmetry is preserved throughout. This would of course result in an error abort when the first LT step is carried out, breaking the linear symmetry.
The here specified symmetry (NOSYM: no symmetry at all) is not the true symmetry of the complete path C(s) but a subgroup. It is always allowed to specify a lower symmetry than the actually present symmetry. Such may be necessary (for instance when the true symmetry cannot be handled by adf) or in special cases required for reasons of analysis. Generally speaking, however, we recommend to use the highest symmetry possible (given the case at hand and taking into account the symmetries recognizable by ADF) to boost performance.
Convergence thresholds in the geometry block are set less tight than the defaults: we need only a reasonable estimate of the path, but no highly converged geometries.
At the end of the run the tape21 result file is saved and renamed t21.LT to serve as restart file for the follow-up calculation.
$ADFBIN/adf <<eor Title HCN Linear Transit NoPrint SFO,Frag,Functions,Computation Restart t21.LT Fragments N t21.N C t21.C H t21.H End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.3 0 0 3 H 1 2 0 1.0 th 0 End symmetry NOSYM Integration 6.0 6.0 Geometry Branch Old LinearTransit 10 Converge Grad=3e-2, Rad=3e-2, Angle=2 END Geovar th 180 0 End End Input eor rm TAPE21 logfile
From the restart file, supplied with the key restart, the program reads off that the first 4 points of the LT path have been done already and the scan is continued with LT point #5. The same path definition is supplied again, including the original starting values for the coordinates. The actual starting coordinates (for LT point #5) are read from the restart file. The input values, however, serve to define and verify consistency of the defined LT path and must therefore be supplied correctly.
The key noprint is used to suppress major parts of standard output: all information pertaining to the sfo analysis, all build-from-fragments information, and the lists of elementary functions in the basis sets and fit sets.
From the results of the Linear Transit run we can sketch the energy barrier that H passes over when going from one side of the molecule to the other. This yields a reasonable guess for the Transition State.
To check that the so-obtained estimate is adequate we compute the frequencies in that geometry: one of them should be imaginary.
Apart from serving as a check that the TS estimate is not too bad, the computed Hessian will also serve in the follow-up calculation to obtain the true TS.
$ADFBIN/adf <<eor Title HCN Frequencies in LT max (approx), moderate precision NoPrint SFO,Frag,Functions,Computation Integration 6.0 6.0 Fragments N t21.N C t21.C H t21.H End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Geometry Frequencies End End Input eor mv TAPE21 t21.Freq
Inspection of the output file shows that one of the frequencies is imaginary, as expected (printed as negative), signalling the proximity of the Transition State.
The TAPE21 result file of the calculation is renamed and saved. Later we will use it as a 'restart' file for a TS search, namely to supply the computed Hessian as the initial 'guess' of the Hessian in the (TS) optimization run.
Now carry out the Transition State search, starting from the lt-derived guess.
In this first attempt to find the TS, no use is made of the tape21 result file from the Frequencies run. That will be done in the next calculation.
$ADFBIN/adf <<eor Title HCN Transition State, automatic initial Hessian NoPrint SFO,Frag,Functions,Computation Integration 6.0 6.0 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor rm TAPE21 logfile
The TS-search run type is specified in the geometryblock.
No symmetry is specified; the program determines the symmetry to be C(s) and consequently carries out the ts search in that symmetry.
$ADFBIN/adf <<eor Title HCN Transition State, initial Hessian from Freq run NoPrint SFO,Frag,Functions,Computation Restart t21.Freq Save TAPE13 Integration 6.0 6.0 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor mv TAPE13 t13.TS rm TAPE21 logfile
The CheckPoint file TAPE13, at normal termination automatically deleted by the program, is here saved, using the SAVE key. TAPE13 is as good a restart file as TAPE21 is, but it is a lot smaller. TAPE21 contains a large amount of information for analysis purposes, while TAPE13 contains essentially only restart-type data.
The input is identical to the previous one, except for the restart file. This is used here to provide the Hessian computed in the Frequencies run as the start-up Hessian for the ts optimization. At the same time the atomic coordinates are read off from the restart file and override the values in the input file. This latter aspect could have been suppressed; see the User's Guide for using the restart key.
Finally the ts search where one coordinate is kept frozen, to illustrate a constrained optimization.
$ADFBIN/adf <<eor Title HCN constrained TS search NoPrint SFO,Frag,Functions,Computation Restart t21.Freq Integration 6.0 6.0 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 rNC 0 0 3 H 1 2 0 1.223 70 0 End GeoVar rNC=1.186 F End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor rm TAPE21 logfile rm t21.Freq
The geovar key specifies that the nc distance, rNC has the initial value 1.15 and remains frozen ('F').
The fact that the optimization is now carried out in a different subspace of atomic coordinates does not prevent us from using the t21.Freq restart file to supply the initial Hessian.
The IRC calculation is split in three steps, to illustrate the Restart facility applied to the IRC functionality.
In the first only a few points are computed, along one of the two paths leading from the TS to the adjacent minima. Since no explicit directives are given in the input to specify the direction of the first path, the so-called 'forward' path is taken. The definition of which is 'forward' and which is 'backward' is in fact quite arbitrary and is determined by the program. See the User's Guide for details.
The saved TAPE13 file from one of the TS calculations is used as restart file. This provides (a) the optimized coordinates of the TS as starting point, (b) the initial Hessian to guide the point-by-point optimizations along the IRC path, and (c) the eigenvector of the lowest Hessian eigenvalue to define the initial direction of the IRC path.
The TAPE13 file from this partial IRC scan is saved to serve as start-up file for the next calculations, which will continue the IRC scan.
In the Geometry key block, the run type is set to IRC and the 'Points' option is used to limit the number of IRC points to compute.
$ADFBIN/adf << eor Title HCN IRC partial path (forward) NoPrint SFO,Frag,Functions, Computation Integration 6.0 6.0 Restart t13.TS Save TAPE13 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=5 End End Input eor mv TAPE13 t13.IRC_1 rm TAPE21 logfile
The IRC is continued in the next calculation, using the TAPE13 file from the previous one as restart file. From this file, the program reads the IRC path information computed sofar. By default, it would continue on the 'forward' path, since that was not yet finished. However, in the Geometry key block, we now specify not only that a limited number of points is to be computed in this run (5 again), but we instruct the program also to compute only points on the 'backward' path.
$ADFBIN/adf << eor Title HCN IRC partial part (backward) NoPrint SFO,Frag,Functions, Computation Restart t13.IRC_1 Save TAPE13 Integration 6.0 6.0 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 END Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=5 Backward End End Input eor mv TAPE13 t13.IRC_2 rm TAPE21 logfile
In the third IRC run, the IRC scan is finished. We start with the TAPE13 file from the previous run and set a maximum of 70 IRC points to compute (which turns out to be sufficient for the complete IRC scan). The program starts on the forward path, continuing where the first (not the previous) had stopped after 5 points, completes the forward path, and then continues on the backward path, starting where the second IRC run had stopped. Both paths are finished and a summary of the path characteristics is printed in the final part of the output.
$ADFBIN/adf << eor Title HCN IRC completion NoPrint SFO,Frag,Functions, Computation Restart t13.IRC_2 Integration 6.0 6.0 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=70 End End Input eor
Sample directory: adf/HCN_CINEB/
This example demonstrates the use of the Nudged Elastic Band method in ADF for finding a transition state of the HCN isomerisation reaction. A shell script used to run the example calculation is shown below:
$ADFBIN/adf <<eor
TITLE Test of the CI-NEB method
SYMMETRY C(S)
NOPRINT SCF SFO
UNITS
length Angstrom
angle Degree
END
ATOMS
1.C 0.000000 0.000000 0.000000
2.N XN 0.000000 0.000000
3.H XH YH 0.000000
END
GEOVAR
XN 1.180 1.163
XH 2.196 1.831 1.006 0.105 -0.718 -1.078
YH 0.000 0.799 1.122 1.163 0.813 0.000
END
BASIS
END
GEOMETRY
CINEB 9
iterations 150
OPTIM selected
converge grad=0.001
nebspring 1 0.06
END
integration 4.0
SCF
Convergence 0.00000001
END
eor
A few important points to note:
Sample directory: adf/TS_C2H6/
Frequently when searching for a transition state, one needs an accurate second derivatives matrix, a Hessian. An exact Hessian may be obtained analytically but this may be very expensive for large molecules. In such cases it may be beneficial to calculate Hessian matrix elements only for atoms directly involved in the reaction for which a transition state is sought for. The rest of the Hessian can then be approximated using a cheaper method.
In this example, a saddle point of the ethane internal rotation around C-C bond is found. In principle, only hydrogen atoms contribute to the normal mode we are interested in. Therefore we calculate a partial Hessian matrix including hydrogen atoms only. For this purpose, the AnalyticalFreq block key is used. In this block, a NUC keyword is added specifying that the second derivatives are calculated for atom 3 (and its symmetry-equivalents) only. Note that the Hessian matrix elements between symmetry-equivalent atoms, for example between 3,H and 4.H are also calculated. The rest of the matrix is estimated using the default method.
$ADFBIN/adf <<eor TITLE Ethane transition state search using partial Hessian ATOMS 1 C 0.000000000000 0.000000000000 0.767685465031 2 C 0.000000000000 0.000000000000 -0.767685465031 3 H 0.964354016767 0.347635559279 1.177128271450 4 H -0.181115782790 -1.008972856410 1.177128271450 5 H -0.783238233981 0.661337297125 1.177128271450 6 H -0.500471876676 0.894626767091 -1.177128271450 7 H -0.524533568868 -0.880734742626 -1.177128271450 8 H 1.025005445540 -0.013892024465 -1.177128271450 END BASIS type DZ core Large END AnalyticalFreq NUC 3 End INTEGRATION 5.0 5.0 5.0 eor
After the Hessian is calculated, the resulting TAPE21 file is used for a subsequent transition state search:
mv TAPE21 ethane-frq.t21
$ADFBIN/adf <<eor
TITLE Ethane transition state search using partial Hessian
ATOMS
1 C 0.000000000000 0.000000000000 0.767685465031
2 C 0.000000000000 0.000000000000 -0.767685465031
3 H 0.964354016767 0.347635559279 1.177128271450
4 H -0.181115782790 -1.008972856410 1.177128271450
5 H -0.783238233981 0.661337297125 1.177128271450
6 H -0.500471876676 0.894626767091 -1.177128271450
7 H -0.524533568868 -0.880734742626 -1.177128271450
8 H 1.025005445540 -0.013892024465 -1.177128271450
END
Fragments
H t21.H
C t21.C
END
GEOMETRY
smooth conservepoints
TransitionState mode=1
optim All Cartesian
iterations 30
step rad=0.15
hessupd BOFILL
converge e=1.0e-4 grad=1.0e-3 rad=1.0e-3
END
RESTART ethane-frq.t21
INTEGRATION 5.0 5.0 5.0
eor
Important note: care should be taken to specify correct mode in the TransitionState keyword. Because a significant part of the Hessian will not be calculated exactly, it is possible that it will have more than one negative eigenvalue, in which case the one we are interested in may not be the first one. In such a case, one needs to specify the correct mode number in the TransitionState keyword.
Sample directory: adf/TS_CH4_HgCl2/
A ZORA scalar relativistic Transition State calculation.
$ADFBIN/adf << eor TITLE Transition State: CH4 + HgCl2 ↔ CH3HgCl + HCl noprint sfo,frag print atdist GEOMETRY TransitionState END TSRC DIST 1 5 1.0 END Relativistic scalar ZORA Basis Type TZP Core small End ATOMS C 0.049484 0.042994 0.000000 H -0.068980 0.638928 -0.915972 H -0.068980 0.638928 0.915972 H -0.841513 -0.626342 0.000000 H 0.555494 -1.148227 0.000000 Hg 2.303289 -0.007233 0.000000 Cl 4.429752 0.776056 0.000000 Cl 1.342057 -2.676083 0.000000 END endinput eor
For the density-functional the Local Density approximation is used (no GGA corrections).
At each geometry cycle the interatomic distance matrix is printed (print atdist).
The initial geometry is a reasonable but not very accuracte estimate of the Transition State. The program needs quite a few cycles to converge, which is rather typical for TS searches: they are a lot more tricky and fail more often than a simple minimization. The TSRC key is used to specify a reaction coordinate along which the transition state is sought for. This feature is especially useful when an accurate Hessian is not available.
Sample directory: adf/TSRC_SN2/
With the TSRC key one can specify a reaction coordinate along which the transition state is sought for. This feature is especially useful when an accurate Hessian is not available.
This example tries to find the TS for the SN2 reaction of F- + CH3Cl ⇔ CH3F + Cl-
$ADFBIN/adf << eor Title TransitionState search for Sn2 reaction of F- + CH3Cl ANALYTICALFREQ END XC LDA VWN GGA OPBE END TSRC dist 1 5 1.0 dist 1 6 -1.0 END ATOMS C 0.000000 0.000000 0.000000 H -0.530807 0.919384693 0.012892 H -0.530807 -0.919384693 0.012892 H 1.061614 0.000000 0.012892 Cl 0.000000 0.000000 -2.124300 F 0.000000 0.000000 2.019100 END Geometry TS End BASIS type TZ2P core NONE createoutput none END INTEGRATION 6.0 6.0 Charge -1 endinput eor
Sample directory: adf/LT_constraint/
The LINEARCONSTRAINTS keyword allows geometry optimizations (old branch) with constraints defined by arbitrary linear combinations of (internal) coordinates. The constraint has to be satisfied at the start of the geometry optimization. Note that the before ADF2008.01 the key LINEARCONSTRAINTS was called CONSTRAINT.
Example for bond length constraint, where at the start of the linear transit rOH1=R1=1.0, and rOH2=R2=1.5, such that (-1.0)*R1+(1.0)*R2=0.5, and in the final geometry -R1+R2=0.0 (Reactcoord 0.5 0.0)
$ADFBIN/adf << eor title constraint keyword XC GGA Becke Perdew END Geometry Branch Old LinearTransit 6 End Integration 5.0 Atoms Internal O 0 0 0 H 1 0 0 R1 H 1 2 0 R2 109.9 End GeoVar R1 =1. R2 =1.5 End LinearConstraints ReactCoord 0.5 0.0 R1 -1.0 R2 1.0 SubEnd End Basis Type DZP End end input eor
Sample directory: adf/Transit_H2O/
In ADF2008.01 a transit calculation option has been added in the new optimization branch. This is capable of performing both linear transits, and non-linear transits, and is the default when the LINEARTRANSIT or TRANSIT sub-block is included in the 'Geometry' block.
The new transit code works differently to the old: the transit is represented as a sequence of constrained optimizations. A 'Constraints' block is used to delineate the constraints applied at each stage of the transit.
Non-linear transits are possible, and can even be combined with linear transits in other coordinates. To perform a non-linear transit in a particular coordinate, explicit values must be given.
$ADFBIN/adf << eor
Title WATER Transit (non-linear), with the new optimizer branch
Atoms
O 0.000000 0.000000 0.000000
H 0.000000 -0.689440 -0.578509
H 0.000000 0.689440 -0.578509
End
Symmetry NOSYM
Constraints
dist 1 2 0.8 1.0 1.25 1.5
angle 2 1 3 start=100.0 end=120.0
End
Basis
Type SZ
Core Large
End
Integration 4.0 2.0
Geometry
Transit 4
Optim Deloc
Converge 0.0001
End
End Input
In the example above, 4 values are given for the distance between atoms 1 and 2. This distance constraint will be applied simultaneously with the linear transit constraints for the angle, with other degrees of freedom optimized at each stage of the transit.
Finally, it should be pointed out that 'partial constraints' are used by default in the transit calculations. These constraints are not required to be fully met at each intermediate geometry, but are fully met at the converged geometries. You can use fully converged constraints by supplying the FULLCONVERGE option to the 'Constraints' subblock of the 'Geometry' block (not to be confused with the 'Constraints' block at root level).
Sample directory: adf/Energy_H2O/
If the TOTALENERGY is included the total energy will be calculated.
This example performs single point runs for H2O with PBE/DZP with frozen cores and all-electron and B3LYP/DZP with all-electron and HARTREEFOCK/DZP with all-electron The tests run in C(2v) symmetry. Integration accuracy is 6.0 which should give total energies accurate at least up to 10-4 atomic units. The key EXACTDENSITY is used for higher accuracy of the results.
First example:
$ADFBIN/adf <<eor Title H2O PBE/DZP (frozen core) single point calculation ATOMS O 0.00000 0.00000 0.00000 H 0.00000 -0.68944 -0.57851 H 0.00000 0.68944 -0.57851 END BASIS Type DZP Core Small END XC GGA PBE END INTEGRATION 6.0 EXACTDENSITY TOTALENERGY eor
Note that only energy difference comparisons are meaningful. These are the only energies that play a role in chemistry of course, and for this one does not need total energies.
Sample directory: adf/SD_CrNH3_6/
The computation of multiplet states corresponding to an open-shell system can be carried out with ADF by first computing the 'Average-of-Configuration' (aoc) state, where all orbitals in the open shell are degenerate and equally occupied. This computation is spin-restricted and serves as a fragment file for the multiplet run, where then different occupation numbers are assigned to the various orbitals in the open shell. The corresponding energies are computed in the field of the aoc, which is achieved by not iterating the self-consistency equations to convergence but only computing the orbitals in the initial field.
Since ADF requires that all symmetry-partners in an irreducible representation (irrep) have equal occupations, the multiplet calculation, where such orbitals are not equally occupied, must be carried out in a formally lower pointgroup symmetry. The pointgroup to select and the appropriate occupation numbers to apply must be worked out by the user 'on paper' in advance. An auxiliary program asf, developed by the group of Claude Daul in Fribourg can be used to determine which calculations are needed, and how to compute the multiplet energies from the results. See the discussion of Multiplet energies in the Theory document.
The script starts with the 'creation' of the required basic atoms, N, H, Cr using a fair basis set quality.
The next step is the computation of the ammonia fragment NH3. This is not a crucial step here: the multiplet state computation can equally well be carried out by not using any intermediate compound fragments. However, it illustrates once more how a bigger molecule can be built up from smaller, but not trivial fragments.
$ADFBIN/adf <<eor title AMMONIA NOPRINT sfo,frag,functions define xH=0.95522523 yH=xH*sqrt(3)/2 zH=0.3711068 end atoms N 0 0 0 H -xH 0 zH H xH/2 -yH zH H xH/2 yH zH end Basis Type TZP Core Small End symmetry C(3V) endinput eor mv TAPE21 t21.NH3
The input of the atomic coordinates uses expressions, in this case to enforce exact symmetry relations that would otherwise require 14-digit input values or some inaccuracy. The symmetry specification is redundant: the program would also find it by itself.
The next step is to compute the reference state, with respect to which we will later compute the multiplet states. The reference state is the so-called 'Average-of-configuration' (aoc) state. The result file (TAPE21) of this calculation will be used as a fragment file.
$ADFBIN/adf <<eor title Cr(NH3)6 : Average-of-Configuration run COMMENT using NH3-fragments END symmetry D(3d) scf iterations 25 mix 0.15 end atoms Cr 0.000000 0.000000 0.000000 N 0.000000 1.714643 1.212436 f=NH3/1 H 0.000000 1.466154 2.206635 f=NH3/1 H -0.827250 2.293404 1.036727 f=NH3/1 H 0.827250 2.293404 1.036727 f=NH3/1 N -1.484924 -0.857321 1.212436 f=NH3/2 H -1.269726 -0.733077 2.206635 f=NH3/2 H -1.572521 -1.863121 1.036727 f=NH3/2 H -2.399771 -0.430282 1.036727 f=NH3/2 N 1.484924 -0.857321 1.212436 f=NH3/3 H 1.269726 -0.733077 2.206635 f=NH3/3 H 2.399771 -0.430282 1.036727 f=NH3/3 H 1.572521 -1.863121 1.036727 f=NH3/3 N 0.000000 -1.714643 -1.212436 f=NH3/4 H 0.000000 -1.466154 -2.206635 f=NH3/4 H 0.827250 -2.293404 -1.036727 f=NH3/4 H -0.827250 -2.293404 -1.036727 f=NH3/4 N 1.484924 0.857321 -1.212436 f=NH3/5 H 1.269726 0.733077 -2.206635 f=NH3/5 H 1.572521 1.863121 -1.036727 f=NH3/5 H 2.399771 0.430282 -1.036727 f=NH3/5 N -1.484924 0.857321 -1.212436 f=NH3/6 H -1.269726 0.733077 -2.206635 f=NH3/6 H -2.399771 0.430282 -1.036727 f=NH3/6 H -1.572521 1.863121 -1.036727 f=NH3/6 H -1.572521 1.863121 -1.036727 f=NH3/6 end fragments Cr t21.Cr NH3 t21.NH3 end occupations A1.G 8.75 A2.G 2 E1.G 16 1.5 0.75 A1.U 2 A2.U 8 E1.U 20 END end input eor mv TAPE21 t21.CrA6ES
Occupation numbers are specified, to make certain what the reference state is that we will start from in the subsequent calculations. The result file TAPE21 is saved to serve as fragment file in the subsequent calculations.
Now, we proceed with the multiplet calculations. In the example they are combined in one single run, but they could also be evaluated in separate runs. For each calculation it is required to:
a) Use the aoc TAPE21 file as fragment file
b) Choose which molecular orbitals in the open shell to occupy: select the appropriate pointgroup symmetry and the UnRestricted key if necessary and specify the occupation numbers, using the irreducible representations of the selected point group.
The results are one-determinant calculations, which must then, later, be combined analytically to obtain the required multiplet energy values.
$ADFBIN/adf <<eor
title Cr(NH3)6 : SlaterDeterminants run
NOPRINT frag
symmetry C(I) ! lower symmetry
scf
iterations 0
end
atoms
Cr 0.000000 0.000000 0.000000 f=CrA6
N 0.000000 1.714643 1.212436 f=CrA6
H 0.000000 1.466154 2.206635 f=CrA6
H -0.827250 2.293404 1.036727 f=CrA6
H 0.827250 2.293404 1.036727 f=CrA6
N -1.484924 -0.857321 1.212436 f=CrA6
H -1.269726 -0.733077 2.206635 f=CrA6
H -1.572521 -1.863121 1.036727 f=CrA6
H -2.399771 -0.430282 1.036727 f=CrA6
N 1.484924 -0.857321 1.212436 f=CrA6
H 1.269726 -0.733077 2.206635 f=CrA6
H 2.399771 -0.430282 1.036727 f=CrA6
H 1.572521 -1.863121 1.036727 f=CrA6
N 0.000000 -1.714643 -1.212436 f=CrA6
H 0.000000 -1.466154 -2.206635 f=CrA6
H 0.827250 -2.293404 -1.036727 f=CrA6
H -0.827250 -2.293404 -1.036727 f=CrA6
N 1.484924 0.857321 -1.212436 f=CrA6
H 1.269726 0.733077 -2.206635 f=CrA6
H 1.572521 1.863121 -1.036727 f=CrA6
H 2.399771 0.430282 -1.036727 f=CrA6
N -1.484924 0.857321 -1.212436 f=CrA6
H -1.269726 0.733077 -2.206635 f=CrA6
H -2.399771 0.430282 -1.036727 f=CrA6
H -1.572521 1.863121 -1.036727 f=CrA6
end
fragments
CrA6 t21.CrA6ES
end
UnRestricted
SlaterDeterminants
Check AOC
A1.g 4 0.375 // 4 0.375
A2.g 1 // 1
E1.g:1 4 0.375 0.1875 // 4 0.375 0.1875
E1.g:2 4 0.375 0.1875 // 4 0.375 0.1875
A1.u 1//1
A2.u 4//4
E1.u:1 5//5
E1.u:2 5//5
SUBEND
State1
A1.g 4 1 // 4 1
A2.g 1 // 1
E1.g:1 4 0 0 // 4 0 1
E1.g:2 4 0 0 // 4 0 0
A1.u 1//1
A2.u 4//4
E1.u:1 5//5
E1.u:2 5//5
SUBEND
State2
A1.g 4 1 // 4 1
A2.g 1 // 1
E1.g:1 4 0 0 // 4 1 0
E1.g:2 4 0 0 // 4 0 0
A1.u 1//1
A2.u 4//4
E1.u:1 5//5
E1.u:2 5//5
SUBEND
State3
A1.g 4 1 // 4 1
A2.g 1 // 1
E1.g:1 4 0 1 // 4 0 0
E1.g:2 4 0 0 // 4 0 0
A1.u 1//1
A2.u 4//4
E1.u:1 5//5
E1.u:2 5//5
SUBEND
end
end input
eor
The SlaterDeterminants block may contain any number of sub blocks, each starting with an (arbitrary) title record, followed by a set of occupation numbers and closed by a SubEnd record. Each such subkey block specifies a single one-determinant-state calculation. All occupation numbers must reference the irreps of the specified pointgroup symmetry, C(I) in the example, and must be just a reassignment of the electrons that are equally distributed over the corresponding degenerate irreps in the reference aoc calculation.
The so-obtained energies of the one-determinant states can now be combined to calculate the desired multiplet energies. See the Theory document and the adf User's Guide.
Note carefully that in the calculation of the SingleDeterminants, the scf procedure is prevented to cycle to convergence by setting the subkey Iterations to zero in the SCF data block.
Sample directory: adf/CuH+_S-squared/
Example calculates expectation value of S2 (< S2 >) of CuH+ in various symmetries, using unrestricted density functional theory. Last example in this example file calculates this value in the case there are more beta electrons than alpha electrons.
$ADFBIN/adf << eor Title calculate expectation value of S-squared ATOMS Z-Matrix Cu 0 0 0 H 1 0 0 1.463 END CHARGE 1.0 -1.0 Unrestricted FRAGMENTS H t21.H Cu t21.Cu END endinput eor
Sample directory: adf/ModStPot_N2+/
This calculation illustrates:
$ADFBIN/adf <<eor title N2+ hole localization atoms N 0 0 -2.0 N 0 0 2.0 end Basis Type DZP Core Small End symmetry C(lin) ! allow symmetry breaking unrestricted Occupations keeporbitals=3 & ! keeporbitals: let the density relax a bit, then fix the MO occups sigma 3 // 1 0 1 pi 2 // 2 end CHARGE 1 1 ! this duplicates info from "OCCUPATIONS" (check) modifystartpotential ! to break the symmetry in the start-up potential N/1 0.5 0.5 N/2 4 1 end end input eor
The purpose of this run is to compute the N2+ ion, with the hole localized on one of the atoms. In a very small system like N2+ this is a tricky thing to do. The program has a tendency towards the symmetric solution, with the hole delocalized. A few trial runs, just putting a net +1 charge into the system, will reveal that clearly.
To achieve the desired situation we apply the key modifystartpotential to break the symmetry of the initial potential. A potential is generated as if the electronic cloud in the second N fragment is spin-polarized in a ratio 4:1 (this precise value is not very relevant), which achieves that initially a non-symmetric solution is obtained. The symmetry must be specified, lest the program determine and use the higher symmetry from the nuclear frame. This would prevent any symmetry breaking altogether.
Next, in order to prevent that the system relaxes to the symmetric situation, we apply the keeporbitals option of the occupations key. This fixes the occupied orbitals in the sense that in each scf cycle the program will try to keep the electrons in orbitals that resemble the previously occupied orbitals as much as possible.
The key modifystartpotential here demonstrated has a more relevant and less unstable application in larger systems. See the User's Guide for references.
Sample directory: adf/Fe4S4_BrokenSymm/
This calculation shows a spin-flip restart feature that allows to exchange alpha and beta fit coefficients for selected atoms upon restart. First the high spin configuration with 8 more α-electrons than beta-electrons is calculated (Sz=4). Next the broken spin-symmetry configuration is calculated (Sz=0), using the subkey spinflip in the restart key. In this case the spin will be flipped for iron atoms 1 and 2.
$ADFBIN/adf <<eor TITLE Fe4S4 High-spin configuration ATOMS Fe -0.000000000000 -1.256142548900 0.888226914500 Fe 0.000000000000 1.256142548900 0.888226914500 Fe -1.256142548900 0.000000000000 -0.888226914500 Fe 1.256142548900 -0.000000000000 -0.888226914500 S -1.845393493800 0.000000000000 1.304890253400 S 1.845393493800 -0.000000000000 1.304890253400 S -0.000000000000 -1.845393493800 -1.304890253400 S 0.000000000000 1.845393493800 -1.304890253400 END Symmetry C(2v) CHARGE 0.0 8.0 UNRESTRICTED BASIS type DZ core Large createoutput None END XC GGA OPBE END end input eor mv TAPE21 Fe4S4-high-spin.t21 $ADFBIN/adf <<eor TITLE Fe4S4 LOW-spin configuration Restart Fe4S4-high-spin.t21 & ! Make sure atoms specified in the SpinFlip keyword are symmetry-equivalent SpinFlip 1 2 End ATOMS Fe -0.000000000000 -1.256142548900 0.888226914500 Fe 0.000000000000 1.256142548900 0.888226914500 Fe -1.256142548900 0.000000000000 -0.888226914500 Fe 1.256142548900 -0.000000000000 -0.888226914500 S -1.845393493800 0.000000000000 1.304890253400 S 1.845393493800 -0.000000000000 1.304890253400 S -0.000000000000 -1.845393493800 -1.304890253400 S 0.000000000000 1.845393493800 -1.304890253400 END Symmetry C(2v) CHARGE 0.0 0.0 UNRESTRICTED BASIS type DZ core Large createoutput None END XC GGA OPBE END eor
Sample directory: adf/CEBE_NNO/
ADF is well suited for calculating Core Electron Binding Energies (CEBEs). In this example it is shown how one can differentiate between the 1s CEBEs of the two non-equivalent nitrogen atoms in N2O, using a delta-SCF technique. It starts with a regular calculation that has the purpose of preparing a reference TAPE21 file for the NNO molecule, which will later be useful in the energy analysis. The result file is saved to t21.NNO.
The same GGA functional is specified throughout the run. The amount of output is reduced by using some print keys.
The prepare the nitrogen atom with a core hole (restricted) will be used as a fragment later. This enables selection of where the core hole should be.
$ADFBIN/adf -n1 << eor title N atom core hole ATOMS N 0.0 0.0 0.0 end Basis Type TZ2P Core None end xc gradients pw86x pw91c end occupations s 1 2 p 3 end end input eor mv TAPE21 t21.N_ch
Now perform the restricted ground state molecule for analysis later. The TAPE21 result file is saved.
$ADFBIN/adf << eor title NNO noprint sfofragpop fragsfo xc gradients pw86x pw91c end ATOMS N 0.0 0.0 -1.1284 N 0.0 0.0 0.0 O 0.0 0.0 1.1841 end Basis Type TZ2P Core None createoutput None end end input eor mv TAPE21 t21.NNO
Next follow two sets of almost identical calculations in which a 1s electron is removed from one or the other N atom (please note that the deepest s level is associated with the 1s of the oxygen atom). The molecular NNO result file is used as fragment. An unrestricted calculation is done and a positive charge is specified. The final result file for the molecule with the core hole is saved. Then another calculation is done to conveniently obtain the energy with respect to the normal molecule. This is repeated for a core hole on the other N atom.
$ADFBIN/adf <<eor unrestricted unrestricted core hole noprint sfofragpop fragsfo ATOMS N 0.0 0.0 -1.1284 f=N_ch N 0.0 0.0 0.0 O 0.0 0.0 1.1841 end xc gradients pw86x pw91c end Basis Type TZ2P Core None createoutput None end fragments N_ch t21.N_ch end charge 1 1 occupation sigma 1 1 1 4 // 1 0 1 4 pi 4 // 4 end end input eor mv TAPE21 t21.NNO.unr1
In the second calculation the result file of one of the unrestricted NNO
calculations is used as restart file, which ensures that the hole stays at its
place, because the starting density is already correct. The result file t21.NNO
for the normal NNO calculation is specified as fragment to serve as an energy
reference. The final Bonding Energy printed by ADF indicates what the CEBE is.
However, please check
Chong, D.P.
Accurate DFT Calculation of Core-Electron Binding Energies
in Reviews in Modern Quantum Chemistry, A Celebration of the Contributions of R.G. Parr,
edited by K.D. Sen (World Scientific Publishing Co., Singapore), 1106-1139 (2002)
for more detailed information on Core-Electron Binding Energies.
This reference also contain infomration on empirical corrections that may have to be made on the final numbers.
$ADFBIN/adf <<eor title NNO unr. core hole noprint sfofragpop fragsfo xc gradients pw86x pw91c end restart t21.NNO.unr1 ATOMS N 0.0 0.0 -1.1284 f=NNO N 0.0 0.0 0.0 f=NNO O 0.0 0.0 1.1841 f=NNO end fragments NNO t21.NNO end unrestricted charge 1 1 occupation sigma 1 1 1 4 // 1 0 1 4 pi 4 // 4 end end input eor
Similarly, one could easily have prepared an oxygen with a core hole and determined the CEBE of the oxygen 1s atom.
Sample directory: adf/Freq_NH3/
Summary:
Computation of frequencies by symmetric displacements. The assumed equilibrium input structure should be given in Cartesian coordinates.
The symmetry is determined automatically by the program as C(3v), from the input coordinates. During the calculation first symmetric atomic displacements are constructed. The number of such displacements in each irreducible representation corresponds to the number of frequencies with the corresponding symmetry. All displaced geometries within one representation have the same symmetry, which enables us to use it to speed up the computation significantly. Another advantage of having the same symmetry is that the numerical integration data can be reused efficiently (see SMOOTH option) thus reducing the level of numerical noise in gradients and force constant matrix.
$ADFBIN/adf <<eor title NH3 frequencies in symmetric displacements atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Basis Type TZP Core Small End geometry frequencies Symm end thermo T=300,400 integration 5.0 end input eor
Computation of frequencies by Cartesian displacements. The assumed equilibrium input structure is given in internal coordinates. A dummy atom is used for a convenient definition of the Z-matrix such that it reflects the pointgroup symmetry C(3v).
$ADFBIN/adf <<eor title NH3 frequencies define rNH=1.02 theta=112 phi=120 end atoms Z-matrix XX 0 0 0 N 1 0 0 1.0 H 2 1 0 rNH theta H 2 1 3 rNH theta phi H 2 1 4 rNH theta phi end Basis Type TZP Core Small End geometry optim cartesian frequencies end thermo T=300,400 integration 5.0 end input eor
The symmetry is determined automatically by the program as C(3v), from the input coordinates. In a Frequencies calculation the symmetry (specified on input or computed internally) is used for analysis and in some cases to speed up the calculation.
The equilibrium coordinate values are supplied as identifiers that are associated with values in the define block.
Unlike using the geovar key, applying the define key does not mean anything in the sense that the various coordinates that refer to the same identifier would be forced to remain equal; it is just a way to display (to the human reader) symmetry in the equilibrium values, to avoid typing errors and to allow an easy adjustment of starting coordinates for another calculation.
Since the atomic coordinates are input in Z-matrix format, the program would by default carry out displacements in internal coordinates to scan the energy surface and hence compute force constants and frequencies. This is overriden by specifying in the geometry block optim cartesian: carry out cartesian displacements.
The key thermo addresses the thermodynamical analysis (only available in a Frequencies calculation, otherwise ignored). The specification 'T=300,400' means that the thermodynamic properties are printed for the temperature range 300-400K, in steps of 10K (default) and for a pressure of 1.0 atmosphere (default).
Frequencies calculations suffer easily from numerical inaccuracies. Therefore, the default numerical integration precision in a Frequencies calculation is much higher than in an ordinary single-point or minimization run. Here we specify the INTEGRATION level to be 5.0 (quite high, but the default for Frequencies is even 6.0).
Rename the TAPE21 result file of the previous calculation so we can restart with other masses. Calculate a different isotope of H, in this case deuterium. It will differ from the original one only in the mass of the nucleus. Repeat the frequency calculation with different fragments. It is important to preserve symmetry at this step so we replace fragment files for ALL H atoms. If you want to replace only one fragment then the original calculation must be performed the same way, with different fragment names.
mv TAPE21 restart.t21 $ADFBIN/adf -n1 <<eor create H M=2.014101779 $ADFRESOURCES/TZP/H eor mv TAPE21 t21.D $ADFBIN/adf <<eor title NH3 frequencies define rNH=1.02 theta=112 phi=120 end atoms Z-matrix XX 0 0 0 N 1 0 0 1.0 H 2 1 0 rNH theta H 2 1 3 rNH theta phi H 2 1 4 rNH theta phi end Fragments N t21.N ! The different isotope mass sits in the next line. H t21.D End geometry optim cartesian frequencies end ! Restart the frequency calculation. ! In fact ADF should perform only one geometry cycle restart restart.t21 thermo T=300,400 integration 5.0 end input eor
Sample directory: adf/Freq_UF6/
Calculation of spin-orbit coupled ZORA gradients for the LDA and GGA functionals
is possible since the ADF2007.01 version.
Summary:
Here only the spin-orbit coupled input file for ADF is given (in the scalar relativistic case change "spinorbit" in "scalar"). The resulting TAPE21 is saved such that it can be used in the frequency calculation.
$ADFBIN/adf <<eor Title UF6 geometry optimization: scalar ZORA integration 5 5 Geometry conv grad=1e-4 End relativistic spinorbit zora Basis Type TZP end ATOMS cartesian 1 U .00000 .00000 .00000 2 F 2.00000 .00000 .00000 3 F -2.00000 .00000 .00000 4 F .00000 2.0000 .00000 5 F .00000 -2.0000 .00000 6 F .00000 .0000 2.00000 7 F .00000 .0000 -2.00000 END end input eor mv TAPE21 UF6.t21
Computation of frequencies by symmetric displacements. The assumed equilibrium input structure should be given in Cartesian coordinates. The calculation starts with the optimized structure read from UF6.t21 (restart file). Again only the spin-orbit coupled input file for ADF is given.
The symmetry is determined automatically by the program as O(H), from the input coordinates. During the calculation first symmetric atomic displacements are constructed. The number of such displacements in each irreducible representation corresponds to the number of frequencies with the corresponding symmetry. All displaced geometries within one representation have the same symmetry, which enables us to use it to speed up the computation significantly. Another advantage of having the same symmetry is that the numerical integration data can be reused efficiently (see SMOOTH option) thus reducing the level of numerical noise in gradients and force constant matrix.
In case of spin-orbit coupling the frequencies can not (yet) be calculated with analytical second derivatives.
$ADFBIN/adf <<eor Title UF6 frequencies and IR intensities: spinorbit ZORA Restart UF6.t21 Geometry Frequencies Symm End integration 5 5 relativistic spinorbit zora Fragments U t21.U F t21.F end ATOMS cartesian 1 U .00000 .00000 .00000 2 F 2.00000 .00000 .00000 3 F -2.00000 .00000 .00000 4 F .00000 2.0000 .00000 5 F .00000 -2.0000 .00000 6 F .00000 .0000 2.00000 7 F .00000 .0000 -2.00000 END end input eor
Sample directory: adf/H2O_HF_freq/
Example shows a Hartree-Fock frequency calculation with an accurate basis set.
The 'DEPENDENCY' key is set to 1e-4. Note that for hybrids and Hartree-Fock the dependency key is always set. The default value in that case is 4e-3. By explicitely setting the 'DEPENDENCY' key we can use a lower value, which is possible in this case. One should check that the results remain reliable if one uses a smaller value for the 'DEPENDENCY' key.
First a geometry optimization is performed.
$ADFBIN/adf << eor
title accurate HF geometry optimization with large QZ4P basis set
basis
type QZ4P
core None
end
adddiffusefit
dependency bas=1e-4
ATOMS
O 0.000000 0.000000 -0.007124
H 0.000000 0.751933 0.556531
H 0.000000 -0.751933 0.556531
END
integration 6 6 6
xc
hartreefock
end
geometry
end
end input
eor
mv TAPE21 t21
Next the frequency calculation is done. A restart is used to pick up the excited state geometry of the previous calculation. The orbitals are not read from the restart file (subkey NoOrb in the restart key). Also the fit coefficients are not read from the restart file (subkey NoSCF in the restart key).
$ADFBIN/adf << eor
title accurate HF frequency calculation with large QZ4P basis set
restart t21 &
noscf
noorb
end
basis
type QZ4P
core None
end
adddiffusefit
dependency bas=1e-4
ATOMS
O 0.000000 0.000000 -0.007124
H 0.000000 0.751933 0.556531
H 0.000000 -0.751933 0.556531
END
integration 6 6 6
xc
hartreefock
end
geometry
frequencies symm
end
end input
eor
Sample directory: adf/EGO_PH2/
Example for an excited state geometry optimization and frequency calculation.
Needed for such excited state optimizations are the key EXCITATIONS (to calculate excitation energies), the key GEOMETRY (to do a geometry optimization) and the key EXCITEDGO (to select for which excitation a geometryy optimization should be performed). The ground state and excited state are open shell.
$ADFBIN/adf <<eor
TITLE PH2 Excited state geometry
atoms
P 0.000000 0.000000 0.0
H 0.7 0.0 0.7
H -0.7 0.0 0.7
end
XC
GGA BP86
END
UNRESTRICTED
CHARGE 0 1
Integration 6.0 6.0 6.0
GEOMETRY
ITERATIONS 50
CONVERGENCE E=0.0001 grad=0.00001
END
SCF
converge 1.0e-9
END
basis
TYPE DZ
CORE NONE
end
excitations
LOWEST 10
onlysinglet
end
EXCITEDGO
STATE B2 1
OUTPUT=1
end
eor
mv TAPE21 PH2.t21
Next the frequencies are calculated of the excited state. A restart is used to pick up the excited state geometry of the previous calculation. Note that in a numerical FREQUENCIES calculation symmetry is turned off except to reduce the number of points calculated. Thus irrespective of the specified point group symmetry the symmetry label A of SYMMETRY NOSYM should be used to select the excited state. Care should be taken to ensure that the correct state is chosen in this frequencies calculation as the excited state number can change when the point group is changed. In this case instead of 'B2 1' one needs to select 'A 1'. Accurate SCF convergence parameters are used.
$ADFBIN/adf <<eor TITLE PH2 Excited state frequencies atoms P 0.000000 0.000000 0.002878 H 1.258230 0.000000 0.655775 H -1.258230 0.000000 0.655775 end XC GGA BP86 END UNRESTRICTED CHARGE 0 1 Integration 6.0 GEOMETRY FREQUENCIES END SCF converge 1.0e-9 END basis TYPE DZ CORE NONE end excitations LOWEST 10 onlysinglet end EXCITEDGO STATE A 1 OUTPUT=2 CPKS EPS=0.000001 end eor
Sample directory: adf/CN_SecDeriv/
The ADF2002.01 version featured analytic second derivatives (SD) for the first time. This initial implementation had severe limitations, both in terms of speed, as in terms of user-friendliness of the output and the number of available options. As of ADF2006.01, the implementation has been significantly improved. More specifically:
Calculation of analytical second derivatives is requested by specifying
AnalyticalFreq End
in the main ADF input.
A high accuracy is specified for the numerical integration to be sure of reliable results. In general, it seems advisable to use high accuracy for heavy nuclei at the moment, whereas default integration accuracy is usually sufficient for light atoms. Further, high integration accuracy is more needed in the atomic spheres than in the rest of the molecule. A cost-effective solution may therefore be to specify a higher integration accuracy in the spheres only (using the accsph subkey of the INTEGRATION keyword).
$ADFBIN/adf << eor title CN atoms N -1.3 0.0 0.0 C 0.0 0.0 0.0 end Basis Type DZ Core None End charge -1 XC LDA Xonly End integration 6.0 AnalyticalFreq End End input eor
After SCF is completed, the energy second derivatives matrix is calculated and analysed, which yields in this case one frequency.
Sample directory: adf/CH4_SecDeriv/
In this example, we use a new feature of adf2006.01: geometry optimization immediately followed by calculation of frequencies. This is done by specifying the Geometry and AnalyticalFreq input blocks in one file.
Note: when using this feature, one should generally set the integration accuracy to a value appropriate for the frequencies calculation, which is about 6.0. In order to save time this recommendation is neglected in this example.
$ADFBIN/adf << eor title CH4 LDA potential Define ZERO = 0.0 RCH = 1.0850 DCH = sqrt(3)*(RCH/3) End Atoms C 0.0 0.0 0.0 H DCH -DCH DCH H DCH DCH -DCH H -DCH DCH DCH H -DCH -DCH -DCH End Basis Type TZP Core None End integration 4.0 Geometry Optim all converge grad=0.0001 End AnalyticalFreq End End input eor
Sample directory: adf/HI_SecDer_ZORA/
The main difference of this example to the previous examples is that a ZORA Hessian is calculated in this example, through the line:
RELATIVISTIC scalar ZORA
Furthermore, the suggestion to use high integration accuracy in the atomic spheres only is shown explicitly here.
$ADFBIN/adf << eor TITLE HI scalar, ZORA, DEFINE R = 1.6090 Z1 = R END ATOMS I 0.0 0.0 0.0 H 0.0 0.0 Z1 END XC LDA Xonly END RELATIVISTIC scalar ZORA Basis Type DZ Core None END integration accint 4.0 accsph 6.0 end AnalyticalFreq End end input eor
Sample directory: adf/MBH_Ethanol/
A frequency calculation is performed using the mobile block Hessian (MBH) method. The coordinates in the ATOMS section should be the partially optimized coordinates (or the fully optimized coordinates would work too). The next input for ADF shows how to perform a frequency calculation with MBH. The flag b=b1 in the ATOMS section adds the label 'b1' to some of the atoms. Only the four atoms labeled 'b1' (CH3) will be considered as a block with fixed internal geometry.
$ADFBIN/adf <<eor TITLE ethanol: second derivatives with MBH approach. CH3 is treated as a rigid block ATOMS 1 C -0.029587 -0.006554 0.008124 b=b1 2 H -0.087498 -0.025163 1.109913 b=b1 3 H 1.027473 -0.056237 -0.302751 b=b1 4 H -0.565305 -0.891154 -0.376242 b=b1 5 C -0.694908 1.238909 -0.501807 b=b2 6 H -0.670258 1.265092 -1.608847 b=b2 7 O -2.069894 1.175059 -0.017251 8 H -0.182335 2.138977 -0.109315 b=b2 9 H -2.586972 1.972802 -0.317216 END SYMMETRY nosym BASIS type DZ core Large CreateOutput None END XC LDA SCF VWN END GEOMETRY frequencies mbh b1 branch new END INTEGRATION 6.0 End input eor
For comparison in this example also a calculation is performed without any restrictions.
$ADFBIN/adf <<eor TITLE ethanol: complete vibrational spectrum, compare with MBH above ATOMS 1 C -0.029587 -0.006554 0.008124 2 H -0.087498 -0.025163 1.109913 3 H 1.027473 -0.056237 -0.302751 4 H -0.565305 -0.891154 -0.376242 5 C -0.694908 1.238909 -0.501807 6 H -0.670258 1.265092 -1.608847 7 O -2.069894 1.175059 -0.017251 8 H -0.182335 2.138977 -0.109315 9 H -2.586972 1.972802 -0.317216 END BASIS type DZ core Large CreateOutput None END XC LDA SCF VWN END AnalyticalFreq End INTEGRATION 5.0 End input eor
Sample directory: adf/MBH_CH4/
A frequency calculation is performed using the mobile block Hessian (MBH)method. The coordinates in the ATOMS section should be the partially optimized coordinates (or the fully optimized coordinates would work too).
Example input how to do a block constraint:
ATOMS
C 0.000000 0.000000 0.000000 b=b1
H 0.634671 0.634671 0.634671 b=b1
H -0.634671 -0.634671 0.634671 b=b1
H -0.634671 0.634671 -0.634671 b=b1
H 0.634671 -0.634671 -0.634671
END
CONSTRAINTS
block b1
END
Such geometry optimization will not be discussed here any further. The next input for ADF shows how to perform a frequency calculation with MBH.
$ADFBIN/adf <<eor
TITLE Methane
BASIS
Type DZ
Core None
END
ATOMS
C 0.000000 0.000000 0.000000 b=b1
H 0.634671 0.634671 0.634671 b=b1
H -0.634671 -0.634671 0.634671 b=b1
H -0.634671 0.634671 -0.634671 b=b1
H 0.634671 -0.634671 -0.634671 b=b2
END
integration 8 8 8
SYMMETRY nosym
GEOMETRY
frequencies disrad=0.001
mbh b1
branch new
END
End input
eor
The flag b=b1 in the ATOMS section adds the label 'b1' to some of the atoms. The four atoms labeled 'b1' will be considered as a block with fixed internal geometry.
In the GEOMETRY section, a Mobile Block Hessian calculation is requested by using the FREQUENCIES and MBH keywords. Here the atoms with label 'b1' are selected to be in the same mobile block. The position/orientation of the block are supposed to be optimized in a preceding partial optimization run. In the vibrational analyis, the block 'b1' is only allowed to vibrate as a whole. The number of resulting modes/frequencies is 3 for the fifth atom plus 6 for the block 'b1' (3 position/3 orientation), resulting in 9 frequencies in total. Since 6 of those frequencies are zero due to translational and rotational invariance of the system, one will find 3 non-zero characteristic frequencies in the output. In practice with ADF not exactly 6 zero's are found, but they are close to zero.
The quality of the frequencies/modes depends largely on the block choice. Best results are obtained when grouping atoms in a block if those atoms are known to form rather rigid structures. For instance, grouping the 11 atoms of benzene side group into a block, will usually result in representative frequencies. In this example the block choice is only illustrative for the methodology.
Sample directory: adf/Freq_NH3_RAMAN/
Summary:
The RamanRange keyword (available since ADF2007.01) can be used to calculate Raman intensities for a range of frequencies only. Using this option is a fast alternative for the existing method of calculating Raman intensities, which is described in the second part of this example.
Two values defining an interval of frequencies to calculate the Raman intensities for. The Raman intensities are calculated by numerical differentiation of the polatizability tensor. Only frequencies frequencies withing the interval that are known to be Raman-active will be included.
$ADFBIN/adf <<eor title NH3 frequencies and calculation of Raman intensities in the range 0-2000 cm-1 atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Basis Type TZP Core Small End AnalyticalFreq end thermo T=300,400 integration 5.0 end input eor mv TAPE21 NH3_freqs.t21 $ADFBIN/adf <<eor title NH3 Raman intensities in the range 0-2000 cm-1 atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Restart NH3_freqs.t21 Fragments H t21.H N t21.N End RamanRange 0.0 2000.0 thermo T=300,400 integration 5.0 end input eor
Raman scattering intensities and depolarization ratios for all molecular vibrations at a certain laser frequency can be calculated in a single run. The run type must be Frequencies and the RESPONSE key is used to specify that Raman intensities are computed.
In this example the static Raman scattering is calculated (ω = 0). This type of calculation is very similar to an IR intensity calculation. In fact, all IR output is automatically generated as well. At all distorted geometries the dipole polarizability tensor is calculated. This is very time-consuming and is only feasible for small molecules.
$ADFBIN/adf <<eor title NH3 frequencies with Raman intensities atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Fragments H t21.H N t21.N End geometry frequencies end response raman end thermo T=300,400 integration 5.0 end input eor
Sample directory: adf/HF_ResonanceRaman/
Example shows a calculation of the Resonance Raman spectrum (RRS) of HF. In this example the RRS is calculated from the geometrical derivatives of the frequency-dependent polarizability, including a finite lifetime.
In the ADF input one then needs to include the subkey FREQUENCIES of the key GEOMETRY (numerical frequencies) and include the subkeys RAMAN and LIFETIME of the key AORESPONSE.
$ADFBIN/adf << eor title HF ao-raman BASIS F DZP/F H DZP/H END GEOMETRY Frequencies END Symmetry NOSYM Atoms H 0.0000 0.0000 0.0000 F 0.0000 0.0000 0.9170 End aoresponse frequency 1 0.52362 Hartree lifetime 0.0034 raman end NOPRINT SFO END INPUT eor
Note that used basis set is too small to get accurate results.
Sample directory: adf/Vibron_RR_uracil/
Example shows a calculation of the Resonance Raman spectrum (RRS) of uracil. In this example the RRS is calculated using the excited-state gradient.
A frequency restart file 'restart.freq' is used as input in the resonance Raman calculation. This restart file is the TAPE21 of a frequency calculation of the runfile 'restart.freq.run'.
First the to ASCII dumped TAPE21 'restart.freq.ascii' is undumped again to make a binary file.
cp $ADFHOME/examples/Test/e_Vibron_RR_uracil/restart.freq.ascii . $ADFBIN/udmpkf < restart.freq.ascii restart.freq
Next the resonance Raman calculation is performed by setting the 'VIBRON' subkey in the 'GEOMETRY' block key, including both the 'EXCITATION' block key and the 'VIBRON' block key. These are the only differences with the frequency run where only the 'FREQUENCIES' subkey was set in the 'GEOMETRY' block key., and the 'EXCITATION' and 'VIBRON' block key were not set.
$ADFBIN/adf << eor Title Input generated by modco EPRINT SFO NOEIG NOOVL NOORBPOP SCF NOPOP END NOPRINT BAS FUNCTIONS UNITS length angstrom angle degree END ATOMS N -0.0147481688 -0.0251586720 0.0000000000 C -0.0263429706 1.3809974655 0.0000000000 N 1.2556533768 1.9305098959 0.0000000000 C 2.5041083561 1.2440596334 0.0000000000 C 2.3755611578 -0.2074475201 0.0000000000 C 1.1446314693 -0.7882184482 0.0000000000 H -0.9346804118 -0.4675883900 0.0000000000 O -1.0845317554 2.0515533614 0.0000000000 H 1.3029888073 2.9549419374 0.0000000000 O 3.5819185026 1.8899458170 0.0000000000 H 3.2859343437 -0.7987226158 0.0000000000 H 0.9976482662 -1.8650665505 0.0000000000 END BASIS type DZ core NONE END XC GGA Becke88 Perdew86 END SYMMETRY tol=0.001 GEOMETRY VIBRON END SCF iterations 50 converge 1.0e-6 1.0e-6 mixing 0.2 lshift 0.0 diis n=10 ok=0.5 cyc=5 cx=5.0 cxx=10.0 END EXCITATION ONLYSING LOWEST 5 END MBLOCKBIG VIBRON NMTAPE restart.freq RESRAMAN STPSIZ 0.1 DOMODES 11 13 16 17 DSCHEME ELCHAR END INTEGRATION 4.0 4.0 END INPUT
Note that used basis set is too small to get accurate results.
Example shows a calculation of normal vibrational Raman optical activity.
In the ADF input one then needs to include the subkey FREQUENCIES of the key GEOMETRY (numerical frequencies) and include the subkey VROA of the key AORESPONSE. A laser frequency need to be added. The other keys in AORESPONSE are recommende by the author of the implementation of VROA in ADF. For accuracy reasons 'INTEGRATION 6' and EXACTDENSITY are used.
$ADFBIN/adf << eor
title VROA TEST H2O2
basis
type TZP
core None
end
ATOMS
1.O -0.750254 -0.034490 0.015133
2.O 0.750254 0.034490 0.015133
3.H -0.943532 0.744006 0.580040
4.H 0.943532 -0.744006 0.580040
END
xc
gga BLYP
end
GEOMETRY
frequencies
END
exactdensity
INTEGRATION 6.0 6.0 6.0
aoresponse
NEWPOLCODE
VROA
scf converge 1d-6 iterations 100
frequency 1 5145 Angstrom
ALDA
FitAOderiv
EL_DIPOLE_EL_DIPOLE VELOCITY
EL_DIPOLE_EL_QUADRUPOLE VELOCITY
EL_DIPOLE_MAG_DIPOLE VELOCITY
end
END INPUT
eor
Note that used basis set is too small to get accurate results. Better is the use larger basis sets, like one of the even tempered basis sets (for example Type ET/ET-QZ3P-1DIFFUSE), or use augmented basis sets (for example Type AUG/ATZ2P).
Sample directory: adf/VROA_RESO/
Example shows a calculation of resonance vibrational Raman optical activity (resonance VROA).
In the ADF input one then needs to include the subkey FREQUENCIES of the key GEOMETRY (numerical frequencies) and include the subkeys VROA and LIFETIME of the key AORESPONSE. A laser frequency need to be added. The other keys in AORESPONSE are recommende by the author of the implementation of VROA in ADF. For accuracy reasons 'INTEGRATION 6' and EXACTDENSITY are used.
$ADFBIN/adf << eor
title VROA TEST H2O2
basis
type TZP
core None
end
ATOMS
1.O -0.750254 -0.034490 0.015133
2.O 0.750254 0.034490 0.015133
3.H -0.943532 0.744006 0.580040
4.H 0.943532 -0.744006 0.580040
END
xc
gga BLYP
end
GEOMETRY
frequencies
END
exactdensity
INTEGRATION 6.0 6.0 6.0
aoresponse
NEWPOLCODE
VROA
scf converge 1d-6 iterations 100
frequency 1 5.15462 eV
lifetime 0.0037
ALDA
FitAOderiv
EL_DIPOLE_EL_DIPOLE VELOCITY
EL_DIPOLE_EL_QUADRUPOLE VELOCITY
EL_DIPOLE_MAG_DIPOLE VELOCITY
end
END INPUT
eor
Note that used basis set is too small to get accurate results. Better is the use larger basis sets, like one of the even tempered basis sets (for example Type ET/ET-QZ3P-1DIFFUSE), or use augmented basis sets (for example Type AUG/ATZ2P).
Sample directory: adf/VCD_COG_NHDT/
Analytical frequencies with subsequent calculation of vibrational circular dichroism (VCD)
The VCD keyword (available since ADF2007.01) can be used to calculate VCD spectra. It is important to note that the VCD keyword only works in combination with the keys AnalyticalFreq and symmetry NOSYM.
Recomended is use to use high accuracy for the geometry optimization which one needs to do before the frequency calculation. This simple example is an NHDT molecule, which is NH3 where one hydrogen atom is replaced with deuterium and another with tritium.
First the atoms are created, next the molecule is calculated.
$ADFBIN/adf -n1 <<eor create H q=1 m=2.014101778 file=$ADFRESOURCES/TZP/H XC gga Becke Perdew end end input eor mv TAPE21 t21.D $ADFBIN/adf -n1 <<eor create H q=1 m=3.01604927 file=$ADFRESOURCES/TZP/H XC gga Becke Perdew end end input eor mv TAPE21 t21.T $ADFBIN/adf -n1 <<eor create N file=$ADFRESOURCES/TZP/N XC gga Becke Perdew end end input eor mv TAPE21 t21.H $ADFBIN/adf -n1 <<eor create N file=$ADFRESOURCES/TZP/N XC gga Becke Perdew end end input eor mv TAPE21 t21.N
Next the molecule is calculated.
$ADFBIN/adf <<eor
Title Single Point calc.
Atoms
N 0.000000 0.000000 0.010272
H -0.471582 -0.816803 0.407861
H.D 0.943163 0.000000 0.407861
H.T -0.471582 0.816803 0.407861
End
Symmetry NOSYM
xc
GGA BP86
end
Fragments
N t21.N
H t21.H
H.D t21.D
H.T t21.T
End
Integration 7.0
AnalyticalFreq
End
VCD
end input
eor
Sample directory: adf/FranckCondon_NO2/
As an example of a Franck-Condon calculation, lets look at the transition of NO2 to NO2-. NO2 is a small molecule with only three vibrational modes. Putting an extra electron on the molecule will cause a big displacement, resulting in large electron-phonon couplings.
In general, the larger the molecule, the smaller the displacement and hence the electron-phonon couplings and Franck-Condon factors. Moreover, larger molecules have more vibrational modes, meaning that the already smaller displacement will generally be smeared out over more modes, resulting in an additional decrease in electron-phonon couplings. This is fortunate, since the number of Franck-Condon factors increases factorially with the number of vibrational modes, making it prohibitively expensive to take more than a few vibrational quanta into account for most molecules.
In order to calculate the Franck-Condon factors for Nitrite and Nitrogen dioxide, the equilibrium positions of the nuclei and the vibrational modes have to be obtained (the geometry optimizations are not shown here):
$ADFBIN/adf << eor
TITLE Nitrogen dioxide
ATOMS
N 0.000000 0.000000 -0.016179
O 0.000000 1.098646 -0.492918
O 0.000000 -1.098646 -0.492918
END
BASIS
CORE NONE
TYPE DZP
END
XC
LDA SCF VWN
END
ANALYTICALFREQ
END
eor
mv TAPE21 NO2.t21
rm t21.* logfile
$ADFBIN/adf << eor
TITLE Nitrite
ATOMS
N 0.000000 0.000000 0.093662
O 0.000000 1.120366 -0.540999
O 0.000000 -1.120366 -0.540999
END
CHARGE -1.0 1.0
UNRESTRICTED
BASIS
CORE NONE
TYPE DZP
END
XC
LDA SCF VWN
END
ANALYTICALFREQ
END
eor
mv TAPE21 NO2-.t21
rm t21.* logfile
This runscript produces two TAPE21 files containing the frequencies and the normal modes for both charge states. Lets first look at the ground state to ground state overlap:
$ADFBIN/fcf << eor
STATES NO2.t21 NO2-.t21
QUANTA 0 0
TRANSLATE
ROTATE
eor
Here, zero vibrational quanta are specified for both charge states, which corresponds to the vibrational ground state. Looking at the standard output, we see for NO2:
| Frequency (cm-1) | λ (dimensionless) |
| 1072.490460 | 1.216108 |
| 1434.990571 | 1.873915 |
| 1875.876562 | 0.000000 |
And for NO2-:
| Frequency (cm-1) | λ (dimensionless) |
| 816.952242 | 0.594643 |
| 1264.390562 | 2.071319 |
| 1314.362101 | 0.000000 |
Both states have two vibrational modes with a significant electron-phonon coupling. The ground state to ground state Franck-Condon factor is therefore expected to be quite small. And indeed, looking at the output, we see that it is 0.7944250686*10-2, less than one percent of the total.
Since NO2 has only three vibrational modes, many quanta can be included, and this indeed turns out to be necessary. Setting the maximum number of quanta at 20 results in 1771 permutations for both states and a total of 3136441 Franck-Condon factors. Even with so many factors, the average sum is still only 0.5196635779. Including one extra vibrational quanta results in an additional 960135 Franck-Condon factors, but an average sum of only 0.5280010614, i.e. less than a percent more. This one percent is smeared out over so many factors that their individual contributions become negligible.
Sample directory: adf/Au2_Resp/
A calculation of response properties of the Au2 dimer, with ZORA relativistic corrections
$ADFBIN/adf << eor Title Au2, Response Properties XC GGA LB94 End Relativistic Scalar ZORA Atoms Au 0.0 0.0 1.236 Au 0.0 0.0 -1.236 End Basis Au ZORA/DZ/Au.4d End Symmetry D(8h) Excitations Lowest 10 TOLERANCE 1d-10 End Response AllComponents End End Input eor
In the response module infinite symmetries cannot be handled (see the User's Guide), so we specify a lower subgroup in the input file, here D(8h).
In this sample run the LB94 potential is used. The implementation implies that the XC potential is evaluated from the exact charge density, rather than the cheaper and faster fitted density (as is the case for other XC functionals). This means that the computation times are longer. In a small molecule like Au2 this hardly shows, but in larger molecules the differences may be more significant. Note that the LB94 is a model potential, thus the calculated bond energy has not so much meaning. If the BASIS key is used and the moldel potential LB94 in the create run of the atoms the BP86 functional is used.
Excitation energies are computed, in principle the lowest 10 in each irrep of the symmetry (see, however, the User's Guide).
Sample directory: adf/CN_unr_exci/
Calculation of the excitation energies of the open shell molecule CN
$ADFBIN/adf << eor Title excitation energies of CN Atoms C .0000 .0000 .0000 N .0000 .0000 1.1718 End unrestricted charge 0 1 excitations lowest 20 end Basis Type AUG/ADZP End End input eor
In this example, the lowest 20 excitation energies of CN are calculated in a spin-unrestricted TDDFT calculation. In the MO → MO transitions part for the excitations of the output file, the spin of each molecular orbitals are also specified to help assign the spin state of the excited states. The transitions are always from α spin-orbital to α spin-orbital or from β spin-orbital to β spin-orbital.
Next the same example for CN is given with the Tamm-Dancoff approximation (TDA) approximation (including TDA in the input). Due to this approximation the calculated excitation energies will not be exactly the same as in the first example.
The third calculation is the calculation of spin-flip excitation energies for CN. Again these energies will not be exactly the same as in the first example. For open-shell molecules, spin-flip transition can result in transition to the ground state with a different Sz value, while the symmetry of the transition density is A1 (Σ+ for linear molecules). The excitation energy of this transition should be zero and this can be used to test the reliability of spin-flip TDDFT. Indeed the calculation of the spin-flip excition energies of CN shows one value which is close to zero and has a transition density of Σ+ symmetry.
$ADFBIN/adf << eor Title spin-flip excitation energies (TDA) of CN Atoms C .0000 .0000 .0000 N .0000 .0000 1.1718 End unrestricted charge 0 1 excitations lowest 20 end SFTDDFT TDA FORCELDA Basis Type AUG/ADZP End End input eor
For accuracy reasons one may need to increase the INTEGRATION accuracy, for example to 6. Note that the basis set is still far from complete, the ATZ2P is better. If one uses the augmented basis sets fFor accuracy reasons one can increase the fit set using a QZ4P Pfit set, for example.
Basis Type AUG/ATZ2P FitType ZORA/QZ4P End INTEGRATION 6
Sample directory: adf/SiH2_spinflip/
Calculation of the spin-flip excitation energies of the open shell molecule SiH2
$ADFBIN/adf << eor Title spin-flip excitation energies of SiH2 Atoms Zmatrix Si 0 0 0 H 1 0 0 1.5145 H 1 2 0 1.5145 92.68 End excitations lowest 20 end unrestricted charge 0 2 SFTDDFT FORCEALDA TDA Basis Type TZ2P End End input eor
In this example, the lowest 20 spin-flip excitation energies of SiH2 are calculated in a spin-unrestricted TDDFT calculation.
In this case an excited state is used as reference, which means that there can also be a negative excitation energy, which is indeed the case. The electron configuration used in the SCF is (a1)1(b1)1, with Sz=1, thus a 3B1 state, which is an excited state. The 1A1 state with electron configuration (a1)2 is lower in energy, and is the ground state.
There is also an excited 1A1 state with electron configuration (b1)2. The transition from the ground 1A1 state to the excited 1A1 state is an excitation from the electron configuration (a1)2 to (b1)2. This transition is actually a double excitation, which means that some double excitations can be reached using spin-flip TDDFT with carefully selected reference states.
In the MO → MO transitions part for the excitations of the output file, the spin of each molecular orbitals are also specified to help assign the spin state of the excited states. Note that in these spin-flip calculations the transitions are always from α spin-orbital to β or from β spin-orbital to α spin-orbital.
For open-shell molecules, spin-flip transition can result in transition to the ground state with a different Sz value, while the symmetry of the transition density is A1. The excitation energy of this transition should be zero and this can be used to test the reliability of spin-flip TDDFT. Indeed the calculation of the spin-flip excition energies of SiH2 shows one value which is close to zero and has a transition density of A1 symmetry.
The 1A1 state with electron configuration (a1)2 can also be used in the calculation of the excitation energies. This is a closed shell configuration, in which case we do not need the spin-flip method.
$ADFBIN/adf << eor Title excitation energies of SiH2 Atoms Zmatrix Si 0 0 0 H 1 0 0 1.5145 H 1 2 0 1.5145 92.68 End excitations lowest 20 end Basis Type TZ2P End End input eor
The transition from the ground 1A1 state to the excited 1A1 state, which is an excitation from the electron configuration (a1)2 to (b1)2, can not be reached in this calculation, since it has mainly double excitation character. Of course, other excited 1A1 states can be reached.
Sample directory: adf/N2_TDHF/
Calculation of the excitation energies of N2 using time-dependent Hartree-Fock (TDHF). It also shows the possibility to use the Tamm-Dancoff approximation (TDA). This examples consists of 4 calculations:
The results will be inaccurate due to small basis set. The key ADDDIFFUSEFIT is required for a more accurate fit of the density.
$ADFBIN/adf << eor Atoms N 0 0 0 N 0 0 1.0977 End XC hartreefock end dependency bas=1e-4 adddiffusefit Basis Type DZ Core None End integration 10 excitations lowest 5 end End Input eor
In case of spin-orbit coupling one needs to include the key RELATIVISTIC with argument ZORA SPINORBIT. In practice one needs to calculate more excitations if one includes spin-orbit coupling, since the singlet-singlet and singlet-triplet excitations are not calculated separately, but will be treated simultaneously, since they may mix.
excitations lowest 20 end relativistic spinorbit zora
The Tamm-Dancoff approximation (TDA) will be used if one includes the key TDA. The calculation is then effectively a CIS calculation.
TDA
Sample directory: adf/TiCl4_CoreExci/
Calculation of the 2p Ti and 2p Cl core excitation energies of TiCl4
The state selective method (key SELECTEXCITATION) can be used to calculate core excitation energies. The use of the key SELECTEXCITATION is similar as the use of the key MODIFYEXCITATION. However, the key SELECTEXCITATION can not be used in combination with the key MODIFYEXCITATION. In the state selective method (key SELECTEXCITATION) the one-electron excited state configuration space remains complete, whereas it is reduced in case the scheme with the MODIFYEXCITATION key.
First an example with the key MODIFYEXCITATION.
$ADFBIN/adf << eor Title TiCl4 TD-DFT scalar relativistic 2p Ti core excitations Units LENGTH BOHR End Atoms Ti 0. 0. 0. Cl 2.36754 2.36754 2.36754 Cl -2.36754 -2.36754 2.36754 Cl 2.36754 -2.36754 -2.36754 Cl -2.36754 2.36754 -2.36754 End SYMMETRY T(D) EPRINT eigval 1000 1000 End XC GGA LB94 End relativistic scalar zora ModifyExcitation UseOccupied T2 2 SubEnd UseScaledZORA END Excitation ONLYSING Davidson & T2 12 SubEnd End Basis Type DZ Core None End end input eor
In this example, the 12 lowest singlet-singlet excitation energies of T2 symmetry are calculated, the dipole allowed excitations. This can also be achieved using the ALLOWED subkey in the key Excitation.
In this example only excitations from the 2t2-orbital are included (see the key MODIFYEXCITATION), an almost pure 2p core orbital of titanium. The orbital energies of the uninteresting other occupied orbitals are artificially changed to a large negative value (-1d6 hartree).
In the second example the 2p Cl core excitation energies of TiCl4 are calculated. The difference between the first example in this one is mainly the use of the key MODIFYEXCITATION:
ModifyExcitation UseOccRange -8.0 -6.0 UseScaledZORA END
In this example only excitations from occupied orbitals are considered which have orbital energies between -8 and -6 hartree, namely the 5a1, 1e, 1t1, 4t2, and 5t2 orbitals, which are almost pure 2p core orbitals of chlorine. The orbital energies of the uninteresting other occupied orbitals are again artificially changed to a large negative value (-1d6 hartree).
Another possibility is the use of the subkey OccVirtRange:
ModifyExcitation UseOccVirtRange 7.0 100.0 UseScaledZORA END
Similarly one can use the key SELECTEXCITATION.
SelectExcitation UseOccupied T2 2 SubEnd UseScaledZORA END
SelectExcitation UseOccVirtRange 7.0 100.0 UseScaledZORA END
Sample directories: adf/Ne_exciso/ and adf/Ne_CoreExci/
Calculation of the (core) excitation energies of Ne including spin-orbit coupling.
The state selective method (key SELECTEXCITATION) can be used to reduce the computational costs of, for example, core excitation energies. In this scheme a guess vector for the orbital transition has to be provided. It should be used in combination with the Davidson method to calculate excitation energies. An overlap criterion is used to follow the wanted eigenvector. This key can also be used in case of spin-orbit coupling. The use of the key SELECTEXCITATION is similar as the use of the key MODIFYEXCITATION. However, the key SELECTEXCITATION can not be used in combination with the key MODIFYEXCITATION. In the state selective method (key SELECTEXCITATION) the one-electron excited state configuration space remains complete, whereas it is reduced in case the scheme with the MODIFYEXCITATION key.
$ADFBIN/adf << eor
Title Ne
Atoms
Ne .0000 .0000 0.0000
End
Basis
Type QZ4P
End
relativistic scalar zora
symmetry d(8h)
integration 6.0
xc
model SAOP
end
Excitations
lowest 10
End
ModifyExcitation
UseOccupied
A1.g 1
SubEnd
UseScaledZORA
END
End input
eor
mv TAPE21 Frag.t21
rm logfile TAPE21
$ADFBIN/adf << eor
Title Ne spin-orbit
Atoms
Ne .0000 .0000 0.0000 f=Frag
End
relativistic spinorbit zora
symmetry d(8h)
xc
model SAOP
end
integration 6.0
Excitations
lowest 12
End
SelectExcitation
UseOccupied
E1/2.g 1
SubEnd
UseScaledZORA
END
Fragments
Frag Frag.t21
End
STCONTRIB
End input
eor
The difference between the core excitation calculation and the standard excitation is the extra subkey MODIFYEXCITATION or SELECTEXCITATION in the core excitation calculation (in italic).
ADF can not handle ATOM and linear symmetries in excitation calculations. Therefore a subsymmetry is used, in this case symmetry d(8h).
A relatively large QZ4P basis set is used, which is still insufficient for excitations to Rydberg-like orbitals, one needs more diffuse functions.
The key STCONTRIB is used, which will give a composition of the spin-orbit coupled excitation in terms of singlet-singlet and singlet-triplet scalar relativistic excitations. In order to use the key STCONTRIB the scalar relativistic fragment should be the complete molecule.
In this case the key MODIFYEXCITATION or SELECTEXCITATION takes care that only excitations from the occupied 1s-orbital (spinor) are included. In symmetry d(8H) the 1s-orbital (spinor) is of A1.g (E1/2.g) symmetry.
Sample directory: adf/AgI_asoexcit/
Calculation of the excitation energies of AgI including spin-orbit coupling in a perturbative way.
$ADFBIN/adf << eor Title [AgI] Atoms Ag .0000 .0000 2.5446 I .0000 .0000 0.0000 End relativistic scalar zora symmetry C(7v) integration 6.0 Excitations lowest 60 End SOPERT Basis Type TZ2P Core None End eor
ADF can not handle ATOM and linear symmetries in excitation calculations. Therefore a subsymmetry is used, in this case symmetry C(7v).
A relatively small TZ2P basis set is used, which is not sufficient for excitations to Rydberg-like orbitals, one needs more diffuse functions.
The key SOPERT is included in scalar relativistic ZORA calculations of excitation energies. First scalar relativistic TDDFT calculations are performed to determine the lowest 60 singlet-singlet and singlet-triplet excited states and the spin-orbit coupling operator is applied to these single-group excited states to obtain the excitation energies with spin-orbit coupling effects included.
Sample directories: adf/EGO_N2/ and adf/EGO_N2_EIGENF/
Two examples for a singlet excited state geometry optimization for N2. In the second example eigenvector following is used.
First example for a singlet excited state geometry optimization, hybrid functional used.
Needed for such excited state optimizations are the key EXCITATIONS (to calculate excitation energies), the key GEOMETRY (to do a geometry optimization) and the key EXCITEDGO (to select for which excitation a geometryy optimization should be performed). In this case a hybrid functional B3LYP is used.
$ADFBIN/adf <<eor
TITLE N2 Excited state geometry
atoms
N 0.0 0.0 -0.7
N 0.0 0.0 0.7
end
XC
HYBRID B3LYP
END
Integration 5
GEOMETRY
ITERATIONS=30
CONVERGE grad=0.0001
END
basis
TYPE DZ
CORE NONE
end
excitations
LOWEST 10
onlysinglet
end
EXCITEDGO
STATE S-.u 1
OUTPUT=2
end
ALLPOINTS
eor
Second example is for a singlet excited state geometry optimization with eigenvector following (subkeyword EIGENFOLLOW of key EXCITEDGO), GGA functional used.
$ADFBIN/adf <<eor
TITLE N2 Eigenvector follow. Lowest state at the starting geometry is not the lowest at min.
atoms
N 0.0 0.0 -0.55
N 0.0 0.0 0.55
end
XC
GGA Becke Perdew
END
Integration 5
GEOMETRY
ITERATIONS=30
CONVERGE grad=0.0001
END
basis
TYPE DZ
CORE NONE
end
excitations
LOWEST 10
onlysinglet
end
EXCITEDGO
STATE A 1
OUTPUT=2
EIGENFOLLOW
end
SYMMETRY NOSYM
eor
Sample directory: adf/EGO_CH2_sf/
Example for a spin-flip excited state geometry optimization with a triplet reference, and a frequency calculation afterwards.
Needed for such excited state optimizations are the key EXCITATIONS (to calculate excitation energies), the key GEOMETRY (to do a geometry optimization) and the key EXCITEDGO (to select for which excitation a geometryy optimization should be performed). In this case spin-flip excitations are calculated.
$ADFBIN/adf <<eor
TITLE CH2 Excited state geometry with triplet reference and spin-flip excitation
atoms
C 0.000000 0.000000 0.0
H 0.7 0.0 0.7
H -0.7 0.0 0.7
end
UNRESTRICTED
CHARGE 0 2
GEOMETRY
ITERATIONS 50
CONVERGENCE E=0.0001 grad=0.0001
END
basis
TYPE DZP
CORE NONE
end
excitations
LOWEST 10
end
TDA
SFTDDFT
FORCEALDA
EXCITEDGO
STATE B2 1
OUTPUT=2
end
SYMMETRY C(2V)
eor
mv TAPE21 CH2.t21
Next the frequencies are calculated of the excited state. A restart is used to pick up the excited state geometry of the previous calculation. Note that in a numerical FREQUENCIES calculation symmetry is turned off except to reduce the number of points calculated. Thus irrespective of the specified point group symmetry the symmetry label A of SYMMETRY NOSYM should be used to select the excited state. Care should be taken to ensure that the correct state is chosen in this frequencies calculation as the excited state number can change when the point group is changed. In this case instead of 'B2 1' one needs to select 'A 2'.
$ADFBIN/adf <<eor
TITLE CH2 Excited state frequencies with triplet reference and spin-flip excitation
atoms
C 0.000000 0.000000 0.0
H 0.7 0.0 0.7
H -0.7 0.0 0.7
end
RESTART CH2.t21
UNRESTRICTED
CHARGE 0 2
GEOMETRY
FREQUENCIES
END
basis
TYPE DZP
CORE NONE
end
excitations
LOWEST 10
end
TDA
SFTDDFT
FORCEALDA
EXCITEDGO
STATE A 2
OUTPUT=2
end
SYMMETRY C(2V)
eor
If the subkey CDSPECTRUM is included in the key EXCITATIONS, the rotatory strength is calculated for the calculated excitations, in order to calculate the CD (Circular Dichroism) spectrum. Only useful for chiral molecules.
With the VELOCITY keyword also the dipole-velocity representation of the rotatory strength is calculated.
Note: results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the CD implementation
Do not use less strict convergence criteria than default, better to use tighter criteria. The approximations in the evaluation of the integrals one makes with the linear scaling techniques are effectively switched off by setting LINEARSCALING 100 (recommended to use this).
Usage:
$ADFBIN/adf <<eor TITLE dimethyloxirane excitations + CD COMMENT results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the CD implementation END XC gga becke perdew END Basis Type DZP Core Small End ATOMS O 0.000129 1.141417 0.000023 C -0.597040 -0.094320 0.428262 C 0.596952 -0.094328 -0.428223 H -0.442927 -0.302650 1.487698 H 0.442944 -0.302474 -1.487720 C -1.978779 -0.386617 -0.093924 H -2.723275 0.220579 0.429114 H -2.043506 -0.157697 -1.159810 H -2.236045 -1.439970 0.055144 C 1.978716 -0.386693 0.093893 H 2.236030 -1.439985 -0.055498 H 2.723156 0.220701 -0.429005 H 2.043497 -0.158088 1.159845 END linearscaling 100 excitations cdspectrum onlysinglet velocity lowest 10 end END INPUT eor
Sample directory: adf/Twist_Ethene_TDDFT/
If the subkey CDSPECTRUM is included in the key EXCITATIONS, the rotatory strength is calculated for the calculated excitations, in order to calculate the CD (Circular Dichroism) spectrum. Only useful for chiral molecules. A hybrid functional is used in this case.
With the VELOCITY keyword also the dipole-velocity representation of the rotatory strength is calculated.
Note: results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the CD implementation
Do not use less strict convergence criteria than default, better to use tighter criteria. The approximations in the evaluation of the integrals one makes with the linear scaling techniques are effectively switched off by setting LINEARSCALING 100 (recommended to use this).
Usage:
$ADFBIN/adf <<eor TITLE twisted ethene COMMENT purpose of this job is to provide a test case for the CD implementation with hybrid functionals, and to test the "velocity" keyword END XC hybrid PBE0 END Basis Type TZ2P Core Small End ATOMS C 0.000000 0.000000 0.000000 C -0.000000 -0.000000 1.350716 H 0.565703 -0.745428 -0.540179 H -0.565703 0.745428 -0.540179 H 0.035837 -0.935093 1.890895 H -0.035837 0.935093 1.890895 END linearscaling 99 excitations cdspectrum onlysinglet velocity lowest 20 end END INPUT eor
Sample directory: adf/H2O_MCD/
Example for the calculation of magnetic circular dichroism (MCD). If the subkey MCD is included in the key EXCITATIONS the MCD parameters of the calculated excitations are calculated (A and B terms). The keys RELATIVISTIC ZORA and SOMCD are required for a calculation of temperature-dependent C terms. The key ALLPOINTS is required for an MCD calculation (if the molecule has symmetry).
$ADFBIN/adf <<eor title water MCD atoms O 0 0 0 H 0 0 1 H 0 1 0 end BASIS TYPE DZP end ALLPOINTS SOMCD UNRESTRICTED CHARGE 1 1 RELATIVISTIC ZORA excitations lowest 20 onlysinglet mcd NMCDTERM=5 end end input eor
Sample directory: adf/H2O_MCD_ZFS/
Example for the calculation of magnetic circular dichroism (MCD) including zero-field splitting (ZFS). If the subkey MCD is included in the key EXCITATIONS the MCD parameters of the calculated excitations are calculated (A and B terms). The keys RELATIVISTIC ZORA and SOMCD are required for a calculation of temperature-dependent C terms. The key ALLPOINTS is required for an MCD calculation (if the molecule has symmetry). For zero-field splitting (ZFS) the electron spin S ≥ 1.
$ADFBIN/adf <<eor title water MCD with zero-field zplitting atoms O 0 0 0 H 0 0 1 H 0 1 0 end BASIS TYPE DZP end ALLPOINTS SOMCD UNRESTRICTED CHARGE 2 2 RELATIVISTIC ZORA ZFS excitations lowest 20 onlysinglet mcd NMCDTERM=5 NODIRECT NTEMP=2 TMIN=1.0 TMAX=300.0 NBFIELD=2 BMIN=1.0 BMAX=8.0 end end input eor
Sample directory: adf/Hyperpol/
This sample illustrates the computation of (hyper) polarizability tensors for the He atom and the H2 molecule.
The symmetry is specified, because the Response module in ADF cannot yet handle the infinite symmetries ATOM, C(lin), D(lin).
$ADFBIN/adf <<EOR Title expt geometrie H2(VII),VWN noprint sfo,frag,functions Symmetry C(8v) Atoms H 0 0 -0.37305 H 0 0 0.37305 End Fragments H t21.H7 End Response HyperPol 0.03 DynaHyp AllComponents End EField 0 0 0.001 end input EOR
The Response data block specifies (AllComponents) that not only the (default) zz-dipole polarizability is to be computed, but the complete tensor. The subkey HyperPol instructs the program to compute hyperpolarizabilities and not only polarizabilities. The DynaHyp subkey implies that the frequency-dependent (hyper)polarizability is calculated. In that case the main laser frequency has to be specified, in hartree units, after the HyperPol subkey.
Only the first hyperpolarizability has been implemented in ADF. Some information on second hyperpolarizabilities can be obtained from the calculation of the first one in a finite field (EFIELD).
In similar fashion the frequency-dependent hyperpolarizability is computed for He, but only the zzz-component because now the AllComponents subkey is omitted.
$ADFBIN/adf <<EOR Title hyperpolarizability He with the LB94 potential noprint sfo,frag,functions Atoms He 0 0 0 End XC GGA LB94 END Fragments He t21.He8 End Response HyperPol 0.07 DynaHyp End integration 5.0 EField 0 0 0.001 Symmetry C(8v) end input EOR
Sample directory: adf/AgI_SO_Pol/
Example shows an frequency-dependent ZORA calculation of complex linear response of AgI including spin-orbit coupling.
The polarizability is calculated using the AORESPONSE key, in which spin-orbit coupling is taken into account. In this case a spin-restricted calculation is required, but, unlike the rest of AORESPONSE, also NOSYM. In this example lifetime effects are included.
$ADFBIN/adf << eor TITLE AgI,SO-ZORA basis type TZ2P core None end ATOMS Ag 0.0 0.0 1.2723 I 0.0 0.0 -1.2723 END symmetry nosym integration 6.0 Relativistic zora spinorbit aoresponse scf frequency 1 0.085 HARTREE lifetime 0.007 ALDA EL_DIPOLE_EL_DIPOLE end END INPUT eor
Sample directory: adf/Disper_HF/
General dispersion coefficients (beyond de dipole-dipole C6 interaction coefficient) are computed with the auxiliary program DISPER. It uses two output files from previous ADF Response calculations. In the example, the two ADF runs are one and the same and the relevant TENSOR output file is used twice.
$ADFBIN/adf <<EOR title Van der Waals coefficients HF atoms H 0 0 -0.8708056087 F 0 0 0.04619439132 end Basis Type DZP Core Small End symmetry C(8v) RESPONSE MAXWAALS 8 VANDERWAALS 7 ALLTENSOR ALLCOMPONENTS END end input EOR
Polarizabilities are computed at 7 (imaginary) frequencies between 0 and infinity. The program determines internally the actual frequency values in this range to use. The user only specifies the number of them, thereby determining the precision of, in fact, a numerical integration over the zero-infinity frequency range. A value of 7 is rather low.
MaxWaals determines that not only the C6 but also C7 and C8 coefficients are computed. A value higher than 8 would not be recommended, because the available basis sets would be inadequate for higher coefficients.
In DISPER calculations the preparatory Response calculation must use the AllTensor and AllComponents subkeys.
The calculation produces a file TENSOR. The subsequent DISPER run uses two such files. In this example, both are taken from the same ADF run, copying the TENSOR file to, respectively, tensorA and tensorB. These names are prescribed for a DISPER calculation.
cp TENSOR tensorA cp TENSOR tensorB $ADFBIN/disper -n1 << eor eor
The DISPER program needs no other input than just the files tensorA and tensorB, which must both be present as local files. Note the '-n1' flag: this enforces that a single-node (non-parallel) run is performed. The current implementation does not support parallelization of DISPER, because the kid processes may not have the (local to the master!) files tensorA and tensorB.
Sample directory: adf/DMO_ORD/
If the subkey OPTICALROTATION is included in the key RESPONSE, the (frequency dependent) optical rotation is calculated.
Note: results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the ORD implementation
Do not use less strict convergence criteria than default, better to use tighter criteria. The approximations in the evaluation of the integrals one makes with the linear scaling techniques are effectively switched off by setting LINEARSCALING 100 (recommended to use this).
Usage:
$ADFBIN/adf <<eor TITLE dimethyloxirane excitations + ORD COMMENT results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the ORD implementation END XC gga becke perdew END Basis Type DZP Core Small End ATOMS O 0.000129 1.141417 0.000023 C -0.597040 -0.094320 0.428262 C 0.596952 -0.094328 -0.428223 H -0.442927 -0.302650 1.487698 H 0.442944 -0.302474 -1.487720 C -1.978779 -0.386617 -0.093924 H -2.723275 0.220579 0.429114 H -2.043506 -0.157697 -1.159810 H -2.236045 -1.439970 0.055144 C 1.978716 -0.386693 0.093893 H 2.236030 -1.439985 -0.055498 H 2.723156 0.220701 -0.429005 H 2.043497 -0.158088 1.159845 END linearscaling 100 response allcomponents opticalrotation end END INPUT eor
Sample directory: adf/DMO_ORD_aoresponse/
If the subkey OPTICALROTATION is included in the key AORESPONSE, the (frequency dependent) optical rotation is calculated. In this example lifetime effects are included. This test example consists of two ORD calculations: one with and one without the velocity gauge.
Note: results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the ORD implementation
Usage:
$ADFBIN/adf <<eor TITLE dimethyloxirane, ORD COMMENT results will be physically meaningless due to small basis set. purpose of this job is to provide a test case for the ORD implementation END XC gga becke perdew END Basis Type DZP Core Small End ATOMS O 0.000129 1.141417 0.000023 C -0.597040 -0.094320 0.428262 C 0.596952 -0.094328 -0.428223 H -0.442927 -0.302650 1.487698 H 0.442944 -0.302474 -1.487720 C -1.978779 -0.386617 -0.093924 H -2.723275 0.220579 0.429114 H -2.043506 -0.157697 -1.159810 H -2.236045 -1.439970 0.055144 C 1.978716 -0.386693 0.093893 H 2.236030 -1.439985 -0.055498 H 2.723156 0.220701 -0.429005 H 2.043497 -0.158088 1.159845 END allpoints aoresponse ALDA opticalrotation frequency 1 5893 angstrom scf iter 20 lifetime 0.007 end END INPUT eor
In the second example the subkey OPTICALROTATION of the key AORESPONSE is changed into VELOCIYORD:
aoresponse ALDA VelocityORD frequency 1 5893 angstrom scf iter 20 lifetime 0.007 end
Sample directory: adf/H2O_Verdet/
Specify the subkey MAGOPTROT in the AORESPONSE key to calculate the Verdet constant.
$ADFBIN/adf <<eor title water basis type TZP core None end atoms O 0.000000 0.134692 0.000000 H 0.869763 -0.538741 0.000000 H -0.869763 -0.538794 0.000000 end symmetry nosym allpoints integration 6.0 linearscaling 99 xc lda vwn gga revPBE end aoresponse scf converge 1d-6 iterations 25 frequency 1 0.088558 Hartree ALDA giao magoptrot end end input eor
Sample directory: adf/DampedVerdet/
Specify the subkey MAGOPTROT in the AORESPONSE key to calculate the Verdet constant. Here it is specified together with the LIFETIME key, such that the real and imaginary part of the damped Verdet constant will be calculated.
$ADFBIN/adf <<eor title Propene ATOMS C 0.867000 1.441800 3.000000 C 0.849400 2.777300 3.000000 C 2.115500 0.591200 3.000000 H -0.088300 0.909000 3.000000 H -0.085900 3.336500 3.000000 H 1.772400 3.363200 3.000000 H 2.737100 0.793300 2.115200 H 1.876900 -0.479100 3.000000 H 2.737100 0.793300 3.884800 END basis type DZP core None end symmetry nosym allpoints INTEGRATION 5.0 XC Model SAOP END noprint sfo aoresponse scf converge 1d-5 iterations 25 frequency 1 0.2 Hartree ALDA giao magoptrot lifetime 0.007 end end input eor
Sample directory: adf/H2O_magnet/
Basic example for a magnetizability calculation.
One should set iterations=0 for STATIC magnetizability. If one does not use SYMMETRY NOSYM, one should set use ALLPOINTS for correct results in AORESPONSE. If a line starts with :: it will be skipped during the reading of the input.
$ADFBIN/adf <<eor title H2O magnetizability test basis type DZP core None end units length bohr angle degree end atoms zmatrix O 0 0 0 0. 0. 0. H 1 0 0 1.808846 0. 0. H 1 2 0 1.808846 104.52 0. end linearscaling 99 tails xc lda gga revPBE end Comment New optiond fro AOResponse below End :: symmetry nosym allpoints :: needed for correct results in AOResponse AOResponse ALDA magneticpert :: needed for magnetizability scf iterations 0 converge 1e-3 :: set iterations=0 for STATIC magnetizability End end input eor
Sample directory: adf/H2O_TD_magnet/
Example for time-dependent magnetizability with GIAOs (Gauge including atomic orbitals).
$ADFBIN/adf <<eor basis type TZP core None end ATOMS O 0.000000 0.134692 0.000000 H 0.869763 -0.538741 0.000000 H -0.869763 -0.538794 0.000000 END symmetry nosym allpoints integration 6.0 linearscaling 99 xc lda vwn gga revPBE end aoresponse scf conv 1d-6 iter 25 frequency 1 5893 Angstrom giao ALDA magneticpert FitAOderiv end END INPUT eor
Sample directory: adf/C2H4_TDCDFT/
Calculation of excitation energies and response properties of C2H4, with the VK functional, thus using time-dependent current-density-functional theory.
$ADFBIN/adf << eor title C2H4 excitation energy calculation with the VK functional ATOMS 1. C 0.000000 0.000000 0.666318 2. C 0.000000 0.000000 -0.666318 3. H 0.000000 0.928431 1.239388 4. H 0.000000 -0.928431 1.239388 5. H 0.000000 0.928431 -1.239388 6. H 0.000000 -0.928431 -1.239388 END BASIS C $ADFHOME/atomicdata/ET/ET-pVQZ/C H $ADFHOME/atomicdata/ET/ET-pVQZ/H END EXCITATIONS END CURRENTRESPONSE END endinput eor
$ADFBIN/adf << eor title C2H4 response calculation with the VK functional ATOMS 1. C 0.000000 0.000000 0.666318 2. C 0.000000 0.000000 -0.666318 3. H 0.000000 0.928431 1.239388 4. H 0.000000 -0.928431 1.239388 5. H 0.000000 0.928431 -1.239388 6. H 0.000000 -0.928431 -1.239388 END BASIS TYPE TZ2P END RESPONSE ALLCOMPONENTS END CURRENTRESPONSE END endinput eor
Sample directories: adf/HBr/ and adf/HBr_SO/
Computation of the NMR chemical shifts for HBr. The second sample uses spin-orbit relativistic corrections.
$ADFBIN/adf << eor TITLE HBr non-relativistic ATOMS 1. H .0000 .0000 .0000 2. Br .0000 .0000 1.4140 End Basis Type DZ Core Large End XC GGA Becke Perdew End End input eor
The TAPE21 result file of ADF must be present under that name for the NMR calculation
mv t21.nmr TAPE21
The NMR program uses only one input (block) key NMR, currently. The subkeys specify what output is produced (OUT) and for which Nuclei the NMR data are computed and printed (NUC). See the User's Guide and the Utilities document for more details.
$ADFBIN/nmr << eor NMR Out TENS Nuc 1 2 End eor
The second run is like the first, except that it uses relativistic corrections, including Spin-Orbit terms. This implies that NOSYM symmetry must be used in the ADF calculation: the NMR program cannot handle symmetry calculations in combination with spin-orbit terms and will stop with an error message if you try to do so.
$ADFBIN/adf << eor TITLE HBr relativistic spinorbit ZORA Atoms 1. H .0000 .0000 .0000 2. Br .0000 .0000 1.4140 End Basis Type DZ Core Large End Symmetry NoSYM XC GGA Becke Perdew End Relativistic SpinOrbit ZORA End input eor rm t12.rel $ADFBIN/nmr << eor NMR U1K BEST OUT TENS NUC 1 2 End eor
Sample directories: adf/HgMeBr_pnr/ (non-relativistic), adf/HgMeBr_psc/ (Pauli scalar relativistic), adf/HgMeBr_zso/ (ZORA relativistic and Spin-Orbit terms included)
NMR data are computed for the 1st and 3rd nucleus only. The UIK subkey is used to indicate that certain terms are to be included in the 'U-matrix', which goes into the first-order change of the MO's due to the applied magnetic field. See the documentation for more information.
The 'BEST' specification means THe best (recommended) options for each relativistic option are included for this sub key. In a non-relativistic run it has no meaning. In a spin-orbit run it would include the ZORA Spin-Orbit terms for a ZORA calculation.
$ADFBIN/nmr << eor NMR U1K BEST NUC 1 3 END eor
The other two calculations are similar, apart from the specification of the applicable relativistic features.
Sample directory: adf/CH4_SAOP/
Computation of the NMR chemical shifts for CH4, with the model potential SAOP.
Important: use SAVE TAPE10. This is necessary for SAOP, since the nmr program does not know about SAOP or other model potentials. On TAPE10 the SCF potential is written, which is read in by the nmr program.
Note: For SAOP one needs an all-electron basis set
$ADFBIN/adf << eor xc model saop end Define RCH = 1.085 XCH = sqrt(3)*(RCH/3) End Atoms C 0 0 0 H XCH XCH XCH H XCH -XCH -XCH H -XCH XCH -XCH H -XCH -XCH XCH End Basis Type TZ2P Core None End save TAPE10 End Input eor $ADFBIN/nmr << eor NMR Out TENS Nuc 1 2 End eor
Sample directory: adf/NMR_NICS/
The NMR program enables the calculation of so-called nucleus-independent chemical shifts (NICS). More details are available in the Properties Programs User's Guide.
In the ADF run, the Efield key is used to define points charges with zero charge. The GHOSTS key in the nmr program then basically copies this block. For the interpretation of the results we refer to the literature.
... Efield 3.0 4.0 5.0 0.0 1.0 2.0 3.0 0.0 End ... eor
Example input for the NMR program.
$ADFBIN/nmr << eor
NMR
Out Iso Tens
GHOSTS
3.0 4.0 5.0
1.0 2.0 3.0
SUBEND
END
END INPUT
eor
Sample directory: adf/NMR_B3LYP/
This example shows how to do hybrid calculation of NMR chemical shifts.
One needs of course a hybrid functional in the XC block key in ADF. One should also SAVE TAPE10, such that it is an input file in the nmr module.
$ADFBIN/adf << eor title PF3-NMR-B3LYP basis type DZP core None end Define RPF = 1.641314 AXPF = 119.702107 End Atoms P 0.00000000 0.00000000 1.00000000 F -0.71283358 1.23466398 1.81325568 F -0.71283358 -1.23466398 1.81325568 F 1.42566716 0.00000000 1.81325568 End integration 6.0 noprint sfo xc hybrid B3LYP end save TAPE10 end input eor rm logfile $ADFBIN/nmr << eor NMR Out TENS Nuc 1 2 SCF 1.d-4 End eor
Next the same calculation is performed with the scalar relativistic ZORA Hamiltonian. In that case one should include in the ADF calculation.
RELATIVISTIC SCALAR ZORA
In the last example spin-orbit coupling is included. Symmetry should be NOSYM.
symmetry nosym RELATIVISTIC SPINORBIT ZORA
In the input for the nmr module one can add the key ZSOAO2007 to approximate the effect of spin on the nucleus in the spin-orbit coupled calculations.
A calculation of NMR nuclear spin-spin coupling constants (NSCCs).
As explained in the ADF manual, the quality of a calculation for spin-spin coupling constants, using the program 'CPL', depends largely on the preceding ADF calculation, which produces the Kohn-Sham orbitals and orbital energies, used as a starting point.
One of the quality-determining factors is the chosen basis set. It should be sufficiently flexible near the nucleus. Although the all-electron basis TZ2P is chosen in this example, it is recommendable to add more functions to the basis and fit sets near the nucleus in case of heavy elements. One could start from a ZORA/QZ basis for example.
The integration accuracy in the ADF calculation is chosen such that the region near the nuclei is described relatively more accurately than the rest of the molecule.
INTEGRATION accint 5 accsph 6 end
The NOSYM symmetry currently needs to be specified in ADF to enable the CPL program to work correctly.
The first call to cpl is as follows:
$ADFBIN/cpl <<eor nmrcoupling dso pso sd scf convergence 1e-7 nuclei 1 2 3 4 nuclei 3 4 end endinput eor
The CPL program can run in parallel.
The specification of what needs to be calculated is given in the nmrcoupling block key.
In this first example, the SD subkey is left out, as this would lead to a very strong increase in the required CPU time. The SD subkey is included in the second CPL run. That subkey controls the calculation of the so-called spin-dipole term.
The subkeys dso and pso specify that, respectively, the diamagnetic and paramagnetic orbital terms will be calculated. The often dominant Fermi contact term (FC) is calculated by default and therefore does not have to be specified explicitly.
The scf convergence subkey, in this context, refers to the convergence for the solution of the coupled-perturbed Kohn-sham equations which need to be solved to obtain to spin-spin couplings.
The lines
nuclei 1 2 3 4 nuclei 3 4
that one coupled-perturbed Kohn-Sham calculation is performed where nucleus number 1 (according to the ordering in the ADF output) is the perturbing nucleus, and nuclei 2, 3, and 4 are the perturbed nuclei, and another coupled-perturbed Kohn-Sham calculation is performed where nucleus 3 is the perturbing nucleus and nucleus 4 is the perturbed nucleus.
The second CPL run also includes the spin-dipole (SD) term, through the SD subkey.
The output of the CPL program first contains a lot of general information, a summary of the specified input, and then produces the desired numbers:
It prints separately the different contributions (FC, DSO, PSO, SD) if specified in input and sums them up to a total number. Experimental NSCCs between two nuclei A and B are usually reported as J(A,B) in Hertz. From a computational point of view, the so-called reduced NSCCs K(A,B) are more convenient for comparisons. CPL outputs both. In this example, the Fermi-contact term is indeed dominant.
The first part of the output refers to the line
nuclei 1 2 3 4
then the same thing is done for the second similar line where nucleus 3 is the perturbing nucleus.
The output for the second CPL run looks very similar, but now the SD term is added to the Fermi contact term, resulting in much longer execution times.
The CPL program also enables calculations using scalar relativistic effects (ZORA) and/or spin-orbit effects.
Schematically, this requires the following changes to the input file with respect to a regular spin-orbit calculation and a nonrelativistic CPL calculation:
Sample directory adf/CPL_HF_hybrid
A calculation of NMR nuclear spin-spin coupling constants (NSCCs) for the hybrid PBE0.
The hybrid PBE0 is chosen as exchange-correlation potential in the ADF calculation. The key 'usespcode' is required for consistency reasons of the PBE0 implementation in ADF and the kernel that is used in the 'CPL' program, that calculates NMR spin-spin coupling constants. Symmetry should be NOSYM. The basis sets used are specially optimized all-electron basis sets for NMR spin-spin coupling calculations (in the directory $ADFHOME/atomicdata/ZORA/jcpl), which have extra tight functions, compared to a default ADF basis set. The integration accuracy is extra high in the core region.
$ADFBIN/adf <<eor
UNITS
length Angstrom
angle Degree
END
:: experimental bond length
ATOMS
F 0.0000 0.0000 0.0000
H 0.0000 0.0000 0.9170
END
BASIS
Type ZORA/jcpl
Core None
END
usespcode
XC
hybrid PBE0
END
SYMMETRY nosym
SCF
converge 1.0e-8 1.0e-6
END
INTEGRATION
accint 6.0
accsph 7.5
end
scf
converge 1e-8 1e-7
end
end input
eor
The first call to cpl is as follows:
$ADFBIN/cpl <<eor gga nmrcoupling dso pso scf convergence 1e-6 iterations 20 nuclei 1 2 end endinput eor
The key 'gga' is included to use the first-order GGA potential instead of the first-order VWN potential. The Hartree-Fock part of the kernel is included automatically if a hybrid potential is used in the ADF calculation.
The second CPL run also includes the spin-dipole (SD) term, through the SD subkey, which is much more time-consuming.
sd
Sample directory: adf/PbH4_finitenuc/
Example for a finite nucleus calculation and the calculation of NMR spin-spin coupling constants.
One of the quality-determining factors for the calculation of NMR coupling constants is the chosen basis set, especially one needs enough tight s functions. If one has a large enough basis set in the core region one can see an effect of using a finite size of the nucleus instead of a point nucleus, especially for heavy nuclei. Such large basis sets can be found for some elements in $ADFRESOURCES/ZORA/jcpl, which are basis sets especially designed for NMR spin-spin coupling calculations. In this example first a basis set for Pb is made which has many tight s functions.
cat < eor > $PWD/Pb
Pb Basis TZ2P - all steep, plus high exponent fcts.
BASIS
1s 39358.2303582207
1s 20080.7297746024
1s 10245.2702931645
1s 5227.17872100229
1s 2666.92791887872
1s 1360.67750963200
1s 694.223219200000
1s 354.195520000000
1s 180.712000000000
2p 39358.2303582207
2p 20080.7297746024
2p 10245.2702931645
2p 5227.17872100229
2p 2666.92791887872
2p 1360.67750963200
2p 694.223219200000
2p 354.195520000000
2p 180.712000000000
1S 92.200
2S 106.000
2S 40.200
3S 51.150
3S 21.850
4S 13.400
4S 9.550
5S 7.650
5S 5.200
6S 3.700
6S 2.400
6S 1.580
2P 196.000
2P 70.200
2P 38.300
3P 21.750
3P 15.900
4P 14.500
4P 9.900
5P 6.600
5P 4.200
6P 2.700
6P 1.600
6P .960
3D 29.400
3D 18.800
4D 13.000
4D 8.450
5D 6.000
5D 3.650
5D 2.200
4F 17.560
4F 9.950
4F 5.750
6D 1.600
5F 2.500
END
CORE 0 0 0 0
END
.....
eor
This large basis set for Pb is used in ADF calculations on PbH4 and the calculation of the NMR spin-spin coupling constants.
$ADFBIN/adf <<eor
UNITS
length Angstrom
angle Degree
END
NuclearModel Gaussian
Print Nuclei
ATOMS
Pb 0.000000 0.000000 0.000000
H -1.023703 1.023703 1.023703
H 1.023703 -1.023703 1.023703
H -1.023703 -1.023703 -1.023703
H 1.023703 1.023703 -1.023703
END
BASIS
type DZP
core None
Pb $PWD/Pb
END
XC
gga PBE
END
usespcode
SCF
iterations 100
converge 1.0e-8 1.0e-6
END
INTEGRATION
accint 5.0
accsph 10.0
end
relativistic scalar zora
end input
eor
$ADFBIN/cpl <<eor
gga
Print Nuclei
nmrcoupling
scf convergence 1e-5 iterations 25
nuclei 1 2
end
endinput
eor
Sample directory: adf/ESR_HfV/
For the ESR g-tensor and D-tensor (zero-field splitting, ZFS) the effect of spin-orbit coupling is important. For the ESR A-tensor and Q-tensor (EFG) spin-orbit coupling is less important.
In this example first spin-orbit coupling is taken into account perturbatively. Next spin-orbit coupling is taken into account self-consistent, using the COLLINEAR keyword.
Note that an all-electron calculation is carried out. This is relevant for the computation of the A-tensor, the nuclear magnetic dipole hyperfine interaction, where an accurate value of the spin-polarization density at the nucleus is important. For the g-tensor this plays a minor role.
In the first ADF calculation the A-tensor (block key ESR) is calculated without the effect of spin-orbit coupling included. The zero-field splitting (key ZFS) is calculated by including spin-orbit coupling perturbatively.
$ADFBIN/adf << eor Atoms Hf 0.0 0.0 0.0 V 0.0 0.0 2.033 End ESR END Unrestricted Symmetry NoSym Charge 0 3 Basis Type TZ2P Core None End ZFS QTENS integration 6 Relativistic Scalar ZORA SAVE TAPE21 TAPE10 end input eor cp TAPE21 hfv.t21 cp TAPE10 hfv.t10
In the next calculation the module nmr calculates the g-tensor (subkey GFACTORS) using spin-orbit coupling and the external magnetic field as perturbation.
$ADFBIN/nmr << eor nmr gfactors u1k best calc all noscl out iso tens end end input
The module cpl can calculate the A-tensor (key HYPERFINE) using spin-orbit coupling and the nuclear magnetic field as perturbation. Note that one needs to set the SCF convergence criterium to a small value.
$ADFBIN/cpl << eor hyperfine atoms 1 2 :: calculates A-tensor for atom 1 and 2, input order SCF Converge 1e-7 end end input eor
ADF can calculate the g-tensor and A-tensor (block key ESR) using only the nuclear or external magnetic field as perturbation, since spin-orbit coupling can be taken into account self-consistently. However, in this case, degenerate perturbation theory is used. The collinear approximation is used (and symmetry NOSYM) to account for spin-polarization effects.
$ADFBIN/adf << eor Atoms Hf 0.0 0.0 0.0 V 0.0 0.0 2.033 End ESR END QTENS collinear Unrestricted Symmetry NoSym Basis Type TZ2P Core None End integration 6 Relativistic Spinorbit ZORA end input
Sample directory: adf/VO_collinear/
The ESR parameters of VO are calculated with the collinear approximation for unrestricted Spin-Orbit coupled calculations. In this example the VO-molecule has three unpaired electrons.
You calculate Electron Spin Resonance properties with the keywords ESR and QTENS. ESR is a block-type key and is used to compute the G-tensor or the Nuclear Magnetic Dipole Hyperfine interaction. QTENS is a simple key and invokes the computation of the Nuclear Electric Quadrupole Hyperfine interaction.
Proper usage of the key ESR requires that you do one of the following:
(a) A Spin-Orbit calculation, spin-restricted, with exactly
one unpaired electron, or
(b) A Spin-Orbit calculation, spin-unrestricted in the collinear approximation, or
(c) No Spin-Orbit terms and spin-unrestricted.
In case (a) and (b) you obtain the G-tensor. In case (b) and (c) you get the Magnetic Dipole Hyperfine interaction.
Note: in case (a) the program also prints a Magnetic Dipole
Hyperfine interaction data, but these have then been computed without the terms
from the spin-density at the nucleus.
Note: in case (b) and (c) one can have more than one unpaired electron.
Note: in case (b) one has to use symmetry NOSYM.
Two calculations are performed:
After the preliminary calculations (DIRAC, to get the relativistic TAPE12 file with relativistic potentials, and the Create runs), we first calculate the Dipole Hyperfine interaction: a spin-unrestricted calculation without Spin-Orbit coupling.
Note that one has to use ALLPOINTS in the calculation for a linear molecule to get results for the nuclear magnetic dipole hyperfine interaction.
$ADFBIN/adf << eor Atoms V 0 0 0 O 0 0 1.589 End XC GGA Becke Perdew End esr end qtens allpoints unrestricted charge 0 3 Relativistic Scalar ZORA CorePotentials t12.rel & V 1 O 2 End integration 5 Fragments V t21.V O t21.O End End input eor
Then a spin-orbit coupled spin-unrestricted calculation is performed using the collinear approximation. Note that symmetry NOSYM is used.
$ADFBIN/adf << eor Atoms V 0 0 0 O 0 0 1.589 End XC GGA Becke Perdew End esr end qtens symmetry nosym unrestricted collinear Relativistic Spinorbit ZORA CorePotentials t12.rel & V 1 O 2 End integration 5 Fragments V t21.V O t21.O End End input eor
Sample directory: adf/ESR_HgF_2der/
This example calculates the ESR g-tensor and A-tensor for HgF. In this example first spin-orbit coupling is taken into account perturbatively.
Note that an all-electron calculation is carried out. This is relevant for the computation of the A-tensor, the nuclear magnetic dipole hyperfine interaction, where an accurate value of the spin-polarization density at the nucleus is important. For the g-tensor this plays a minor role.
In the first example the module nmr calculates the g-tensor (subkey GFACTORS) using spin-orbit coupling and the external magnetic field as perturbation.
$ADFBIN/adf << eor Atoms Hg 0.0 0.0 0.0 F 0.0 0.0 2.804 End Unrestricted Symmetry NoSym Charge 0 1 Basis Type TZ2P Core None End usespcode XC GGA PBE End integration 7 Relativistic Scalar ZORA SCF Iterations 500 Converge 1e-7 1e-7 END SAVE TAPE21 TAPE10 end input eor $ADFBIN/nmr << eor nmr gfactors u1k best calc all noscl out iso tens end end input eor
In the second example the module cpl calculates the A-tensor (key HYPERFINE) using spin-orbit coupling and the nuclear magnetic field as perturbation. Note that one needs to set the SCF convergence criterium to a small value. For an accurate calculation of the A-tensor one needs a very large basis set in the core region (especially tight s-functions), especially for heavy nuclei. If one has such a large basis set in the core region, one can also see an effect of using a finite size of the nucleus instead of a point nucleus. Such large basis sets can be found for some elements in $ADFRESOURCES/ZORA/jcpl, which are basis sets especially designed for ESR A-tensor and NMR spin-spin coupling calculations.
$ADFBIN/adf << eor Atoms Hg 0.0 0.0 0.0 F 0.0 0.0 2.804 End Unrestricted Symmetry NoSym Charge 0 1 Basis Hg ZORA/jcpl/Hg F ZORA/jcpl/F End usespcode XC GGA PBE End integration 7 NuclearModel Gaussian Relativistic Scalar ZORA SCF Iterations 500 Converge 1e-7 1e-7 END end input
$ADFBIN/cpl << eor gga hyperfine atoms 1 2 :: calculates A-tensor for atom 1 and 2, input order SCF Converge 1e-7 end end input eor
Sample directory: adf/ESR_TiF3/
You calculate Electron Spin Resonance properties with the keywords ESR and QTENS. ESR is a block-type key and is used to compute the G-tensor or the Nuclear Magnetic Dipole Hyperfine interaction. QTENS is a simple key and invokes the computation of the Nuclear Electric Quadrupole Hyperfine interaction.
Proper usage of the key ESR requires that you do one of the following:
(a) A Spin-Orbit calculation, spin-restricted, with exactly
one unpaired electron, or
(b) A Spin-Orbit calculation, spin-unrestricted in the collinear approximation, or
(c) No Spin-Orbit terms and spin-unrestricted.
In case (a) and (b) you obtain the G-tensor. In case (b) and (c) you get the Magnetic Dipole Hyperfine interaction.
Note: in case (a) the program also prints a Magnetic Dipole
Hyperfine interaction data, but these have then been computed without the terms
from the spin-density at the nucleus.
Note: in case (b) and (c) one can have more than one unpaired electron.
Note: in case (b) one has to use symmetry NOSYM.
Five calculations are performed:
First a scalar relativistic spin-restricted calculation is performed. The TAPE21 of this calculation is saved as a fragment in the next spin-unrestricted calculation, using only 1 SCF iteration, which is a way to get the scalar relativistic spin-restricted open shell result for the magnetic dipole hyperfine interaction.
$ADFBIN/adf << eor title TiF3 scalar relativistic restricted noprint sfo,frag,functions Define RTIF = 1.780 RY = RTIF*SQRT(3)/2 End Atoms Ti 0 0 0 F RTIF 0 0 F -RTIF/2 RY 0 F -RTIF/2 -RY 0 End Basis Type TZ2P Core None End XC GGA Becke Perdew End relativistic scalar zora End input eor mv TAPE21 t21.TiF3 rm logfile $ADFBIN/adf << eor title TiF3 scalar relativistic open shell restricted noprint sfo,frag,functions Define RTIF = 1.780 RY = RTIF*SQRT(3)/2 End ESR End qtens Atoms Ti 0 0 0 f=TiF3 F RTIF 0 0 f=TiF3 F -RTIF/2 RY 0 f=TiF3 F -RTIF/2 -RY 0 f=TiF3 End Fragments TiF3 t21.TiF3 End XC GGA Becke Perdew End charge 0 1 unrestricted scf iter 0 End relativistic scalar zora End input eor
Next a spin-unrestricted SCF calculation is performed to get the scalar relativistic spin-unrestricted result for the magnetic dipole hyperfine interaction.
$ADFBIN/adf << eor title TiF3 relativistic open shell unrestricted noprint sfo,frag,functions Define RTIF = 1.780 RY = RTIF*SQRT(3)/2 End ESR End qtens Atoms Ti 0 0 0 f=TiF3 F RTIF 0 0 f=TiF3 F -RTIF/2 RY 0 f=TiF3 F -RTIF/2 -RY 0 f=TiF3 End Fragments TiF3 t21.TiF3 End XC GGA Becke Perdew End charge 0 1 unrestricted relativistic scalar zora End input eor
Then, for the same molecule, we compute the G-tensor in a Spin-Orbit run (spin-restricted).
The here-computed and printed Dipole Hyperfine interaction misses the terms from the spin-density at the nucleus: compare with the outcomes from the first calculation.
In each of the calculations, the QTENS key invokes the computation of the Electric Quadrupole Hyperfine interaction.
Note that an all-electron calculation is carried out. This is relevant for the computation of the A-tensor, the nuclear magnetic dipole hyperfine interaction, where an accurate value of the spin-polarization density at the nucleus is important. For the G-tensor (and also for the Q-tensor) this plays a minor role, but for reasons of consistency both calculations use the same basis set and (absence of) frozen core.
$ADFBIN/adf << eor title TiF3 relativistic spinorbit open shell restricted noprint sfo,frag,functions Define RTIF = 1.780 RY = RTIF*SQRT(3)/2 End ESR End qtens Atoms Ti 0 0 0 f=TiF3 F RTIF 0 0 f=TiF3 F -RTIF/2 RY 0 f=TiF3 F -RTIF/2 -RY 0 f=TiF3 End Fragments TiF3 t21.TiF3 End XC GGA Becke Perdew End relativistic spinorbit zora End input eor
Finally a spin-orbit coupled spin-unrestricted calculation is performed using the collinear approximation. Note that symmetry NOSYM is used.
$ADFBIN/adf << eor title TiF3 relativistic spinorbit open shell unrestricted collinear noprint sfo,frag,functions Define RTIF = 1.780 RY = RTIF*SQRT(3)/2 End ESR End qtens symmetry nosym unrestricted collinear Atoms Ti 0 0 0 f=TiF3 F RTIF 0 0 f=TiF3 F -RTIF/2 RY 0 f=TiF3 F -RTIF/2 -RY 0 f=TiF3 End Fragments TiF3 t21.TiF3 End XC GGA Becke Perdew End relativistic spinorbit zora End input eor
The zero-field splitting (ZFS) can be calculated for open shell molecules with electron spin S ≥ 1, using the key ZFS.
Only the spin-orbit contribution to ZFS is evaluated. Can be used in combination with LDA and GGAs. RELATIVISTIC ZORA is also required.
$ADFBIN/adf << eor TITLE NH Zero-field splitting RELATIVISTIC ZORA ATOMS N 0.000000 0.000000 0.007767 H 0.000000 0.000000 -1.043445 END BASIS TYPE DZP END XC GGA BECKE PERDEW END CHARGE 0.0 2.0 UNRESTRICTED ZFS INTEGRATION 5.0 5.0 5.0 ENDINPUT eor
Sample directory: adf/Mossbauer/
By default in ADF the electron density at the nuclei is calculated, no input key is required. The electron density at the nuclei could be used for the interpretation of isomer shifts in Mössbauer spectroscopy. The absolute electron density at a nucleus heavily depends on the accuracy of the basis set in the core region of this nucleus, especially if relativistic effects are included. Important is to use the same basis set, same exchange correlation functional, same integration accuracy, if electron densities at nuclei in different molecules are compared. For the calculation of Mössbauer quadrupole splittings the key QTENS is required. For 57Fe quadrupole splittings will be written in units of mm/s, used in Mössbauer spectroscopy Example shows a calculation on ferrocene with a non-relativistic, and two with a scalar relativistic ZORA Hamiltonian using a different all electron basis set.
First a non-relativistic calculation.
$ADFBIN/adf << eor
title ferrocene
Atoms
FE 0.000000 0.000000 0.000000
C 1.215650 0.000000 1.600813
C 0.375656 -1.156152 1.600813
C -0.983481 -0.714541 1.600813
C -0.983481 0.714541 1.600813
C 0.375656 1.156152 1.600813
C 1.215650 0.000000 -1.600813
C 0.375656 1.156152 -1.600813
C -0.983481 0.714541 -1.600813
C -0.983481 -0.714541 -1.600813
C 0.375656 -1.156152 -1.600813
H 2.310827 0.000000 1.629796
H 0.714085 -2.197727 1.629796
H -1.869498 -1.358270 1.629796
H -1.869498 1.358270 1.629796
H 0.714085 2.197727 1.629796
H 2.310827 0.000000 -1.629796
H 0.714085 2.197727 -1.629796
H -1.869498 1.358270 -1.629796
H -1.869498 -1.358270 -1.629796
H 0.714085 -2.197727 -1.629796
End
xc
gga blyp
end
Basis
Type TZP
Core none
End
qtens
Integration 6
exactdensity
End Input
eor
Next the scalar relativistic ZORA calculations. ADF will also calculate the quadrupole splittings including the small component density, also called SR ZORA-4. The only difference is the RELATIVISTIC keyword:
relativistic scalar zora
Next a scalar relativistic calculation is performed with a much larger basis set (QZ4P) in the core region. Changing the basis set will have a large effect on the electron density at the nucleus and a noticeable effect on the calculated quadrupole splittings.
Basis
Type QZ4P
Core none
End
relativistic scalar zora
Sample directory: adf/AT_transferintegrals/
ADF can can calculate charge transfer integrals, that are needed in approximate methods that model charge transport properties. The molecular system typically should be build from 2 fragments. In this example charge transfer integrals are calculated between Adenine and Thymine. First these two molecules are calculated. In the fragment calculation full symmetry can be used.
$ADFBIN/adf <<eor TITLE Adenine fragment ATOMS 1 N 0.000000000000 0.656191000000 4.473450000000 2 C 0.000000000000 1.850911000000 5.098850000000 3 N 0.000000000000 2.094911000000 6.411070000000 4 C 0.000000000000 0.951291000000 7.115010000000 5 C 0.000000000000 -0.355699000000 6.611740000000 6 C 0.000000000000 -0.487619000000 5.203330000000 7 N 0.000000000000 0.791131000000 8.484350000000 8 C 0.000000000000 -0.567649000000 8.729290000000 9 N 0.000000000000 -1.292469000000 7.631450000000 10 N 0.000000000000 -1.672349000000 4.572610000000 11 H 0.000000000000 2.715551000000 4.433920000000 12 H 0.000000000000 1.540301000000 9.166150000000 13 H 0.000000000000 -0.961519000000 9.739820000000 14 H 0.000000000000 -2.515699000000 5.129900000000 15 H 0.000000000000 -1.718459000000 3.541030000000 END BASIS type DZ core None END eor mv TAPE21 Adenine.t21 $ADFBIN/adf <<eor TITLE Thymine fragment ATOMS 1 N 0.000000000000 0.617991000000 1.666040000000 2 C 0.000000000000 1.851251000000 1.046260000000 3 N 0.000000000000 1.768641000000 -0.347380000000 4 C 0.000000000000 0.582611000000 -1.042160000000 5 C 0.000000000000 -0.621999000000 -0.417040000000 6 C 0.000000000000 -0.627269000000 1.045880000000 7 O 0.000000000000 -1.670479000000 1.720780000000 8 O 0.000000000000 2.924531000000 1.636600000000 9 C 0.000000000000 -1.937039000000 -1.138130000000 10 H 0.000000000000 0.635221000000 2.733380000000 11 H 0.000000000000 2.660141000000 -0.830100000000 12 H 0.000000000000 0.676731000000 -2.127100000000 13 H 0.880180000000 -2.533409000000 -0.860650000000 14 H 0.000000000000 -1.793509000000 -2.225780000000 15 H -0.880180000000 -2.533409000000 -0.860650000000 END BASIS type DZ core None END eor mv TAPE21 Thymine.t21
Next the the base pair is calculated that consists of Adenine and Thymine. To calculate the charge transfer integrals, spatial overlap integrals and site energies, include the key TRANSFERINTEGRALS in the input for ADF. Symmetry NOSYM should be used.
$ADFBIN/adf <<eor
TITLE AT
ATOMS
1 N 0.000000000000 0.656191000000 4.473450000000 f=Adenine
2 C 0.000000000000 1.850911000000 5.098850000000 f=Adenine
3 N 0.000000000000 2.094911000000 6.411070000000 f=Adenine
4 C 0.000000000000 0.951291000000 7.115010000000 f=Adenine
5 C 0.000000000000 -0.355699000000 6.611740000000 f=Adenine
6 C 0.000000000000 -0.487619000000 5.203330000000 f=Adenine
7 N 0.000000000000 0.791131000000 8.484350000000 f=Adenine
8 C 0.000000000000 -0.567649000000 8.729290000000 f=Adenine
9 N 0.000000000000 -1.292469000000 7.631450000000 f=Adenine
10 N 0.000000000000 -1.672349000000 4.572610000000 f=Adenine
11 H 0.000000000000 2.715551000000 4.433920000000 f=Adenine
12 H 0.000000000000 1.540301000000 9.166150000000 f=Adenine
13 H 0.000000000000 -0.961519000000 9.739820000000 f=Adenine
14 H 0.000000000000 -2.515699000000 5.129900000000 f=Adenine
15 H 0.000000000000 -1.718459000000 3.541030000000 f=Adenine
16 N 0.000000000000 0.617991000000 1.666040000000 f=Thymine
17 C 0.000000000000 1.851251000000 1.046260000000 f=Thymine
18 N 0.000000000000 1.768641000000 -0.347380000000 f=Thymine
19 C 0.000000000000 0.582611000000 -1.042160000000 f=Thymine
20 C 0.000000000000 -0.621999000000 -0.417040000000 f=Thymine
21 C 0.000000000000 -0.627269000000 1.045880000000 f=Thymine
22 O 0.000000000000 -1.670479000000 1.720780000000 f=Thymine
23 O 0.000000000000 2.924531000000 1.636600000000 f=Thymine
24 C 0.000000000000 -1.937039000000 -1.138130000000 f=Thymine
25 H 0.000000000000 0.635221000000 2.733380000000 f=Thymine
26 H 0.000000000000 2.660141000000 -0.830100000000 f=Thymine
27 H 0.000000000000 0.676731000000 -2.127100000000 f=Thymine
28 H 0.880180000000 -2.533409000000 -0.860650000000 f=Thymine
29 H 0.000000000000 -1.793509000000 -2.225780000000 f=Thymine
30 H -0.880180000000 -2.533409000000 -0.860650000000 f=Thymine
END
Fragments
Adenine Adenine.t21
Thymine Thymine.t21
end
SYMMETRY NOSYM
TRANSFERINTEGRALS
eor
After the calculation has finished in the output one will find the charge transfer (overlap integrals and site energies) that are needed to calculate hole mobility or electron mobility calculations.
Overlap integral (hole) HOMO fragment 1 - HOMO fragment 2: 0.00000 Site energy (hole) HOMO fragment 1 (eV): -6.89282 Site energy (hole) HOMO fragment 2 (eV): -6.49009 Charge transfer integral (hole) HOMO fragment 1 - HOMO fragment 2 (eV): 0.00000 Overlap integral (electron) LUMO fragment 1 - LUMO fragment 2: 0.00396 Site energy (electron) LUMO fragment 1 (eV): -2.23423 Site energy (electron) LUMO fragment 2 (eV): -2.64264 Charge transfer integral (electron) LUMO fragment 1 - LUMO fragment 2 (eV): -0.04543
Sample directory: adf/green_Al/
As an example of a non-self-consistent Green's function calculation, we will look at the density of states (DOS) and transmission of an infinite 1D chain of Aluminum atoms.
First we need to perform a single-point calculation with ADF on a principal layer, consisting, in this case, of four atoms. Since bulk Aluminum has an FCC structure with a lattice constant of 4.05 Angstrom, the nearest neighbor distance is approximately 2.83 Angstrom. green requires SYMMETRY NOSYM, so we have the following input file for the principal layer:
$ADFBIN/adf << eor
TITLE Principal layer
ATOMS
Al -4.290000 0.000000 0.000000
Al -1.430000 0.000000 0.000000
Al 1.430000 0.000000 0.000000
Al 4.290000 0.000000 0.000000
END
SYMMETRY NOSYM
BASIS
Type DZP
Core Large
CreateOutput None
END
eor
mv TAPE21 layer.t21
The bulk contact geometry consists of three principal layers:
$ADFBIN/adf << eor
TITLE Bulk
ATOMS
Al -15.730000 0.000000 0.000000 f=left
Al -12.870000 0.000000 0.000000 f=left
Al -10.010000 0.000000 0.000000 f=left
Al -7.150000 0.000000 0.000000 f=left
Al -4.290000 0.000000 0.000000 f=center
Al -1.430000 0.000000 0.000000 f=center
Al 1.430000 0.000000 0.000000 f=center
Al 4.290000 0.000000 0.000000 f=center
Al 7.150000 0.000000 0.000000 f=right
Al 10.010000 0.000000 0.000000 f=right
Al 12.870000 0.000000 0.000000 f=right
Al 15.730000 0.000000 0.000000 f=right
END
SYMMETRY NOSYM
FRAGMENTS
left layer.t21
center layer.t21
right layer.t21
END
SCF
Iterations 100
END
eor
mv TAPE21 bulk.t21
Notice that we have increased the number of SCF iterations. The combination of SYMMETRY NOSYM with a 1D chain of metal atoms generally leads to convergence problems. This is the main reason why the principal layer consists of only four atoms. Fortunately, for larger 3D contacts, the convergence is generally better.
From the bulk TAPE21 file green can calculate the self-energies of the left and right contacts. As discussed in the introduction, the self-energy of the left contact needs the center and right fragments of the bulk calculation, and the self-energy of the right contact needs the center and left fragments. Since we need a self-energy matrix for every energy for which we want to calculate the DOS and transmission, already here we have to specify the energy range. We take 1000 points between -0.4 and 0 Hartree.
$ADFBIN/green << eor
SURFACE bulk.t21
FRAGMENTS center right
END
EPS -0.4 0 1000
ETA 1e-6
eor
mv SURFACE left.kf
$ADFBIN/green << eor
SURFACE bulk.t21
FRAGMENTS center left
END
EPS -0.4 0 1000
ETA 1e-6
eor
mv SURFACE right.kf
Since we want to calculate the DOS and transmission of bare aluminum, we can reuse the bulk.t21 file for the extended molecule. We couple the left self-energy to the "left" fragment and the right self-energy to the "right" fragment in bulk.t21. Since we performed restricted ADF calculations, there is no difference between spin-A and spin-B and we can omit spin-B from the calculation.
$ADFBIN/green << eor
DOS bulk.t21
TRANS bulk.t21
EPS -0.4 0 1000
ETA 1e-6
LEFT left.kf
FRAGMENT left
END
RIGHT right.kf
FRAGMENT right
END
NOSAVE DOS_B, TRANS_B
eor
The resulting DOS and transmission are shown in the following figure:
As would be expected for a 1D system, the DOS shows Van Hove singularities at the band edges. Apart from oscillations due to the finite size of the system in ADF, the transmission only reaches integer values. Between approximately -0.35 and -0.15 Hartree, only the sigma channel contributes to the transmission. Above -0.15 Hartree also the two pi channels start to contribute.
Sample directory: adf/green_Al/
run file: green_Au.run
In this example of green, the self-energies are calculated of gold electrodes, the material most often used in molecular electronics. In the example for the Benzenedithiol junction these self-energies will be used to calculate the DOS and transmission of a benzenedithiol junction. The geometry of the electrodes is shown in Fig. 1.
Figure 1: Geometry of the gold contact used in the calculation of the self-energy. The lead consists of two surface layers, left (red) and right (blue), and a bulk layer (green). Each principal layer in turn consists of three atomic layers. This should be sufficient to ensure that the Hamiltonian of the central (green) layer is a bulk Hamiltonian.
Each principal layer contains 3x3x3=27 gold atoms. For the calculation of the self-energies three principal layers are needed, and therefore 81 gold atoms in total. To keep the runtimes manageable it is therefore important to choose the basis set as small as possible. For transport calculations, a DZ basis with a large frozen core is generally sufficient. Unfortunately, even with the largest frozen core (Au.4f), the basis set for Au still contains 19 electrons. A significant speedup can be obtained by limiting this to 11 electrons (only the outer d and s shells). Be advised that even with this reduction the total runtime of calculation can be long.
To facilitate the calculation of the electrodes, first a gold atom fragment will be calculated with the smallest possible basis. The sample directory contains the required Au.5p and Au.5p.dirac files. Note that for gold relativistic effects are important. Therefore RELATIVISTIC Scalar ZORA will be used throughout this example.
$ADFBIN/dirac < Au.5p.dirac
mv TAPE12 t12.rel
$ADFBIN/adf -n1 << eor
CREATE Au file=Au.5p
RELATIVISTIC Scalar ZORA
COREPOTENTIALS t12.rel
XC
LDA SCF VWN
END
eor
mv TAPE21 t21.Au
$ADFBIN/adf << eor
TITLE Gold atom
ATOMS
Au 0.000000 0.000000 0.000000
END
RELATIVISTIC Scalar ZORA
FRAGMENTS
Au t21.Au
END
XC
LDA SCF VWN
END
eor
mv TAPE21 Au.t21
A principal layer of gold consists of three atomic layers, which should be sufficient due to the small screening length. An atomic layer contains 3x3=9 atoms in a (111) FCC configuration. This allows one to use the top-, bride-, and hollow-site binding configurations for a molecule. For the following calculations it is necessary to first construct a fragment of a principal layer.
$ADFBIN/adf << eor
TITLE Principal layer
ATOMS
Au -2.355588 -6.662612 0.000000
Au -2.355589 -4.164133 -1.442498
Au -2.355589 -4.164133 1.442498
Au -2.355589 -1.665653 -2.884996
Au -2.355589 -1.665653 0.000000
Au -2.355589 -1.665653 2.884996
Au -2.355589 0.832826 -1.442498
Au -2.355589 0.832826 1.442498
Au -2.355589 3.331306 0.000000
Au 0.000000 -4.996959 0.000000
Au 0.000000 -2.498480 -1.442498
Au 0.000000 -2.498480 1.442498
Au 0.000000 0.000000 -2.884996
Au 0.000000 0.000000 0.000000
Au 0.000000 0.000000 2.884996
Au 0.000000 2.498480 -1.442498
Au 0.000000 2.498480 1.442498
Au 0.000000 4.996959 0.000000
Au 2.355589 -3.331306 0.000000
Au 2.355589 -0.832826 -1.442498
Au 2.355589 -0.832826 1.442498
Au 2.355589 1.665653 -2.884996
Au 2.355589 1.665653 0.000000
Au 2.355589 1.665653 2.884996
Au 2.355589 4.164133 -1.442498
Au 2.355589 4.164133 1.442498
Au 2.355588 6.662612 0.000000
END
SYMMETRY NOSYM
RELATIVISTIC Scalar ZORA
FRAGMENTS
Au Au.t21
END
XC
LDA SCF VWN
END
eor
mv TAPE21 layer.t21
Three principal layers are stacked together to calculate the self-energies (see Fig. 1). The names of the fragments are significant, since one needs to refer to them by name in the calculation of the self-energies.
$ADFBIN/adf << eor
TITLE Bulk gold
ATOMS
Au -9.422355 -11.659571 0.000000 f=left
Au -9.422356 -9.161092 -1.442498 f=left
Au -9.422356 -9.161092 1.442498 f=left
Au -9.422356 -6.662612 -2.884996 f=left
Au -9.422356 -6.662612 0.000000 f=left
Au -9.422356 -6.662612 2.884996 f=left
Au -9.422356 -4.164133 -1.442498 f=left
Au -9.422356 -4.164133 1.442498 f=left
Au -9.422356 -1.665653 0.000000 f=left
Au -7.066767 -9.993918 0.000000 f=left
Au -7.066767 -7.495439 -1.442498 f=left
Au -7.066767 -7.495439 1.442498 f=left
Au -7.066767 -4.996959 -2.884996 f=left
Au -7.066767 -4.996959 0.000000 f=left
Au -7.066767 -4.996959 2.884996 f=left
Au -7.066767 -2.498479 -1.442498 f=left
Au -7.066767 -2.498479 1.442498 f=left
Au -7.066767 0.000000 0.000000 f=left
Au -4.711178 -8.328265 0.000000 f=left
Au -4.711178 -5.829785 -1.442498 f=left
Au -4.711178 -5.829785 1.442498 f=left
Au -4.711178 -3.331306 -2.884996 f=left
Au -4.711178 -3.331306 0.000000 f=left
Au -4.711178 -3.331306 2.884996 f=left
Au -4.711178 -0.832826 -1.442498 f=left
Au -4.711178 -0.832826 1.442498 f=left
Au -4.711179 1.665653 0.000000 f=left
Au -2.355588 -6.662612 0.000000 f=center
Au -2.355589 -4.164133 -1.442498 f=center
Au -2.355589 -4.164133 1.442498 f=center
Au -2.355589 -1.665653 -2.884996 f=center
Au -2.355589 -1.665653 0.000000 f=center
Au -2.355589 -1.665653 2.884996 f=center
Au -2.355589 0.832826 -1.442498 f=center
Au -2.355589 0.832826 1.442498 f=center
Au -2.355589 3.331306 0.000000 f=center
Au 0.000000 -4.996959 0.000000 f=center
Au 0.000000 -2.498480 -1.442498 f=center
Au 0.000000 -2.498480 1.442498 f=center
Au 0.000000 0.000000 -2.884996 f=center
Au 0.000000 0.000000 0.000000 f=center
Au 0.000000 0.000000 2.884996 f=center
Au 0.000000 2.498480 -1.442498 f=center
Au 0.000000 2.498480 1.442498 f=center
Au 0.000000 4.996959 0.000000 f=center
Au 2.355589 -3.331306 0.000000 f=center
Au 2.355589 -0.832826 -1.442498 f=center
Au 2.355589 -0.832826 1.442498 f=center
Au 2.355589 1.665653 -2.884996 f=center
Au 2.355589 1.665653 0.000000 f=center
Au 2.355589 1.665653 2.884996 f=center
Au 2.355589 4.164133 -1.442498 f=center
Au 2.355589 4.164133 1.442498 f=center
Au 2.355588 6.662612 0.000000 f=center
Au 4.711179 -1.665653 0.000000 f=right
Au 4.711178 0.832826 -1.442498 f=right
Au 4.711178 0.832826 1.442498 f=right
Au 4.711178 3.331306 -2.884996 f=right
Au 4.711178 3.331306 0.000000 f=right
Au 4.711178 3.331306 2.884996 f=right
Au 4.711178 5.829785 -1.442498 f=right
Au 4.711178 5.829785 1.442498 f=right
Au 4.711178 8.328265 0.000000 f=right
Au 7.066767 0.000000 0.000000 f=right
Au 7.066767 2.498479 -1.442498 f=right
Au 7.066767 2.498479 1.442498 f=right
Au 7.066767 4.996959 -2.884996 f=right
Au 7.066767 4.996959 0.000000 f=right
Au 7.066767 4.996959 2.884996 f=right
Au 7.066767 7.495439 -1.442498 f=right
Au 7.066767 7.495439 1.442498 f=right
Au 7.066767 9.993918 0.000000 f=right
Au 9.422356 1.665653 0.000000 f=right
Au 9.422356 4.164133 -1.442498 f=right
Au 9.422356 4.164133 1.442498 f=right
Au 9.422356 6.662612 -2.884996 f=right
Au 9.422356 6.662612 0.000000 f=right
Au 9.422356 6.662612 2.884996 f=right
Au 9.422356 9.161092 -1.442498 f=right
Au 9.422356 9.161092 1.442498 f=right
Au 9.422355 11.659571 0.000000 f=right
END
SYMMETRY NOSYM
RELATIVISTIC Scalar ZORA
FRAGMENTS
left layer.t21
center layer.t21
right layer.t21
END
XC
LDA SCF VWN
END
eor
mv TAPE21 bulk.t21
Similar to the other examples, the self-energies of the left and right contacts is calculated for 1000 energy points between -0.5 and 0 Hartree. This results in two keyfiles of approximately 2.5 GB each. Since the self-energies are independent of whatever is placed between the contacts, they can be reused many times.
$ADFBIN/green << eor
SURFACE bulk.t21
FRAGMENTS center right
END
EPS -0.5 0 1000
ETA 1e-6
eor
mv SURFACE left.kf
$ADFBIN/green << eor
SURFACE bulk.t21
FRAGMENTS center left
END
EPS -0.5 0 1000
ETA 1e-6
eor
mv SURFACE right.kf
In order to interpret transmissions calculated with these self-energies, it is necessary to know the location of the Fermi energy. An estimate for the Fermi energy can be obtained from the bulk SCF calculation by taking the average of the HOMO and LUMO energies, which in this case equals -0.195 Hartree.
Usually the self-energies will be used to calculate the transmission of a molecular junction. However, it is instructive to use a principal layer of gold as the "molecule" and study the DOS and transmission of bulk gold.
$ADFBIN/green << eor
DOS bulk.t21
TRANS bulk.t21
EPS -0.5 0 1000
ETA 1e-6
LEFT left.kf
FRAGMENT left
END
RIGHT right.kf
FRAGMENT right
END
NOSAVE DOS_B, TRANS_B
eor
The results are shown in the following figure:
From this figure it can be seen that around the Fermi energy (-0.2 Hartree), both the DOS and the transmission of gold are relatively constant. This feature makes gold an attractive material for electrodes, since one can expect that the transmission of a molecular junction will be dominated by the molecular properties.
Sample directory: adf/green_Al/
run file: green_BDT.run
In this example of green, the DOS and transmission of a benzenedithiol molecule between gold electrodes is calculated. The calculation uses the self-energies obtained in the example for the gold electrodes. Note that this is a relatively expensive calculation.
First a fragment of the isolated molecule is constructed:
$ADFBIN/adf << eor
TITLE Benzenedithiol
ATOMS
C -1.400000 0.000000 0.000000
C -0.700000 0.000000 -1.200000
C -0.700000 0.000000 1.200000
C 0.700000 0.000000 -1.200000
C 0.700000 0.000000 1.200000
C 1.400000 0.000000 0.000000
H -1.200000 0.000000 -2.200000
H -1.200000 0.000000 2.200000
H 1.200000 0.000000 -2.200000
H 1.200000 0.000000 2.200000
S -3.200000 0.000000 0.000000
S 3.200000 0.000000 0.000000
END
SYMMETRY NOSYM
RELATIVISTIC Scalar ZORA
BASIS
type DZP
core Large
createOutput None
END
XC
LDA SCF VWN
END
eor
mv TAPE21 molecule.t21
Next the molecule is sandwiched between the electrodes in the configuration of Fig. 2. For this the fragment of the principal layer obtained in the example for the gold electrodes is needed.
Figure 2: Geometry of the extended molecule used in the calculation of a benzenedithiol junction. The molecule is shown in green, while the left and right contact regions are shown in red and blue, respectively. Note that the red region corresponds to the blue surface layer in Figure 1 of the example for the gold electrodes and vice versa.
$ADFBIN/adf << eor
TITLE Benzenedithiol
ATOMS
Au -9.911177 -6.662612 0.000000 f=left
Au -9.911178 -4.164133 -1.442498 f=left
Au -9.911178 -4.164133 1.442498 f=left
Au -9.911178 -1.665653 -2.884996 f=left
Au -9.911178 -1.665653 0.000000 f=left
Au -9.911178 -1.665653 2.884996 f=left
Au -9.911178 0.832826 -1.442498 f=left
Au -9.911178 0.832826 1.442498 f=left
Au -9.911178 3.331306 0.000000 f=left
Au -7.555589 -4.996959 0.000000 f=left
Au -7.555589 -2.498480 -1.442498 f=left
Au -7.555589 -2.498480 1.442498 f=left
Au -7.555589 0.000000 -2.884996 f=left
Au -7.555589 0.000000 0.000000 f=left
Au -7.555589 0.000000 2.884996 f=left
Au -7.555589 2.498480 -1.442498 f=left
Au -7.555589 2.498480 1.442498 f=left
Au -7.555589 4.996959 0.000000 f=left
Au -5.200000 -3.331306 0.000000 f=left
Au -5.200000 -0.832826 -1.442498 f=left
Au -5.200000 -0.832826 1.442498 f=left
Au -5.200000 1.665653 -2.884996 f=left
Au -5.200000 1.665653 0.000000 f=left
Au -5.200000 1.665653 2.884996 f=left
Au -5.200000 4.164133 -1.442498 f=left
Au -5.200000 4.164133 1.442498 f=left
Au -5.200001 6.662612 0.000000 f=left
C -1.400000 0.000000 0.000000 f=molecule
C -0.700000 0.000000 -1.200000 f=molecule
C -0.700000 0.000000 1.200000 f=molecule
C 0.700000 0.000000 -1.200000 f=molecule
C 0.700000 0.000000 1.200000 f=molecule
C 1.400000 0.000000 0.000000 f=molecule
H -1.200000 0.000000 -2.200000 f=molecule
H -1.200000 0.000000 2.200000 f=molecule
H 1.200000 0.000000 -2.200000 f=molecule
H 1.200000 0.000000 2.200000 f=molecule
S -3.200000 0.000000 0.000000 f=molecule
S 3.200000 0.000000 0.000000 f=molecule
Au 5.200001 -6.662612 0.000000 f=right
Au 5.200000 -4.164133 -1.442498 f=right
Au 5.200000 -4.164133 1.442498 f=right
Au 5.200000 -1.665653 -2.884996 f=right
Au 5.200000 -1.665653 0.000000 f=right
Au 5.200000 -1.665653 2.884996 f=right
Au 5.200000 0.832826 -1.442498 f=right
Au 5.200000 0.832826 1.442498 f=right
Au 5.200000 3.331306 0.000000 f=right
Au 7.555589 -4.996959 0.000000 f=right
Au 7.555589 -2.498480 -1.442498 f=right
Au 7.555589 -2.498480 1.442498 f=right
Au 7.555589 0.000000 -2.884996 f=right
Au 7.555589 0.000000 0.000000 f=right
Au 7.555589 0.000000 2.884996 f=right
Au 7.555589 2.498480 -1.442498 f=right
Au 7.555589 2.498480 1.442498 f=right
Au 7.555589 4.996959 0.000000 f=right
Au 9.911178 -3.331306 0.000000 f=right
Au 9.911178 -0.832826 -1.442498 f=right
Au 9.911178 -0.832826 1.442498 f=right
Au 9.911178 1.665653 -2.884996 f=right
Au 9.911178 1.665653 0.000000 f=right
Au 9.911178 1.665653 2.884996 f=right
Au 9.911178 4.164133 -1.442498 f=right
Au 9.911178 4.164133 1.442498 f=right
Au 9.911177 6.662612 0.000000 f=right
END
SYMMETRY NOSYM
RELATIVISTIC Scalar ZORA
FRAGMENTS
left layer.t21
molecule molecule.t21
right layer.t21
END
XC
LDA SCF VWN
END
eor
mv TAPE21 fock.t21
The DOS and transmission can now be calculated:
$ADFBIN/green << eor
DOS fock.t21
TRANS fock.t21
EPS -0.5 0 1000
ETA 1e-6
LEFT left.kf
FRAGMENT left
END
RIGHT right.kf
FRAGMENT right
END
NOSAVE DOS_B, TRANS_B
eor
The results are shown in the following figure:
The Fermi energy of the electrodes is -0.2 Hartree (see the example for the gold electrodes). This is just above the HOMO of the molecular junction. Consistent with literature, the HOMO and lower orbitals are combined into a broad peak just below the Fermi energy, while the LUMO is much sharper and situated about 2 eV above the Fermi energy.
The current can be calculated by integrating the transmission around the Fermi energy. At low temperatures this means that the differential conductance is equal to the transmission times the quantum of conductance.
Sample directory: adf/Frags_NiCO4/
An illustration of the fragment feature of adf
A transition metal complex is built from a Nickel atom and four CO fragments. The outcomes allows for an analysis (of molecular orbitals and the Bonding energy) in terms of the fragment properties. It is a Single Point calculation. Geometry optimization would not have been possible in this set-up because an optimization requires that only single-atom fragments are used.
The three atoms are created first: C, O, and Ni. For Carbon and Oxygen a type-DZ basis set is used (double-zeta) using the Basis key, while Ni gets a type-TZP basis (triple-zeta plus polarization).
The CO molecule, to serve as a fragment template in Ni(CO)4, is computed from the atomic fragments C and O. The coordinate values (atoms) are in bohr, rather than in Angstrom because the unit-of-length is redefined by the key units with subkey length.
The key scf is used to specify a somewhat tighter convergence criterion than the default, just to illustrate how to do this (normal settings are quite adequate).
The TAPE21 result file is renamed t21.CO.
$ADFBIN/adf <<eor title CO (as fragment for NiCO4) SCF converge 1e-8 end EPRINT SFO eig ovl END units length bohr end atoms C 0 0 0 O 0 0 2.15617844 end Basis Type DZ Core Small End endinput eor mv TAPE21 t21.CO
Starting from ADF2008.01 one needs to include the subkey SFO of the key EPRINT with arguments eig and ovl in order to get the SFO MO coefficients and SFO overlap matrix printed on standard output.
Apart from the title, the input contains comment. This does not specify computational parameters but is only echoed in the output header, similar to the title. Contrary to the title, however, such comments are not preserved, apart from their echo in output and they are not written to TAPE21 or any other result file.
The atomic coordinates (atoms) are given in bohr (Units). To supply the numerical values use is made of user-defined constants (define): xyzC and xyzOx. This is convenient and it prevents typing errors when several coordinate values are identical, in particular when they carry a lot of decimal places.
The Atoms records contain also a specification of the fragments to which the respective atoms belong: four different CO fragments. No fragment is specified for the Ni atom, which implies that it is a fragment on its own.
The numbers at the very left of the records (1 through 9, with (optionally) a period after them), have no relevance. You can set them for ease of reference or counting.
$ADFBIN/adf <<eor title Ni(CO)4, from fragments Ni and CO COMMENT No geometry optimization possible, because not all fragments are single atoms END units length bohr end EPRINT SFO eig ovl END DEFINE xyzC=2.0053211 xyzOx=3.2501913 END atoms 1. Ni 0 0 0 2. C xyzC xyzC xyzC f=CO/1 3. C -xyzC -xyzC xyzC f=CO/2 4. C xyzC -xyzC -xyzC f=CO/3 5. C -xyzC xyzC -xyzC f=CO/4 6. O xyzOx xyzOx xyzOx f=CO/1 7. O -xyzOx -xyzOx xyzOx f=CO/2 8. O xyzOx -xyzOx -xyzOx f=CO/3 9. O -xyzOx xyzOx -xyzOx f=CO/4 end fragments CO t21.CO Ni t21.Ni end endinput eor
Sample directory: adf/Frags_PtCl4H2
The (scalar) ZORA relativistic option formalism) is used because of the presence of the heavy Pt atom. The complex is built from fragments H2 and PtCl42-.
The calculations of the molecule and larger fragments are performed with GGA's.
The two fragments H2 and PtCl42- are first calculated, from which we are going to build the final complex.
$ADFBIN/adf <<eor Title H2 R=1.68a.u. NoPrint sfo,frag,functions Units length bohr End Atoms H 0.0 0.0 0.84 H 0.0 0.0 -0.84 End Basis Type DZP End XC GGA becke perdew End Relativistic Scalar ZORA End Input eor mv TAPE21 t21H2
The result file TAPE21 is renamed and saved to serve as fragment file.
$adf <<eor title PtCl4 (2-) noprint sfo,frag,functions units length bohr end ATOMS Pt 0 0 0 Cl 4.361580 0.000000 0 Cl 0.000000 4.361580 0 Cl -4.361580 0.000000 0 Cl 0.000000 -4.361580 0 end Basis Type DZP Pt DZ/Pt.4d End xc GGA becke perdew end relativistic scalar ZORA charge -2 end input eor mv TAPE21 t21PtCl4
The key charge is used to specify the net total charge. The default for the net total charge is the sum-of-fragment-charges. The fragments (Pt and Cl atoms) have been computed neutrally, but we want to calculate the PtCl4 complex as a 2- ion.
Finally we compute PtCl4H22- from the fragments PtCl42- and H2.
$ADFBIN/adf <<eor title PtCl4 H2 units length bohr end EPRINT SFO eig ovl END integration 4.0 xc GGA becke perdew end relativistic scalar ZORA ATOMS Pt 0 0 0 f=PtCl4 Cl 0.000000 -4.361580 0.00000000 f=PtCl4 Cl 0.000000 4.361580 0.00000000 f=PtCl4 Cl -4.361580 0.000000 0.00000000 f=PtCl4 Cl 4.361580 0.000000 0.00000000 f=PtCl4 H 0.0 0.0 5.58 f=H2 H 0.0 0.0 7.26 f=H2 end fragments PtCl4 PtCl4.t21 H2 H2.t21 end end input eor
Note that, although the key charge is not supplied, the molecule is not neutral: the default charge (that is, omitting the keys charge, occupations) is the sum-of-fragments: the fragments here are H2 and PtCl42-, yielding a net charge for the molecule of minus two.
Note the f= fragment specification in the Atoms block. No fragment-numbering suffix (/n) is required because there is only one fragment of each fragment type.
Sample directory: adf/UnrFrag_H2
This is a small but important example to illustrate what goes into an accurate calculation of the 'true' bond energy of a molecule. The (ADF-specific) problem is that in a straightforward molecular calculation, the bond energy is computed as the energy difference between at the one hand the molecule, and at the other hand the isolated spherically symmetric spin-restricted atoms. The italic-typed features imply that the reference (comparison) state is usually not the physical ground state of the reference system (isolated atoms) and hence the computed energy difference has no direct relation to experimental data. To account for the true atomic ground states, one has to add correction terms. Study this sample carefully to make sure that you fully understand the steps to take and consult the User's Guide for details. See also the this document for a discussion of multiplet states.
See also the example, SD_Cr(NH3)6.
The H2 case consists of a sequence of simple calculations to demonstrate the Unrestricted Fragments option. The energy difference between an unrestricted fragment as it is used in adf and a self-consistent unrestricted fragment is also computed. This turns out to be quite small, confirming that the adf approach, although not formally exact, is adequate for practical purposes.
$ADFBIN/adf <<eor create H file=$ADFRESOURCES/DZP/H end input eor mv TAPE21 t21H $ADFBIN/adf <<eor title H unrestr., not self-consistent (as used in unr.frag. calcs) scf iterations 0 ! prohibit relaxation end unrestricted charge 0 1 ! if not specified up and down electrons ! will both get 0.5 electron: in fact restricted fragments H t21H end atoms H 0 0 0 end endinput eor rm TAPE21 logfile
By setting the scf iterations to zero (a value of one (1) would give the same result) we prevent cycling to self-consistency. The energy of the 'final' one-electron orbitals is consequently computed in the start-up potential, i.e. the field of the restricted (basic) atom, where spin-α and spin-β are equally occupied, in this case by 0.5 electron each. The not-self-consistent, unrestricted H atom is precisely the 'unrestricted' fragment as it can be used in an adf calculation with unrestricted fragments. The fragment file must be the TAPE21 result file from a restricted run, but at start-up you can specify that the Fragment Orbitals are, for purposes of reference and comparison, occupied in an unrestricted way in the final molecule.
A calculation that uses restricted fragments right away computes the bonding energy relative to the restricted fragments. The difference between using restricted and unrestriced fragments is the 'bonding' energy computed in the run above.
$ADFBIN/adf <<eor title H unr. self-consistent from unr.0 unrestricted charge 0 1 fragoccupations H s 1 // 0 subend end Atoms H 0 0 0 end fragments H t21H end end input eor rm TAPE21 logfile
Here we start with the unrestricted fragment and relax to self-consistency. The 'bonding energy', i.e. the relaxation energy, is very small, demonstrating that using non-self-consistent unrestricted fragments involves only a small error (which, moreover, can be computed as shown here).
The key UnRestricted sets the spin-unrestricted mode. The key Charge is used to specify a net total charge of zero and a net total spin polarization by an excess of 1.0 spin-α electrons over spin-β.
$ADFBIN/adf <<eor title H2 restricted, from restricted fragments ATOMS H 0 0 0.375 H 0 0 -0.375 end fragments H t21H end end input eor rm TAPE21 logfile
This is the simplest approach, using restricted fragments. The bonding energy must be corrected because the reference (restricted H atoms, with 0.5 electrons in spin-α and 0.5 in spin-β) is far from the true H-atom ground state: see the previous runs on the single H atom.
$ADFBIN/adf <<eor
title H2 from unrestricted fragments
ATOMS
H.1 0 0 0.375
H.2 0 0 -0.375
end
fragments ! two different fragment types are necessary
! because the two atoms get different FragOccupations
! (see below), while the key FragOc.. addresses
! only fragmentTYPES
H.1 t21H
H.2 t21H
end
charge 0
occupations
sigma 2 ! specify the state (not always
! necessary)
end
fragoccupations
H.1
s 1 // 0
subend
H.2
s 0 // 1
subend
end
modifystartpotential
H.1 1 // 0 ! this helps SCF start-up
H.2 0 // 1 ! but is here not necessary
end
end input
eor
rm TAPE21 logfile
This should be a fair approximation (in the lda model) to the bonding energy of H2 with respect to the unrestricted H atoms. The difference between the bonding energies of this and the previous run should be very close to the energy of the not-self-consistent unrestricted H-atom with respect to the restricted basic atom (calculation #2).
$ADFBIN/adf <<eor title H2 excited ATOMS H 0 .0 0.375 H 0 .0 -0.375 end fragments H t21H end fragoccupations H s 1 // 0 subend end unrestricted charge 0 2 occupations sigma.g 1 // 0 sigma.u 1 // 0 end end input eor
Finally the calculation of an excited state, with respect to unrestricted fragments. The excitation energy is obtained by comparing the energy with the energy of the ground state calculation. This difference compares reasonably, but not accurately, to the difference in one-electron ground state energies of the involved orbitals (Koopman's theorem).
Note that excitation energies can also be calculated with Time-Dependent DFT, using the RESPONSE module of ADF. See related sample runs.
Sample directory: adf/PCCP_Unr_BondEnergy/
This example illustrates advanced usage of the bond energy decomposition scheme used in ADF.
A proper decomposition of an electron-pair bond energy requires specifying opposite spins for the unpaired electrons of the respective radical fragments, which can be done with the input key FragOccupations. The specified alpha- and beta-spin configurations of the radical fragments are shown in the output section B U I L D.
Please note that if one neglects explicitly specifying opposite spins for the unpaired electrons of the fragments, each of them is treated as being half an alpha and half a beta electron and consequently, they enter into a spurious Pauli repulsive interaction. This results, amongst others, into the Pauli repulsion term being too repulsive and the orbital interaction term being too much stabilizing.
The example consists of an analysis of the C-C single bond between two CP radicals in the four-atomic molecule PCCP. The CP fragment calculations used to provide the TAPE21 for the overall PCCP calculation must be done, for technical reasons, in the restricted mode ("cp_fpccp_asr"). The proper spins are then specified in the calculation of the overall molecule using the FragOccupations key ("pccp_fa1_as"). Note that this implies a slight approximation because the bond energy computed in this way refers to the energy difference between closed-shell PCCP and two CP radicals that are described by orbitals from a spin-restricted SCF calculation, which have been given an unrestricted occupation. In other words, the set of alpha- and beta-spin orbitals are identical and the effect of spin polarization is missing. In practice, this leads to minor energy differences with respect to the correct bond energy, that is, the energy difference between closed-shell PCCP and two CP radicals treated in the unrestricted mode, i.e., for which the set of alpha- and beta-spin orbitals are allowed to relax toward different solutions in the SCF procedure. This correction term can be computed directly by carrying out
An unrestricted computation of the CP radical ("cp_fpccp_asu") using the restricted CP radical ("cp_fpccp_asr") as a fragment.
$ADFBIN/adf<<eor TITLE cp_fpccp_asr EPRINT SFO eig ovl END XC GRADIENTS BECKE PERDEW END ATOMS C .0000 .0000 .6681 P .0000 .0000 2.2555 END BASIS Type TZ2P Core Small END integration 5.0 END INPUT eor mv TAPE21 t21cp_fpccp $ADFBIN/adf<<eor TITLE cp_fpccp_asu EPRINT SFO eig ovl END XC GRADIENTS BECKE PERDEW END ATOMS C .0000 .0000 .6681 f=CP P .0000 .0000 2.2555 f=CP END FRAGMENTS CP t21cp_fpccp END UNRESTRICTED CHARGE 0 1 integration 5.0 END INPUT eor rm TAPE21 logfile $ADFBIN/adf<<eor TITLE pccp_fa1_as EPRINT SFO eig ovl ORBPOP 20 20 SUBEND END XC GRADIENTS BECKE PERDEW END ATOMS P .0000 .0000 2.2555 f=CP_A C .0000 .0000 .6681 f=CP_A C .0000 .0000 -.6681 f=CP_B P .0000 .0000 -2.2555 f=CP_B END integration 5.0 FRAGMENTS CP_A t21cp_fpccp CP_B t21cp_fpccp END SYMMETRY C(LIN) FRAGOCCUPATIONS CP_A SIGMA 3//2 PI 2//2 SUBEND CP_B SIGMA 2//3 PI 2//2 SUBEND END END INPUT eor
Sample directory: adf/EDA_meta_gga_hybrid/
This example illustrates the bond energy decomposition scheme using metaGGA or metahybrid or hybrid functionals in ADF.
The first example is straightforward with closed shell atomic fragments: Zn2. The second example has open shell atomic fragments: Cr2, and the extra complication that spin symmetry breaking lowers the energy of the molecule, although the total Sz-value is zero. The third example has open shell atomic fragments, and the molecule is open shell: CrH.
In the first example for Zn2 the metahybrid TPSSh is used. In the bond energy analysis, the bond energy is split in a Pauli repulsion term, a steric interaction, and an orbital interaction.
$ADFBIN/adf<<eor Atoms Zn 0.0 0.0 0.0 Zn 0.0 0.0 3.2 End XC metahybrid TPSSh end Basis Type TZ2P Core None End dependency bas=1e-4 integration 6 End Input eor
In the second example for Cr2 the metaGGA TPSS is used. Since the fragments are open shell, one may want to use unrestricetd fragments, however, this is not possible in ADF. A fair approximation to a computation with unrestricted fragments can be achieved with the key FRAGOCCUPATIONS. You tell ADF that you want to treat the fragments as if they were unrestricted; this causes the program to duplicate the one-electron orbitals of the fragment: one set for spin-α and one set for spin-β. You can then specify occupation numbers for these spin-unrestricted fragments, and occupy spin-α orbitals differently from spin-β orbitals. Especially for the Pauli-repulsion it is important that one chooses the spin-occupations on the different fragments such that they are 'prepared for bonding'.
Of course, the unrestricted fragments that you use in this way, are not self-consistent: different numbers of spin-α and spin-β electrons usually result in different spatial orbitals and different energy eigenvalues for spin-α and spin-β when you go to self-consistency, while here you have spatially identical fragment orbitals. Nevertheless it is often a fair approximation which gives you a considerable extension of analysis possibilities.
Spin-symmetry breaking is enforced by the use of the key ModifyStartPotential in combination with the key key UNRESTRICTED. In the ADF output one can find that there is spin-density on both of the atoms.
$ADFBIN/adf<<eor Atoms Cr.1 0.0 0.0 0.0 Cr.2 0.0 0.0 1.8 End XC metagga TPSS end Basis Type TZ2P Core None End dependency bas=1e-4 integration 6 unrestricted charge 0 0 ModifyStartPotential Cr.1 1 // 0 Cr.2 0 // 1 End FragOccupations Cr.1 S 4 // 3 P 6 // 6 D 5 // 0 SubEnd Cr.2 S 3 // 4 P 6 // 6 D 0 // 5 SubEnd End End Input eor
In order to calculate the effect of self-consistency one should calculate the Cr atom spin-unrestrictedly.
$ADFBIN/adf<<eor Atoms Cr 0.0 0.0 0.0 End XC metagga TPSS end Basis Type TZ2P Core None End integration 6 unrestricted charge 0 6 FragOccupations Cr S 4 // 3 P 6 // 6 D 5 // 0 SubEnd End End Input eor
In the third example for CrH the hybrid B3LYP is used.
$ADFBIN/adf<<eor Atoms Cr 0.0 0.0 0.0 H 0.0 0.0 1.65 End XC hybrid B3LYP end Basis Type TZ2P Core None End dependency bas=1e-4 integration 6 unrestricted charge 0 5 FragOccupations Cr S 4 // 3 P 6 // 6 D 5 // 0 SubEnd H S 0 // 1 SubEnd End End Input eor
In order to calculate the effect of self-consistency of spin-polarization on the atoms one should calculate the Cr and H atom spin-unrestrictedly.
$ADFBIN/adf<<eor Atoms Cr 0.0 0.0 0.0 End XC hybrid B3LYP end Basis Type TZ2P Core None End dependency bas=1e-4 integration 6 unrestricted charge 0 6 FragOccupations Cr S 4 // 3 P 6 // 6 D 5 // 0 SubEnd End End Input eor rm TAPE21 logfile $ADFBIN/adf<<eor Atoms H 0.0 0.0 0.0 End XC hybrid B3LYP end Basis Type TZ2P Core None End dependency bas=1e-4 integration 6 unrestricted charge 0 1 FragOccupations H S 1 // 0 SubEnd End End Input eor
Sample directory: adf/TlH_SO_analysis/
Application of the Spin-Orbit relativistic option (using double-group symmetry) to TlH with a detailed analysis of the spinors in terms of SFOs (Symmetrized Fragment Orbitals).
In order to get the population analysis, one should have one scalar relativistic fragment, which is the whole molecule. The SFOs in this case are the scalar relativistic orbitals, which are already orthonormal, because one has only one fragment which is the whole molecule.
First the relativistic fragment is made, including the create of the atoms:
$ADFBIN/dirac -n1 < $ADFRESOURCES/Dirac/Tl $ADFBIN/dirac -n1 < $ADFRESOURCES/Dirac/H mv TAPE12 t12.rel $ADFBIN/adf <<eor create Tl file=$ADFRESOURCES/ZORA/TZ2P/Tl xc LDA vwn GGA becke perdew end relativistic scalar zora corepotentials t12.rel & Tl 1 H 2 end end input eor mv TAPE21 t21.Tl $ADFBIN/adf <<eor create H file=$ADFRESOURCES/ZORA/TZ2P/H xc LDA vwn GGA becke perdew end relativistic scalar zora corepotentials t12.rel & Tl 1 H 2 end end input eor mv TAPE21 t21.H $ADFBIN/adf <<eor title TlH, scalar relativistic zora integration 6.0 relativistic scalar zora corepotentials t12.rel & Tl 1 H 2 end ATOMS Tl 0.0 0.0 0.0 H 0.0 0.0 1.870 end fragments Tl t21.Tl H t21.H end xc LDA vwn GGA becke perdew end EPRINT SFO eig ovl END end input eor mv TAPE21 t21.TlH
In order to get the population analysis, one should have one scalar relativistic fragment, which is the whole molecule, which is TlH in this case.
$ADFBIN/adf <<eor title TlH from fragment TlH, with SpinOrbit coupling integration 6.0 relativistic spinorbit zora corepotentials t12.rel & Tl 1 H 2 end ATOMS Tl 0.0 0.0 0.0 f=TlH H 0.0 0.0 1.870 f=TlH end fragments TlH t21.TlH end xc LDA vwn GGA becke perdew end EPRINT SFO eig ovl END end input eor
The output gives something like:
=======================
Double group symmetry : *** J1/2 ***
=======================
=== J1/2:1 ===
Spinors expanded in SFOs
....
Spinor: 21 22 23 24
occup: 1.00 1.00 1.00 0.00
------ ---- ---- ---- ----
SFO SIGMA
13.alpha: 0.7614+0.0000i 0.0096+0.0000i 0.0052+0.0000i -0.0006+0.0000i
14.alpha: 0.0154+0.0000i -0.9996+0.0000i 0.0208+0.0000i -0.0077+0.0000i
15.alpha: -0.0146+0.0000i 0.0185+0.0000i 0.9849+0.0000i 0.1625+0.0000i
SFO PI:x
8.beta : 0.4578+0.0000i 0.0091+0.0000i 0.0112+0.0000i 0.0030+0.0000i
9.beta : 0.0005+0.0000i -0.0074+0.0000i -0.1119+0.0000i 0.6910+0.0000i
SFO PI:y
8.beta : 0.0000+0.4578i 0.0000+0.0091i 0.0000+0.0112i 0.0000+0.0030i
9.beta : 0.0000+0.0005i 0.0000-0.0074i 0.0000-0.1119i 0.0000+0.6910i
....
Left out are a lot of small numbers. The meaning is that a spinor
of J_z=1/2 symmetry can have SIGMA and PI character, for example,
the 21st spinor with occupation number 1.0, is approximately
(21 J_z=1/2) = 0.76 (13 SIGMA alpha) + 0.46 (8 PI:x beta) + i 0.46 (8 PI:y beta)
Next in the SFO contributions per spinor the real and imaginary spin alpha part and real and imaginary spin beta part are all summed together to give a percentage of a certain SFO. are summed. For example the 21st spinor has almost 60% (13 SIGMA) character.
SFO contributions (%) per spinor
Spinor: 21 22 23 24
occup: 1.00 1.00 1.00 0.00
------ ---- ---- ---- ----
SFO SIGMA
13: 57.97 0.01 0.00 0.00
14: 0.02 99.92 0.04 0.01
15: 0.02 0.03 97.01 2.64
SFO PI:x
8: 20.96 0.01 0.01 0.00
9: 0.00 0.01 1.25 47.75
SFO PI:y
8: 20.96 0.01 0.01 0.00
9: 0.00 0.01 1.25 47.75
Starting from the ADF2008.01 version in ADF one calculate Bader atomic charges using a grid based method. This is described in this example. Another possibility for Bader's analysis is to use the adf2aim utility such that a third party program Xaim can be used.
With the BADER input key the ADF program will calculate Bader charges (AIM charges) using a grid based method.
$ADFBIN/adf <<eor
Title Calculate Bader analysis for water
Atoms
O 0.000000 0.000000 -0.001658
H -0.769048 0.000000 0.595209
H 0.769048 0.000000 0.595209
End
Basis
Type TZP
Core none
End
Bader
End Input
eor
Next a similar calculation for ferrocene is given, which is not repeated here.
Sample directory: adf/BondOrder/
With the key BONDORDER a bond order analysis is performed based on SFOs. Note that symmetry used in the calculation should be NOSYM. Shown here is only the example for benzene, where the bond orders calculated are with respect to the atomic fragments.
$ADFBIN/adf <<eor title benzene BP/SZ bondorders tol=0.05 define cc=1.38476576 ccc=120.0 dih=0.0 hc=1.07212846 hcc= 120.0 dih2=180.0 end atoms Z-matrix C 0 0 0 C 1 0 0 cc C 2 1 0 cc ccc C 3 2 1 cc ccc dih C 4 3 2 cc ccc dih C 5 4 3 cc ccc dih H 2 1 3 hc hcc dih2 H 3 2 4 hc hcc dih2 H 4 3 5 hc hcc dih2 H 5 4 3 hc hcc dih2 H 6 5 4 hc hcc dih2 H 1 2 3 hc hcc dih2 end basis Type SZ Core None end symmetry NOSYM xc gga becke perdew end bondorder tol=0.05 printall noprint sfo eor
Sample directories: adf/Diimina_NOCV/ and adf/Hplus_CO_etsnocv
Example for calculation of ETS-NOCV for spin-restricted fragments. ETS-NOCV: energy analysis using the Natural Orbitals for Chemical Valence. The ethylene molecule and a Ni-diimina form a complex together. This example will be discussed first. The other example is H+ and CO form together HCO+, this example is similar to the discussed example. All electron basis sets are required.
First the two fragments are calculated.
$ADFBIN/adf << eor Title: et-----Ni-diimina: ethylene run atoms cartesian C -0.430177075 -1.815433265 0.860288229 C -0.363705637 -1.910722338 -0.515633302 H 0.533109934 -2.284970854 -1.016904201 H -1.279922499 -1.884673940 -1.115144723 H -1.389295819 -1.753589602 1.377541080 H 0.440296224 -2.041861443 1.484489314 end basis Type DZP Core Small end symmetry NOSYM xc gga scf becke perdew end endinput eor mv TAPE21 t21.etfrag $ADFBIN/adf << eor Title: et-----Ni-diimina: Ni-diimina run atoms cartesian Ni 0.022615419 0.037783871 0.025751533 N 0.386170317 1.871072585 0.306265538 C 1.612863056 2.248007643 0.148716016 C 2.540686607 1.163409862 -0.183603690 N 1.976290003 0.008161589 -0.301176178 H -0.288333328 2.609667211 0.546869047 H 1.942601454 3.283060847 0.269249237 H 3.613259273 1.338293482 -0.302134814 H 2.621707427 -0.766258151 -0.517479818 H -1.351756655 0.253389698 0.386197419 end charge 1 basis Type DZP Core Small end symmetry NOSYM xc gga scf becke perdew end endinput eor mv TAPE21 t21.Nifrag
Next these fragments are used in the calculation of the full complex. The keys ETSNOCV and 'PRINT etslowdin' are needed in this case to to analyze the bonding in the molecule with respect to the fragments. The symmetry must be NOSYM.
$ADFBIN/adf << eor Title: et-----Ni-diimina run atoms Ni 0.022615419 0.037783871 0.025751533 f=k N 0.386170317 1.871072585 0.306265538 f=k C 1.612863056 2.248007643 0.148716016 f=k C 2.540686607 1.163409862 -0.183603690 f=k N 1.976290003 0.008161589 -0.301176178 f=k H -0.288333328 2.609667211 0.546869047 f=k H 1.942601454 3.283060847 0.269249237 f=k H 3.613259273 1.338293482 -0.302134814 f=k H 2.621707427 -0.766258151 -0.517479818 f=k H -1.351756655 0.253389698 0.386197419 f=k C -0.430177075 -1.815433265 0.860288229 f=m C -0.363705637 -1.910722338 -0.515633302 f=m H 0.533109934 -2.284970854 -1.016904201 f=m H -1.279922499 -1.884673940 -1.115144723 f=m H -1.389295819 -1.753589602 1.377541080 f=m H 0.440296224 -2.041861443 1.484489314 f=m end charge 1 fragments m t21.etfrag k t21.Nifrag end symmetry NOSYM xc gga scf becke perdew end ETSNOCV print etslowdin endinput eor
Next one could do densf calculations, to view the natural orbitals in this method, see also the the documentation for the densf analysis program and the ADF-GUI. Input is the TAPE21 of the molecular calculation.
$ADFBIN/densf << eor GRID MEDIUM NOCV THRESH 0.01 END END INPUT eor mv TAPE41 nocv2.t41 $ADFBIN/densf << eor GRID MEDIUM NOCV ALL END END INPUT eor mv TAPE41 nocv3.t41
Sample directory: adf/NOCV_CrCO5-CH2/
Example for calculation of ETS-NOCV for spin-restricted fragments. ETS-NOCV: energy analysis using the Natural Orbitals for Chemical Valence. The CH2 molecule and Cr(CO)5 are the fragments, which form Cr(CO)5CH2 molecule.
First the two fragments are calculated.
$ADFBIN/adf -n1 << eor
Title CrCO5--[CH2] run from CrCO5 and CH2 closed shell fragments,FULL electron calc.!
atoms cartesian
C -0.429104 1.732058 -0.225052
H 0.407023 2.440417 -0.352323
H -1.385325 2.281354 -0.254124
end
basis
Type DZP
Core None
end
symmetry NOSYM
xc
gga becke perdew
end
endinput
eor
mv TAPE21 t21.CH2
$ADFBIN/adf -n1 << eor
Title [CrCO5] run
atoms cartesian
Cr -0.248053 -0.169062 0.005810
C -0.072963 -2.080685 0.229583
O 0.030811 -3.223220 0.361925
C -0.182894 0.049840 1.909128
O -0.142780 0.212309 3.050403
C -0.299940 -0.409118 -1.894730
O -0.331795 -0.521589 -3.042336
C -2.138631 -0.242152 0.075713
O -3.295036 -0.249916 0.115045
C 1.624487 0.092244 -0.083118
O 2.763411 0.288575 -0.140976
end
basis
Type DZP
Core None
Cr $ADFRESOURCES/TZP/Cr
end
symmetry NOSYM
xc
gga becke perdew
end
endinput
eor
mv TAPE21 t21.Crfragment
Next these fragments are used in the calculation of the full complex. The keys ETSNOCV and 'PRINT etslowdin' are needed in this case to to analyze the bonding in the molecule with respect to the fragments. The symmetry must be NOSYM. Note the '-n1' flag for adf: this enforces that a single-node (non-parallel) run is performed.
$ADFBIN/adf -n1 << eor
Title:CrCO5--[CH2], etsnocv activated by etsnocv and print etslowdin
atoms
C -0.429104 1.732058 -0.225052 f=f1
Cr -0.248053 -0.169062 0.005810 f=f2
C -0.072963 -2.080685 0.229583 f=f2
O 0.030811 -3.223220 0.361925 f=f2
C -0.182894 0.049840 1.909128 f=f2
O -0.142780 0.212309 3.050403 f=f2
C -0.299940 -0.409118 -1.894730 f=f2
O -0.331795 -0.521589 -3.042336 f=f2
C -2.138631 -0.242152 0.075713 f=f2
O -3.295036 -0.249916 0.115045 f=f2
C 1.624487 0.092244 -0.083118 f=f2
O 2.763411 0.288575 -0.140976 f=f2
H 0.407023 2.440417 -0.352323 f=f1
H -1.385325 2.281354 -0.254124 f=f1
end
fragments
f1 t21.CH2
f2 t21.Crfragment
end
symmetry NOSYM
xc
gga becke perdew
end
ETSNOCV RHOKMIN=1.e-3 EKMIN=1.5 ENOCV=0.05
print etslowdin
endinput
eor
Sample directory: adf/CH3_CH3_etsnocv/
Example for calculation of ETS-NOCV for unrestricted fragments. ETS-NOCV: energy analysis using the Natural Orbitals for Chemical Valence. The ethane molecule is built from two methyl radicals
First the two methyl fragments are calculated. The fragments should be spin-restricted.
$ADFBIN/adf << eor
Title CH3-CH3 built from CH3 radicals, FULL electron calc.!
atoms cartesian
C 0.019664 -0.034069 0.009101
H 0.039672 -0.069395 1.109620
H 1.063205 -0.065727 -0.341092
H -0.474230 -0.953693 -0.341621
end
basis
H $ADFRESOURCES/DZP/H
C $ADFRESOURCES/DZP/C
end
symmetry NOSYM
SCF
Iterations 2500
Converge 1E-6
end
xc
gga scf becke perdew
end
endinput
mv TAPE21 t21.frag1
$ADFBIN/adf << eor
Title CH3 radical
atoms cartesian
C -0.703210 1.217999 -0.497874
H -0.723753 1.252869 -1.598316
H -1.746567 1.250049 -0.147169
H -0.208833 2.137544 -0.147653
end
basis
H $ADFRESOURCES/DZP/H
C $ADFRESOURCES/DZP/C
end
symmetry NOSYM
SCF
Iterations 2500
Converge 1E-6
end
xc
gga scf becke perdew
end
endinput
eor
mv TAPE21 t21.frag2
Next these fragments are used in the calculation of the molecule ethane, using the key FRAGOCCUPATIONS to use an unrestricted fragment occupation for the methyl radicals, such that they are prepared for bonding. In the one fragment the singly occupied orbital will be an α-orbital, and in the other fragment the singly occupied orbital will be a β-orbital, such that the calculated Pauli repulsion between the fragments will be small.
The keys ETSNOCV and 'PRINT etslowdin-unrestricted' are needed in this case to to analyze the bonding in a molecule with unpaired electrons in the fragments. The symmetry must be NOSYM.
$ADFBIN/adf << eor
Title: final [CH3]-[CH3], etsnocv activated by etsnocv and etslowdin-unrestricted
atoms
C 0.019664 -0.034069 0.009101 f=f1
H 0.039672 -0.069395 1.109620 f=f1
H 1.063205 -0.065727 -0.341092 f=f1
H -0.474230 -0.953693 -0.341621 f=f1
C -0.703210 1.217999 -0.497874 f=f2
H -0.723753 1.252869 -1.598316 f=f2
H -1.746567 1.250049 -0.147169 f=f2
H -0.208833 2.137544 -0.147653 f=f2
end
fragments
f1 t21.frag1
f2 t21.frag2
end
fragoccupations
f1
A 5 // 4
subend
f2
A 4 // 5
subend
end
symmetry NOSYM
SCF
Iterations 800
Converge 1E-6
end
xc
gga scf becke perdew
end
ETSNOCV RHOKMIN=1.e-3 EKMIN=1.5 ENOCV=0.05
PRINT etslowdin-unrestricted
end input
eor
Next 2 densf calculations, to view the natural orbitals in this method, see also the the documentation for the densf analysis program and the ADF-GUI. Input is the TAPE21 of the molecular calculation.
$ADFBIN/densf << eor GRID MEDIUM NOCV Alpha 1 2 59 60 Beta 1 2 59 60 END END INPUT eor mv TAPE41 nocv1.t41 $ADFBIN/densf << eor GRID MEDIUM NOCV THRESH 0.01 END END INPUT eor mv TAPE41 nocv2.t41 eor
Sample directory: adf/Cntrs_NO2/
This example illustrates using the utility programs cntrs and densf. See the ADF manual for details.
$ADFBIN/adf << eor title NO2 atoms N 0 0 0 O 1.009356 0 0.464189 O -1.009356 0 0.464189 end Basis Type DZ Core Small End unrestricted charge 0 1 endinput eor
After the normal ADF calculation on NO2 has been completed, the utility program densf is executed to generate a TAPE41 file with user-specified items evaluated in a regular, user-specified grid.
The TAPE21 on which densf operates must be present as a local file with name TAPE21.
$ADFBIN/densf << eor density scf ortho frag fitdensity scf ortho frag orbitals scf alpha a1 1 2 a2 1 b1 2 beta a2 1 b1 1 2 b2 1 end coulpot frag ortho scf grid -7.5 -7.5 0.0 51 51 1.0 0.0 0.0 15.00 0.0 1.0 0.0 15.00 end end input eor
The charge density values in the grid are requested for all available types of density: exact and fitted, for the initial (sum-of-fragments), intermediate (orthogonalized fragments, see the ADF User's Guide) and final (SCF) situation.
Several SCF molecular orbitals are computed by specifying their indices in the energy-ordered list (a separate list for each symmetry subspecies).
The coulomb potentials (again: for sum-of-fragments, orthogonalized fragments, and SCF) are generated.
The grid is defined by an 'origin', the numbers of points in all independent grid directions and the direction vectors with the total grid size in each direction separately.
Since there are only two 'numbers-of-points' (51 each) a 2-dimensional grid is generated. 1D and 3D grids are also possible. See the ADF manual for a more detailed survey of the available options.
The result of the densf run is a file tape41 (binary, KF). This contains all computed data. tape41 can be used by cntrs to generate plot data.
$ADFBIN/cntrs << eor scan 0.02 0.05 0.10 0.2 0.5 0.0 -.02 -.05 -.10 -.2 -.5 end dash 0.2 file cont.d SCF%Density_A SCF%Density_B endinput eor
In this example eleven (11) scan values are defined to draw contours for, with a dash length of 0.2 bohr.
An ascii file cont.d will be opened by cntrs on which the specified items (SCF-densities for spin-up and spin-down) will be combined (by default: simply added) into one quantity.
For this quantity the contour lines that correspond to the specified scan values are stored. See the ADF manual for precise specifications and options.
gnuplot << eor set term dumb 100 80 set output "outplot" plot "cont.d" using 1:2 with lines eor cat outplot
The public domain software gnuplot (not included in the adf package) is applied here to display the result from cntrs. The resulting picture on your screen (if you have gnuplot available) looks like
10 *+
*
* ***********************
* ***** *****
* ** ************ **
* ** ********* ******** ********** **
* ** ** *********** **** ********** ** *
* * ***** *** **** **** *
* * * * ********** ********** **** *
8 *+ * * * ** ******* * * **
* * * * ** ** * * *
* * * * * * **** * * * *
* * * * * ** * * * * *** ** *
* * * * * * **** * * * *
* ** ** ** * ** * * *
* * * * ** ******* * * *
* ** * * ********** ********** **** *
* * ** *** ********* ****** **
* ** ** ********* ********* ** *
6 *+ ** *********** ****** ********** **
* *** ********** ***
* ********** **********
* ***********
Sample directory: adf/Cntrs.LocOrb_C2H2/
An illustration of the computation of localized molecular orbitals in C2H2.
The delocalized molecular orbitals as they result from the scf are localized in two different ways. In the first the three σ bonds are recombined only among themselves (no π bonds are mixed in), yielding two equivalent localized CH σ bonds and one localized σ bond. In a second step the localization of the remaining bond (the two π's) is performed, but this produces nothing new since no combination of the two π's is more localized then they are already by themselves.
$ADFBIN/adf << eor title C2H2, localization Sigma and Pi separately Atoms C 0 0 .63 C 0 0 -.63 H 0 0 1.63 H 0 0 -1.63 END Basis Type TZP Core Small End LocOrb alfa 4 5 alfa 1 2 3 END integration 4.0 end input eor
In the first localization cycle the π-orbitals are left out: #4 and #5 in the list of all occupied valence MOs: first 3 MOs of the first irreducible representation (s), then the 2 from the second irrep (π). In the second localization step the first three (meanwhile localized) orbitals are kept aside.
With densf the local orbitals can be computed in a user-defined grid (for plotting purposes). densf requires a file with name TAPE21.
$ADFBIN/densf << eor Grid 0. -5. -5. 100 100 0. 0. 1. 10. 0. 1. 0. 10. End Orbitals Local 1 2 3 4 5 End END INPUT eor
The program cntrs is applied to process the densf result file TAPE41.
$ADFBIN/cntrs << eor SCAN 0.01 0.02 0.04 0.08 0.16 0.32 0.64 1.28 -0.01 -0.02 -0.04 -0.08 -0.16 -0.32 -0.64 -1.28 END file ctr.a1 LocOrb%1 1.00 file ctr.a2 LocOrb%2 1.00 file ctr.a3 LocOrb%3 1.00 file ctr.a4 LocOrb%4 1.00 file ctr.a5 LocOrb%5 1.00 END INPUT eor
Again, gnuplot may be used to display the result on your screen.
$gnuplot << eor set term dumb 100 80 set output "outplot" plot "ctr.a1" using 1:2 with lines plot "ctr.a2" using 1:2 with lines plot "ctr.a3" using 1:2 with lines plot "ctr.a4" using 1:2 with lines plot "ctr.a5" using 1:2 with lines eor cat outplot
This results in 5 pictures, the first one looking like:
****************
**** ****
*** ***
** ***** ***** **
** *** *** **
** *** *** **
** ** ********** **
* ** **** **** **
** * *** *** * **
* ** ** ** ** *
* * *** ********** *** * *
* * * *** *** * * *
* * ** ** ** ** * *
***** * * * * * * * *
** ** * * * * ******** * * * *
* ** * * * ** **** **** ** * * *
* * *** * * *** *** * * ***
** **** * ***** * ** ** * *****
* *** **** ** **** * * * * **** *
** ** ****** ****** * * ****** *****
* * * *********** ************ *************
** * ** ********* ** ** **********
* * ** *** ** ******** * * ******** **
* ** * ** ** * ********* ********* * ** *
* * * *** ** *** ***** ****** ** * **
* * * *** * * ******* ****** ** ** **
* * * * ** * ********* ********* * ** *
* * ** *** ** ******** * * ******** **
** * ** ********** ** ** **********
* * ************ ************ *************
** ** ***** ****** * * ****** ******
* ***** * ***** * ** ** * ***** *
** ** *** * * *** *** * * ***
* ** * * * ** **** **** ** * * *
*** ** * * * * ******** * * * *
***** * * * * * * * *
* * ** ** ** ** * *
* * * *** *** * * *
* * *** ********** *** * *
* ** ** ** ** *
** * *** *** * **
* ** **** **** ** *
** ** ********** ** **
** *** *** **
** *** *** **
** ***** ***** **
** ********** **
*** ***
**** ****
****************
The second illustration of the computation of localized molecular orbitals in C2H2 combines directly all MOs (σ and π). This yields 3 equivalent 'banana' bonds, mixtures of σ and π, and two equivalent pure σ bonds.
$ADFBIN/adf << eor title C2H2 localization without frozen orbitals Atoms C 0 0 .63 C 0 0 -.63 H 0 0 1.63 H 0 0 -1.63 end fragments C t21.C H t21.H end integration 4.0 locorb end end input eor
Sample directory: adf/DOS_Cu4CO/
This sample illustrates the DOS property program to compute density-of-states data, for energy-dependent analysis.
First, the Cu4CO molecule is calculated (ADF), using single-atom fragments.
$ADFBIN/adf <<eor title Cu4CO (3,1) from atoms units length bohr end define rCu=2.784 end atoms 1. Cu rCu 0.0 0.0 2. Cu -rCu/2 rCu*sqrt(3)/2 0.0 3. Cu -rCu/2 -rCu*sqrt(3)/2 0.0 4. Cu 0.0 0.0 -rCu*sqrt(2) 5. C 0.0 0.0 2.65 6. O 0.0 0.0 4.91 end Basis Type TZP Core small end XC GGA PostSCF Becke Perdew END endinput eor
The PostSCF feature in the specification of the XC functional is used: the 'Becke-Perdew' GGA corrections are not included self-consistently but applied to the energy evaluation after the self-consistent LDA solution has been obtained.
The utility program dos requires a file named TAPE21 in the current directory, unless overridden using a TAPE21 keyword (not used in this example).
$ADFBIN/dos << eor file dostxt energyrange npoint=36 e-start=-25 e-end=10 tdos ! Cu 3d partial DOS gpdos a1 14:22 32:34 a2 5:10 e1:1 18:32 37:42 e1:2 18:32 37:42 end ! The same but using BAS gpdos BAS 17:34 57:74 97:114 137:154 end ! The same as above, but using much less complicated input gpdos ATYPE Cu d end ! Overlap PDOS between Cu 3d and CO 2p opdos ATYPE Cu 3d SUBEND ATOM 5 2p ATOM 6 2p end end input eor
Here, the total density of states, as well as various partial densities of states, are computed. You may feed the results found in the dostxt file into a plotting program such as gnuplot. The result is not displayed here. See the ADF manual for more detailed info about the dos program.
One can calculate Bade atomic charges and other Atoms in Molecule properties directly in ADF using a grid based method (from ADF2008.01 onwards). Another possibility for Bader's analysis, an example is described here, is to use the adf2aim utility such that a third party program Xaim can be used.
The ADF utility adf2aim (original name rdt21, now part of the ADF package) developed by Xavi López, Engelber Sans and Carles Bo converts an ADF TAPE21 to WFN format (for Bader analysis).
The WFN file is an input file for the third party program Xaim (see http://www.quimica.urv.es/XAIM for details), which is a graphical user interface to programs that can perform the Bader analysis.
Usage of adf2aim:
$ADFBIN/adf <<eor TITLE HF ATOMS 1. H .0000 .0000 .0000 2. F .0000 .0000 0.917 End Basis End End input eor $ADFBIN/adf2aim TAPE21 echo 'Contents of rdt21.res:' cat rdt21.res echo 'Contents of WFN:' cat WFN
Sample directory: adf/H2O_ADFNBO/
Dr. Autschbach, SCM, and Prof. Weinhold have collaborated to prepare a simple in
put file generator, called adfnbo, for
the GENNBO program of Prof. Weinholds Natural Bond Orbital (NBO) 5.0 package:
http://www.chem.wisc.edu/~nbo5
The GENNBO executable is included in the
ADF distribution and can be enabled via the license file
for all those who buy an NBO manual from either the NBO authors
or from SCM (info@scm.com).
Usage:
$ADFBIN/adf <<eor
Title simple NBO example for water
Atoms Z-Matrix
O 0 0 0
H 1 0 0 0.9
H 1 2 0 0.9 100
End
Basis
CORE NONE
TYPE DZ
End
FULLFOCK
AOMAT2FILE
SAVE TAPE15
SYMMETRY NOSYM
End Input
eor
$ADFBIN/adfnbo <<eor
write
fock
spherical
end input
eor
$ADFBIN/gennbo < FILE47
Note added: recommended is to use the key 'spherical' in the adfnbo input.
A File named FILE47 is generated by adfnbo which is an input file for the general NBO program gennbo. ADF needs to write some data to file, which is done by including these keywords in the adf input file:
FULLFOCK AOMAT2FILE SAVE TAPE15 SYMMETRY NOSYM
This section contains a brief summary of the capabilities of GENNBO, made available by Prof. Weinhold.
GENNBO implements most capabilities of the full NBO 5.0 program suite
as described on the NBO website:
http://www.chem.wisc.edu/~nbo5
These include determination of natural atomic orbitals (NAOs), bond orbitals
(NBOs), and localized MOs (NLMOs), as well as the associated NPA
(atomic charges and orbital populations) and NRT (resonance structures,
weightings, bond orders) valence descriptors, for a wide variety of
uncorrelated and correlated (variational, perturbative, or density
functional) theoretical levels. GENNBO-supported options include
all keywords except those explicitly requiring interactive communication
with the host electronic structure system (viz., $DEL deletions, NEDA,
NCS, NJC). The GENNBO program typically sits conveniently on the
PC desktop, ready to analyze (or re-analyze at will, with altered
options) the final results of a complex ADF calculation performed
on a remote cluster.
GENNBO "communicates" with the original ADF calculation through an archive file (JOB.47 file, preserving all necessary details of the final density) that is initially generated by ADF and subsequently becomes the input file for GENNBO. The .47 file contains a standard $NBO ... $END keylist that can be edited with a standard word processor or text editor to include chosen NBO keyword options, just as though they might have appeared in the original input stream of an interactive ADFNBO run. The stand-alone GENNBO program therefore allows many alternative NBO analysis options to be explored at leisure, without costly re-calculation of the wavefunction.
Sample directory: adf/AlCl3_efgnbo/
Example shows an NBO analysis of an EFG calculation for AlCl3.
In the ADF input one then needs to include the QTENS (EFG calculation) and include the subkey EFG of the key AORESPONSE. A higher integration accuracy in the core region is used (subkey ACCSPH of the key INTEGRATION). Other keywords are necessary because of the NBO analysis afterwards. Note that ADF, ADFNBO, and GENNBO have to run several times.
$ADFBIN/adf << eor Title AlCl3 atoms Al 0.000000 0.000000 -0.237368 Cl 1.808813 0.000000 0.807083 Cl 0.000000 0.000000 -2.326083 Cl -1.808813 0.000000 0.807083 End xc lda vwn gga revPBE end Symmetry NOSYM Integration accint 4.5 accsph 5.5 End Aoresponse efg 1 nbo end qtens save TAPE15 FULLFOCK AOMAT2FILE END INPUT eor $ADFBIN/adfnbo << eor write spherical fock TESTJOB end input eor rm -f adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 $ADFBIN/gennbo < FILE47 $ADFBIN/adfnbo << eor copy spherical fock TESTJOB end input eor $ADFBIN/adfnbo << eor read spherical fock TESTJOB end input eor rm -f adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 rm -f TAPE21 TAPE13 TAPE15 $ADFBIN/adf << eor Title AlCl3 atoms Al 0.000000 0.000000 -0.237368 Cl 1.808813 0.000000 0.807083 Cl 0.000000 0.000000 -2.326083 Cl -1.808813 0.000000 0.807083 End xc lda vwn gga revPBE end Symmetry NOSYM Integration accint 4.5 accsph 5.5 End Basis Type TZP Core none End Aoresponse efg 1 nbo end qtens End Input eor
Sample directory: adf/CH4_nmrnbo/
Example shows an NBO analysis of an NMR shielding calculation for CH4.
Some keywords are necessary because of the NBO analysis afterwards. A higher integration accuracy in the core region is used (subkey ACCSPH of the key INTEGRATION). First the scalar relativistic calculation is performed and the scalar relativistic localized orbitals are made. Next the spin-orbit coupled ADF calculation is done, and a calculation of NMR shielding constants is performed with an analyisis in scalar relativistic localized orbitals.
$ADFBIN/adf << eor ATOMS C 0.0000 0.0000 0.0000 H 0.6316 0.6316 0.6316 H 0.6316 -0.6316 -0.6316 H -0.6316 0.6316 -0.6316 H -0.6316 -0.6316 0.6316 END save TAPE15 FULLFOCK AOMAT2FILE BASIS type DZP core None END INTEGRATION accint 4.5 accsph 5.5 end relativistic scalar zora eor $ADFBIN/adfnbo << eor write spherical fock TESTJOB end input eor rm adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 $ADFBIN/gennbo < FILE47 $ADFBIN/adfnbo << eor copy spherical fock end input eor rm adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 rm TAPE15 TAPE21 logfile $ADFBIN/adf << eor ATOMS C 0.0000 0.0000 0.0000 H 0.6316 0.6316 0.6316 H 0.6316 -0.6316 -0.6316 H -0.6316 0.6316 -0.6316 H -0.6316 -0.6316 0.6316 END BASIS type DZP core None END SYMMETRY nosym INTEGRATION accint 4.5 accsph 5.5 end relativistic spinorbit zora end input eor rm TAPE15 TAPE10 logfile SINFO IINFO $ADFBIN/nmr << eor noscale nmr atoms 2 1 u1k best calc all :: noscl out iso tens end analysis print 0.01 canonical nbo components end end input eor
Sample directory: adf/CPL_CH3OH_NBO/
Example shows an NBO analysis of an NMR spin-spin coupling constants calculation for CH3OH.
Some keywords are necessary because of the NBO analysis afterwards. A higher integration accuracy in the core region is used (subkey ACCSPH of the key INTEGRATION). First the scalar relativistic calculation is performed and the scalar relativistic localized orbitals are made, and a calculation of NMR spin-spin coupling constants is performed with an analyisis in scalar relativistic localized orbitals. Next the spin-orbit coupled ADF calculation is done, and a calculation of NMR spin-spin coupling constants is performed with an analyisis in scalar relativistic localized orbitals.
$ADFBIN/adf << eor ATOMS 1 O 0.151078120000 -0.158942890000 -0.184382010000 2 H 0.762854510000 0.480823600000 0.187867830000 3 C 0.654254930000 -1.481762230000 0.026343630000 4 H 1.616760580000 -1.581906770000 -0.455670800000 5 H -0.035909520000 -2.200223490000 -0.393433960000 6 H 0.761359880000 -1.661537720000 1.087000640000 END save TAPE15 FULLFOCK AOMAT2FILE BASIS type DZP core None END SCF converge 1.0e-8 END SYMMETRY nosym INTEGRATION accint 4.5 accsph 5.5 end relativistic scalar zora end input eor $ADFBIN/adfnbo << eor write spherical fock TESTJOB end input eor rm adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 $ADFBIN/gennbo < FILE47 $ADFBIN/adfnbo << eor copy spherical fock end input eor $ADFBIN/adfnbo << eor spherical fock read end input eor rm adfnbo.37 adfnbo.39 adfnbo.49 adfnbo.48 rm TAPE15 TAPE21 TAPE13 logfile $ADFBIN/cpl << eor nmrcoupling xalpha dso pso sd scf convergence 1e-5 iterations 10 contributions 1e19 nbo nuclei 3 5 6 end endinput eor rm TAPE15 TAPE21 TAPE13 logfile $ADFBIN/adf << eor ATOMS 1 O 0.151078120000 -0.158942890000 -0.184382010000 2 H 0.762854510000 0.480823600000 0.187867830000 3 C 0.654254930000 -1.481762230000 0.026343630000 4 H 1.616760580000 -1.581906770000 -0.455670800000 5 H -0.035909520000 -2.200223490000 -0.393433960000 6 H 0.761359880000 -1.661537720000 1.087000640000 END BASIS type DZP core None END SYMMETRY nosym SCF converge 1.0e-8 END INTEGRATION accint 4.5 accsph 5.5 end relativistic spinorbit zora end input eor rm TAPE15 $ADFBIN/cpl << eor nmrcoupling xalpha dso pso sd scf convergence 1e-5 iterations 10 contributions 1e19 nbo nuclei 3 5 6 end end input eor
Sample directory: adf/BSSE_CrCO6/
A study of the Basis Set Superposition Error (BSSE) in the formation of Cr(CO)6. from CO and Cr(CO)5.
For the BSSE calculation special chemical elements must be created to describe the 'ghost' atoms, which have zero nuclear charge and mass. They do have basis functions (and fit functions), however, and they are used to calculate the lowering of the energy of the system to which the ghost atoms are added, due to the enlargement of the basis by the ghost basis. The ghost atom has the same basis and fit set as the normal element but no nuclear charge and no frozen core. The BASIS key recognizes elements denoted with Gh.atom in the ATOMS key as being ghost atoms. If the basis file specifies a frozen core ADF will treat it as if no frozen core is present.
The following calculations are carried out:
This series of calculations is carried out with basis set DZ.
Finally, the whole thing might be redone with basis set TZP, to show that the BSSE decreases with larger basis.
The calculations for the type DZ basis are contained in the sample script (with input- and output files). Those for type TZP bases can be set up easily and may be done as an exercise.
For the first series of calculations, with basis type DZ, the input files are discussed below and the relevant parts are echoed from the output files where the energy decomposition and the total bond energy are printed.
For the other series, using type TZP basis sets, only a summary of the results is given.
The calculations in this example all use:
For the BSSE calculations we first do the 'normal' calculations of CO and Cr(CO)5, yielding the fragment files t21.CO and t21.CrCO5. The input files for these calculations are not shown here.
For the CO BSSE calculation the standard CO must have been computed first. In the BSSE run a Cr(CO)5 ghost fragment basis set is then added to the 'normal' CO input. Important is the use of the BASIS key. In this case the BASIS key is used for the generation of the ghost atoms, it should have the same definition for the atoms as will be used later for the Cr(CO)5 fragment. The FRAGMENTS key is used for the fragment CO. The energy change (the printed 'bond energy' in the output) is the BSSE.
The input file for the CO-BSSE run is:
title BSSE for CO due to Cr(CO)5 ghost
noprint sfo,frag,functions
atoms
Gh.Cr 0 0 0
Gh.C -1.86 0 0
Gh.C 1.86 0 0
Gh.C 0 1.86 0
Gh.C 0 -1.86 0
Gh.C 0 0 -1.86
Gh.O 3.03 0 0
Gh.O -3.03 0 0
Gh.O 0 3.03 0
Gh.O 0 -3.03 0
Gh.O 0 0 -3.03
C 0 0 1.86 f=CO
O 0 0 3.03 f=CO
end
Basis
Type DZ
Core Small
end
fragments
CO CO.t21
end
symmetry C(4V)
integration 4
endinput
In the output we find in the Bond Energy section:
hartree eV kcal/mol kJ/mol
-------------------- ----------- ---------- -----------
Pauli Repulsion
Kinetic (Delta T^0): 0.000000000000009 0.0000 0.00 0.00
Delta V^Pauli Coulomb: -0.000000000000007 0.0000 0.00 0.00
Delta V^Pauli LDA-XC: -0.000000000000003 0.0000 0.00 0.00
-------------------- ----------- ---------- -----------
Total Pauli Repulsion: -0.000000000000001 0.0000 0.00 0.00
(Total Pauli Repulsion =
Delta E^Pauli in BB paper)
Steric Interaction
Pauli Repulsion (Delta E^Pauli): -0.000000000000001 0.0000 0.00 0.00
Electrostatic Interaction: -0.000000000000017 0.0000 0.00 0.00
(Electrostatic Interaction =
Delta V_elstat in the BB paper)
-------------------- ----------- ---------- -----------
Total Steric Interaction: -0.000000000000018 0.0000 0.00 0.00
(Total Steric Interaction =
Delta E^0 in the BB paper)
Orbital Interactions
A1: -0.001838638722848 -0.0500 -1.15 -4.83
A2: 0.000000000000000 0.0000 0.00 0.00
B1: 0.000000000000000 0.0000 0.00 0.00
B2: 0.000000000000000 0.0000 0.00 0.00
E1: -0.002025936656647 -0.0551 -1.27 -5.32
-------------------- ----------- ---------- -----------
Total Orbital Interactions: -0.003864575379498 -0.1052 -2.43 -10.15
Alternative Decomposition Orb.Int.
Kinetic: -0.056036605580477 -1.5248 -35.16 -147.12
Coulomb: 0.048666195764206 1.3243 30.54 127.77
XC: 0.003505834436773 0.0954 2.20 9.20
-------------------- ----------- ---------- -----------
Total Orbital Interactions: -0.003864575379498 -0.1052 -2.43 -10.15
Residu (E=Steric+OrbInt+Res): -0.000000000000003 0.0000 0.00 0.00
Total Bonding Energy: -0.003864575379519 -0.1052 -2.43 -10.15
Summary of Bonding Energy (energy terms are taken from the energy decomposition above)
======================================================================================
Electrostatic Energy: -0.000000000000017 0.0000 0.00 0.00
Kinetic Energy: -0.056036605580468 -1.5248 -35.16 -147.12
Coulomb (Steric+OrbInt) Energy: 0.048666195764196 1.3243 30.54 127.77
XC Energy: 0.003505834436770 0.0954 2.20 9.20
-------------------- ----------- ---------- -----------
Total Bonding Energy: -0.003864575379519 -0.1052 -2.43 -10.15
The BSSE for CO is computed as 2.43 kcal/mole
In similar fashion the BSSE is computed for Cr(CO)5. In the BSSE run a ghost atoms C and O at the positions they will have in the Cr(CO)6 molecule are added to the normal Cr(CO)5 input:
title BSSE for Cr(CO)5 due to CO ghost noprint sfo,frag,functions atoms Cr 0 0 0 f=CrCO5 C 1.86 0 0 f=CrCO5 C -1.86 0 0 f=CrCO5 C 0 1.86 0 f=CrCO5 C 0 -1.86 0 f=CrCO5 C 0 0 -1.86 f=CrCO5 O 3.03 0 0 f=CrCO5 O -3.03 0 0 f=CrCO5 O 0 3.03 0 f=CrCO5 O 0 -3.03 0 f=CrCO5 O 0 0 -3.03 f=CrCO5 Gh.C 0 0 1.86 Gh.O 0 0 3.03 end Basis Type DZ Core Small end fragments CrCO5 CrCO5.t21 end symmetry C(4v) integration 4 endinput
The Bond Energy result yields 1.93 kcal/mole for the BSSE.
The bonding of CO to Cr(CO)5 is computed in the normal way: from fragments CO and Cr(CO)5. The obtained value for the bond energy can then simply corrected for the two BSSE terms, 4.36 kcal/mole together.
The two BSSE runs (#2 and #4 in the list above) can also be repeated, but now with the core orthogonalization functions omitted from the ghost bases. To to this one can not use the BASIS key, but one needs to explicitely 'create' the ghost atoms. This will not be done here, but only the results will be discussed. One may argue about whether these functions should be included in the ghost basis sets, but since they are very contracted around the ghost nuclei they are not expected to contribute significantly anyway and may then just as well be omitted. This has been explicitly verified by test examples. /p>
The Core Functions (the functions in the valence basis set that serve only for core-orthogonalization, for instance the 1S 5.40 in the Carbon basis set (see the $ADFHOME/atomicdata/DZ/C.1s database file) are removed from the Create data files used for the creation of the ghost atoms.
This yields as BSSE values for CO and Cr(CO)5 respectively 2.33 and 1.87 kcal/mole (compare 2.43 and 1.93 kcal/mole for the case with Core Functions included). The net total effect of including/removing the Core Functions is therefore (2.43-2.33)+(1.93-1.87)=0.16 kcal/mole. This is an order of magnitude smaller than the BSSE effect itself.
BSSE effects should diminish with larger bases and disappear in the limit of a perfect basis. This can be studied by comparing the BSSE for basis DZ, see above, with the BSSE for basis TZP. The procedure is completely similar to the one above and yields:
For the BSSE terms: 0.7 kcal/mole for CO (compare: 2.4 kcal/mole for basis DZ), and 0.6 kcal/mole for Cr(CO)5 (1.9 for basis DZ)
The total BSSE drops from 4.4 kcal/mole in basis DZ to 1.3 in basis TZP.
A systematic study with adf
of the BSSE in metal-carbonyl complexes can be found in
Rosa, A., et al., Basis Set Effects in Density Functional Calculations on the Metal-Ligand and Metal-Metal Bonds of Cr(CO)5-CO and (CO)5.
Journal of Physical Chemistry, 1996, 100: p. 5690-5696.
Sample directory: adf/SCF_Ti2O4/
One can run into SCF convergence problems when calculating certain types of systems. Some of the notorious examples are transition metal oxides and lantanide compounds. Below, several approaches to solving the SCF convergence problem are demonstrated.
NewDIIS keyword
The first approach is to try a new DIIS algorithm, which will probably become default in a future version. The new algorithm is switched on by using the keyword NewDIIS anywhere in the input file:
$ADFBIN/adf << eor Title Ti2O4 SCF aid test (NewDIIS) Atoms Ti 1.730 0.000 0.000 Ti -1.730 0.000 0.000 O 0.000 1.224 0.000 O 0.000 -1.224 0.000 O 3.850 0.000 0.000 O -3.850 0.000 0.000 End XC GGA Becke Perdew End Basis Type DZ Core Small End SCF Iterations 300 End NewDIIS End input eor
Multi-step smearing
Second approach is an extension to the so-called "electron smearing" method. In this method, the electrons are distributed among orbitals around Fermi-level using a pseudo-thermal distribution function. Although the result with fractional occupation number has no physical sense, the method can be used to achieve integer occupation numbers by reducing the smearing parameter step-wise. In the example above, replace the NewDIIS keyword with the following line of text:
Occupations Smear=0.2,0.1,0.07,0.05,0.03,0.02,0.01,0.007,0.005,0.001
A few notes:
Steepest descent method
The third example demonstrates the use of the Occupations Steep= option (see the User's Guide for details). There are two differences from the previous example shown below:
SCF Iterations 300 Mixing 0.05 DIIS N=0 End Occupations Steep=0.5,0.3
One difference is, obviously, in the Occupations keyword. The other difference is more subtle. For stable convergence, it is often essential to switch off DIIS and set the mixing parameter to a low value. Of course, it will make convergence quite (sometimes very) slow. Ultimately you should get either an aufbau configuration or a configuration with exactly degenerate HOMO. In this example, the result is an aufbau solution.
Both methods should, in principle, give the same result, which is the case in this example.
A-DIIS
The fourth example uses the so called A-DIIS method. The A-DIIS method combines the strength of the ARH and DIIS methods. It does not require energy evaluation so it is much cheaper than the ARH and Energy-DIIS methods.
SCF Iterations 300 Mixing 0.05 ADIIS End
Energy-DIIS
The fifth example uses the so called Energy-DIIS method. Please note that similar to ARH and unlike the standard SCF procedure in ADF this method requires energy evaluation at each SCF cycle, which makes it significantly slower compared to energy-free SCF.
SCF Iterations 300 Mixing 0.05 EDIIS End
Augmented Roothaan-Hall
The sixth example uses the Augmented Roothaan-Hall (ARH) method. The basic idea of this method is that the density matrix is optimized directly to minimize the total energy. Important: the ARH method can be used with SYMMETRY NOSYM only.
Symmetry NOSYM SCF Iterations 300 Mixing 0.05 ARH End
LISTi
The seventh example uses the LISTi method. LISTi is very similar to the usual DIIS but typically it performs much better. It is also computationally less expensive and scales better in parallel even though DIIS is rarely a scaling bottleneck.
SCF Iterations 300 LISTi End
Sample directory: adf/Freq_NH3_Scan/
Sometimes spurious imaginary frequencies are calculated where one would expect a very low (nearly zero) frequency. Most frequently this happens when there is a barrier-free rotation of, for example, methyl groups. The SCANFREQ keyword allows one to rescan calculated frequencies in order to find out if they wre calculated accurately.
In this example analytical frequencies are calculated. Next recalculation of certain NH3 frequencies are performed by scanning along normal modes from a restart file. In this calculation the frequencies are calculated numerically with finite displacements using symmetry.
$ADFBIN/adf <<eor title NH3 analytic frequencies atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Basis Type TZP Core Small End AnalyticalFreq End integration 5.0 end input eor mv TAPE21 NH3_anl.t21 $ADFBIN/adf <<eor title Re-calculate NH3 frequencies by scanning along normal modes from a restart file atoms N 0.0000 0.0000 0.0000 H 0.4729 0.8190 0.3821 H -0.9457 0.0000 0.3821 H 0.4729 -0.8190 0.3821 end Fragments N t21.N H t21.H End ScanFreq 0 4000 Restart NH3_anl.t21 integration 5.0 end input eor
Sample directory: adf/BakersetSP/
In this example you will find how to use adfprep to run a particular job (a single point calculation in this case) for all molecules in the Baker set. The molecules are simply xyz files and contain no ADF specific information. adfreport is used to collect the resulting bonding energies.
rm -f runset
for f in $ADFHOME/examples/adf/BakersetSP/Bakerset/*.xyz
do
"$ADFBIN/adfprep" -t SP -i 2.5 -b DZ -c Large -m "$f" -j `basename $f .xyz`>> runset
done
chmod +x runset
./runset
echo Results
ls -t -1 *.t21 | while read f
do
"$ADFBIN/adfreport" "$f" BondingEnergy
done
echo Ready
Sample directory: adf/ConvergenceTestCH4/
In this example you will find how to use adfprep to test convergence of the bonding energy with respect to basis set and integration accuracy. adfreport is used to collect the resulting bonding energies.
rm -f runset
for b in SZ DZ DZP TZP TZ2P QZ4P
do
"$ADFBIN/adfprep" -t "$ADFHOME/examples/adf/ConvergenceTestCH4/methane.adf" \
-b $b -j methane.$b >> runset
done
chmod +x runset
./runset
echo Results
echo Basis set convergence of Bonding Energy, SZ DZ DZP TZP TZ2P QZ4P
for b in SZ DZ DZP TZP TZ2P QZ4P
do
"$ADFBIN/adfreport" "methane.$b.t21" BondingEnergy
done
rm -f runset
for i in 2 3 4 5
do
"$ADFBIN/adfprep" -t "$ADFHOME/examples/adf/ConvergenceTestCH4/methane.adf" \
-b DZP -i $i -j methane.$i >> runset
done
chmod +x runset
./runset
echo Integration convergence of Bonding Energy, 2 3 4 5
for i in 2 3 4 5
do
"$ADFBIN/adfreport" "methane.$i.t21" BondingEnergy
done
echo Ready
| 3DRISM-Glycine [1] | ESR_HfV [1] | Hplus_CO_etsnocv [1] |
| AgI_asoexcit [1] | ESR_HgF_2der [1] | Hyperpol [1] |
| AgI_SO_Pol [1] | ESR_TiF3 [1] | LT_constraint [1] |
| AIM_HF [1] | FDE_Energy_H2O-Ne_unrestricted [1] | MBH_CH4 [1] |
| AlCl3_efgnbo [1] | FDE_Energy_NH3-H2O [1] | MBH_Ethanol [1] |
| AT_transferintegrals [1] | FDE_H2O_128 [1] | methane_dimer_dispersion [1] |
| Au2_Resp [1] | FDE_HeCO2_freezeandthaw [1] | MM_Dispersion [1] |
| Au2_ZORA [1] | FDE_NMR_relax [1] | ModStPot_N2+ [1] |
| AuH_analyse_exciso [1] | Fe4S4_BrokenSymm [1] | Mossbauer [1] |
| Bader [1] | Field.PtCO [1] | N2_TDHF [1] |
| BakersetSP [1] | Frags_NiCO4 [1] | Ne_CoreExci [1] |
| BondOrder [1] | Frags_PtCl4H2 [1] | Ne_exciso [1] |
| BSSE_CrCO6 [1] | FranckCondon_NO2 [1] | NH_ZFS [1] |
| C2H4_TDCDFT [1] | Freq_NH3 [1] | NMR_B3LYP [1] |
| CEBE_NNO [1] | Freq_NH3_RAMAN [1] | NMR_NICS [1] |
| CH3_CH3_etsnocv [1] | Freq_NH3_Scan [1] | NOCV_CrCO5-CH2 [1] |
| CH4_nmrnbo [1] | Freq_UF6 [1] | OH_MetaGGA [1] |
| CH4_SAOP [1] | GO_constraints [1] | PbH4_finitenuc [1] |
| CH4_SecDeriv [1] | GO_FDE_H2O-Li [1] | PCCP_Unr_BondEnergy [1] |
| CN_SecDeriv [1] | GO_FDE_NH3-H2O [1] | pdb2adf [1] |
| CN_unr_exci [1] | GO_Formaldehyde [1] | QMMM_Butane [1] |
| Cntrs.LocOrb_C2H2 [1] | GO_H2O [1] | QMMM_CYT [1] |
| Cntrs_NO2 [1] | GO_LiF_Efield [1] | QMMM_Surface [1] |
| CO_model [1] | GO_restraint [1] | SCF_Ti2O4 [1] |
| ConvergenceTestCH4 [1] | green_Al [1] | SD_CrNH3_6 [1] |
| CPL_C2H2 [1] | green_Au [1] | SiH2_spinflip [1] |
| CPL_CH3OH_NBO [1] | green_BDT [1] | SO_ZORA_Bi2 [1] |
| CPL_HF_hybrid [1] | H2O_ADFNBO [1] | Solv_HCl [1] |
| CuH+_S-squared [1] | H2O_HF_freq [1] | SUBEXCI_dimer [1] |
| DampedVerdet [1] | H2O_magnet [1] | TiCl4_CoreExci [1] |
| DelocalGO_aspirin [1] | H2O_MCD [1] | Tl_noncollinear [1] |
| Diimina_NOCV [1] | H2O_MCD_ZFS [1] | TlH_SO_analysis [1] |
| Disper_HF [1] | H2O_TD_magnet [1] | Transit_H2O [1] |
| DMO_CD [1] | H2O_Verdet [1] | TS_C2H6 [1] |
| DMO_ORD [1] | H2PO_B3LYP [1] | TS_CH4_HgCl2 [1] |
| DMO_ORD_aoresponse [1] | HBr [1] | TSRC_SN2 [1] |
| DOS_Cu4CO [1] | HBr_SO [1] | Twist_Ethene_TDDFT [1] |
| EDA_meta_gga_hybrid [1] | HCN [1] | UnrFrag_H2 [1] |
| Efield.PntQ_N2 [1] | HCN_CINEB [1] | VCD_COG_NHDT [1] |
| EGO_CH2_sf [1] | HF_ResonanceRaman [1] | Vibron_RR_uracil [1] |
| EGO_CH2O_trip_constr [1] | HgMeBr_pnr [1] | VO_collinear [1] |
| EGO_N2 [1] | HgMeBr_psc [1] | VROA [1] |
| EGO_N2_EIGENF [1] | HgMeBr_zso [1] | VROA_RESO [1] |
| EGO_PH2 [1] | HI_EFG [1] | ZORA_GO_AuH [1] |
| Energy_H2O [1] | HI_SecDer_ZORA [1] |
