As appearing at pubs.acs.org/toc/cmatex/32/5 – icing on the cake that was the accepted article (as described at this link and courtesy of VMD and GIMP).

## A Quick Guide To Running DFTB (With Available Parameters) In GAMESS-US

This post comes out of my great appreciation for just how well Yoshio Nishimoto's DFTB (density functional-based tight binding method) implementation in GAMESS-US runs, both as an additional functionality in an already considerable program and in comparison to a few other programs I've worked with to do the same.

Also, the use of unmodified Slater-Koster files in the GAMESS-US implementation is a nice touch.

This all begins at the dftb.org website with the downloading of available Slater-Koster parameter files for available sets of elements. Note that your favorite elements might not yet have parameters – or parameters within any given parameter set – for immediate use. Until others publish new parameter sets and post them somewhere – and if you're an academic – you might consider giving Stefan Grimme's XTB serious consideration.

Additionally! A recent find while looking for new parameter sets was the KIST Integrated Force Field Platform (kiff.vfab.org), providing the complete set of DFTB parameters and up-to-date ReaxFF parameters as well.

## Basic Input Format

A few key points (highlighted in the generic pointer input file) below:

! $scf soscf=.t. fdiff=.f. shift=.f. extrap=.f. damp=.t. diis=.f. $end

$system modio=31 $end

$basis gbasis=dftb $end

! $dftb ndftb=2 dampxh=.t. dampex=4.0 itypmx=0 etemp=300 $end

$dftb ndftb=3 dampxh=.t. dampex=4.0 disp=skhp etemp=300

disppr(1)=0.31,0.386,0.386,0.000,0.000,0.000,0.000,3.5,3.5,3.5,

3.5,3.5,3.5,0.80,0.69,1.382,1.382,1.382,1.064,1.064,1.064,3.8,3.8,

3.8,3.8,3.8,3.8,2.50 $end

$dftb hubder(1)=-0.1857,-0.1492 $end

$dftbsk

C C "/path/to/skfiles/3ob-3-1/C-C.skf"

C Si "/path/to/skfiles/3ob-3-1/C-Si.skf"

C H "/path/to/skfiles/3ob-3-1/C-H.skf"

Si C "/path/to/skfiles/3ob-3-1/Si-C.skf"

Si Si "/path/to/skfiles/3ob-3-1/Si-Si.skf"

Si H "/path/to/skfiles/3ob-3-1/Si-H.skf"

H C "/path/to/skfiles/3ob-3-1/H-C.skf"

H Si "/path/to/skfiles/3ob-3-1/H-Si.skf"

H H "/path/to/skfiles/3ob-3-1/H-H.skf"

$end $data

C1

C 6.0 -2.775607 0.000000 0.000000

...

H 1.0 -2.809590 0.000420 50.594100

$end

1. There have been many improvements to the DFTB code. The myweb.liu.edu/~nmatsuna/gamess manual (among others hosted *not* on the official GAMESS-US website) is among the first that *always* come up when searching out GAMESS-US keywords – and it's several years out of date (esp. for DFTB keywords). The 2019.1 manual has a significantly expanded **$DFTB** section, including calls for including dispersion corrections from the original **$DFT** block.

2. No **$SCF** – as the manual states, converger specifications are pre-defined in the **ITYPMX** keyword, with additional keywords in the **$DFTB** block (esp. **ETEMP**) available to aid in problematic convergence.

3. One big **$DFTBSK** – as addressed on Jan Jensen's initial post about running DFTB in GAMESS-US, GAMESS-US will only read atom pairs needed for the atoms listed in the **$DATA** section. You are fine to simply make an all-encompassing **$DFTBSK** block for all pairs in a parameter folder and write that to the input file.

4. **DISPPR** for **SKHP** – a mild manual issue. For the "Slater–Kirkwood + bond number polarizability dependence" dispersion correction, the current manual states that "For **DISP**=**SKHP**, a set for a species has 14 parameters. The first six are the polarizabilities depending on the number of bonds, and the next six are cutoff length, and the last is atomic charge." You'll note that 6 + 6 + 1 adds up to 13. The 14 numbers needed for each element in this keyword begins with the covalent atomic radius for the element. In the absence of a list for those values, Table S1.1 (freely available) from "DFTB Parameters for the Periodic Table: Part 1, Electronic Structure" is an excellent resource. I did not bother to address the "different covalent radii values for different hybridizations" issue and have seen little on the topic related to DFTB calculations. For the other 13 numbers, the DFTB+ manual contains a set in Appendix E and reproduced below, with notes to consider other publications for more and different versions of the same elements (and note the NOTES column in the Appendix E table for P and S). Local copy of the manual – 2019Aug24_DFTBplus_manual.pdf

Element Polarisability (6 #'s) ---- [Å3]Cutoff (6 #'s) -- [Å]Chrg O 0.560 0.560 0.000 0.000 0.000 0.0003.8 3.8 3.8 3.8 3.8 3.8 3.15 N 1.030 1.030 1.090 1.090 1.090 1.0903.8 3.8 3.8 3.8 3.8 3.8 2.82 C 1.382 1.382 1.382 1.064 1.064 1.0643.8 3.8 3.8 3.8 3.8 3.8 2.50 H 0.386 0.386 0.000 0.000 0.000 0.0003.5 3.5 3.5 3.5 3.5 3.5 0.80 P 1.600 1.600 1.600 1.600 1.600 1.6004.7 4.7 4.7 4.7 4.7 4.7 4.50 S 3.000 3.000 3.000 3.000 3.000 3.0004.7 4.7 4.7 4.7 4.7 4.7 4.80

5. **HUBDER** – If you don't specify Hubbard Derivatives in your input file, GAMESS-US will include them from internal values. Values appear to be the set from the 3ob-3-1 parameter set and are available in the set's README file (reproduced from that file below for available elements).

List of all atomic Hubbard derivatives (atomic units): Br = -0.0573 C = -0.1492 Ca = -0.0340 Cl = -0.0697 F = -0.1623 H = -0.1857 I = -0.0433 K = -0.0339 Mg = -0.02 N = -0.1535 Na = -0.0454 O = -0.1575 P = -0.14 S = -0.11 Zn = -0.03

6. **HUBDER** order – the order for these numbers is as per the FIRST appearance of the element in the **$DATA** block. This is ieasy to see if you don't specity the **HUBDER** values in the input file and simply let GAMESS-US write out the "HUBBARD DERIVATIVES" block in your DFTB3 run.

7. **$DFT** Dispersion – according to the manual, you can now use more the traditional dispersion corrections in **$DFTB** calculations. In the **$DFTB** block, use **DISP=DFT,** then modify your **$DFT** section as you otherwise would – and include **DC=.T. **in your **$DFT** block.

8. DFTB Block Summary in the output file – below for a DFTB3 run of a CH-containing molecule using the 3ob-3-1 SK parameters:

********************************************************** ** DENSITY-FUNCTIONAL TIGHT-BINDING (DFTB) CALCULATION ** ********************************************************** WRITTEN BY YOSHIO NISHIMOTO (NAGOYA UNIVERSITY) NUMBER OF SPECIES: 2 SPECIES 1 : H WITH S ORBITAL SPECIES 2 : C WITH S+P ORBITALS $DFTB OPTIONS ------------- NDFTB = 3 SCC = T DFTB3 = T SRSCC = F DAMPXH = T DAMPEX = 4.00 DISP =NONE ITYPMX = 0 ETEMP = 0.00 MODESD = 0 MODGAM = 8 PRTORB = F --- SCC CALCULATION --- INCLUDE 3RD ORDER CORRECTION SETTING HUBBARD DERIVATIVES FOR DFTB3 H : -0.14920 (USER DEFINED) C : -0.18570 (USER DEFINED) USE X-H DAMPING: 4.00000 BROYDEN'S CHARGE MIXING 4 SLATER-KOSTER FILES WILL BE READ START READING SLATER-KOSTER FILES 1 1 (H - H ) = /home/ec2-user/3ob-3-1/H-H.skf 1 2 (H - C ) = /home/ec2-user/3ob-3-1/H-C.skf 2 1 (C - H ) = /home/ec2-user/3ob-3-1/C-H.skf 2 2 (C - C ) = /home/ec2-user/3ob-3-1/C-C.skf

## DFTB2

Starting with the less keyword-involved run, below is a completely generic DFTB2 input file with no funny business in the **$DFTB** block to aid in convergence (my first test would be to change **ETEMP** to something higher, which has often done the trick when a run won't complete its first SCF cycle in 200 steps).

$contrl runtyp=optimize icharg=0 maxit=200 $end $system mwords=100 $end $system modio=31 $end $basis gbasis=dftb $end $dftb ndftb=2 dampxh=.t. dampex=4.0 itypmx=0 etemp=0 $end $data C1 C 6.0 -2.775607 0.000000 0.000000 … H 1.0 -2.809590 0.000420 50.594100 $end $dftbsk C C "/path/to/dftb/params/pbc-0-3/C-C.skf" C H "/path/to/dftb/params/pbc-0-3/C-H.skf" H C "/path/to/dftb/params/pbc-0-3/H-C.skf" H H "/path/to/dftb/params/pbc-0-3/H-H.skf" $end

## DFTB3

A more involved DFTB3 input file, including (1) the **SKHP** dispersion-correction and carriage returns to show that all 14 numbers are there – in order of appearance – for the atoms in **$DATA**, (2) the **$SCF** block as a comment that is ignored by DFTB calculations, (3) a reordering of all keywords above the **$DATA** block (because why not?), (4) the **HUBDER** values in order of appearance in the **$DATA** block, and (5) the Si-containing pairs in the **$DFTBSK** block just to show that their presence isn't an issue for the run.

! $scf soscf=.t. fdiff=.f. shift=.f. extrap=.f. damp=.t. diis=.f. $end $system modio=31 $end $basis gbasis=dftb $end $dftb ndftb=3 dampxh=.t. dampex=4.0 disp=skhp etemp=300 disppr(1)= 0.31,0.386,0.386,0.000,0.000,0.000,0.000,3.5,3.5,3.5,3.5,3.5,3.5,0.80, 0.69,1.382,1.382,1.382,1.064,1.064,1.064,3.8,3.8,3.8,3.8,3.8,3.8,2.50 $end $dftb hubder(1)=-0.1857,-0.1492 $end $dftbsk C C "/path/to/skfiles/3ob-3-1/C-C.skf" C Si "/path/to/skfiles/3ob-3-1/C-Si.skf" C H "/path/to/skfiles/3ob-3-1/C-H.skf" Si C "/path/to/skfiles/3ob-3-1/Si-C.skf" Si Si "/path/to/skfiles/3ob-3-1/Si-Si.skf" Si H "/path/to/skfiles/3ob-3-1/Si-H.skf" H C "/path/to/skfiles/3ob-3-1/H-C.skf" H Si "/path/to/skfiles/3ob-3-1/H-Si.skf" H H "/path/to/skfiles/3ob-3-1/H-H.skf" $end $data C1 H 1.0 -2.775607 0.000000 0.000000 C 6.0 -3.775607 0.000000 0.000000 .... H 1.0 -2.809590 0.000420 50.594100 $end

And, just to keep you from bouncing between tabs, the $DFTB section from the 2019.1 version of the manual (obviously check there at some point for changes and improvements).

$DFTB group (relevant for GBASIS=DFTB) Density-functional tight-binding (DFTB) is turned on by selecting GBASIS=DFTB in $BASIS. $DFTB controls optional parameters for a DFTB calculation. DFTB is formulated in a two-center approximation utilizing implicitly a minimal pseudoatomic orbital basis set with corresponding, pretabulated one- and two-center integrals. Because of this, many properties (for instances, multipoles higher than dipoles) and many options are ignored or not available in the current implementations of DFTB. DFTB also uses an independent SCF driver (SCF in DFTB is also called SCC, see below), so most SCF options are not available for DFTB. Only SCFTYP=RHF and UHF are implemented. SCFTYP=ROHF is available, only when all SPNCST values are zero. DFTB does not explicitly use symmetry (C1 throughout) since integrals are never computed during the calculations. Slater-Koster tables are only defined for spherical functions (5d) so DFTB sets ISPHER=1. Most $GUESS options do not work for DFTB (DFTB does not use initial orbitals in the usual sense). Other than the default (METHOD=HUCKEL, which is ignored), only METHOD=MOREAD works (note that SCC-DFTB can use initial charges on atoms, derived from the orbitals). RUNTYP=OPTIMIZE, HESSIAN and RAMAN are available for full (non-FMO) DFTB and FMO-DFTB. Excited state calculations for full DFTB may be performed through the standard (linear- response) time-dependent formalism (only closed shell). PCM can be used for both ground and excited state calculations, and energy and gradient can be evaluated. In DFTB calculation, the atom type is determined by its name, not its nuclear charge as elsewhere in GAMESS. The nuclear charge (the second column in $DATA) is used only in population analysis, but not in SCF. DFTB uses a notion of "species", which means an atomic type. The species are numbered according to the order in which atoms appear in $DATA. For instances, in water there are two species, O and H. An atomic type of each species needs MAXANG, which for most but not all atoms is set automatically. NDFTB order of the Taylor expansion of the total energy around a reference density in the DFTB model. = 1 NCC-DFTB, also called DFTB1. NCC stands for non-charge-consistent, i.e., no explicity charge-charge interaction term is included in the energy calculation. = 2 SCC-DFTB, also called DFTB2. SCC means a self-charge-consistent approach, and SCC implies that SCF iterations are carried out that converge monopolar charges towards self-consistency. = 3 DFTB3, including 3rd order correction using Hubbard derivatives (HUBDER). In order to reproduce the published DFTB3 approach, it is necessary to also specify DAMPXH=.TRUE. to add other terms. Gaus, M. et al. J. Chem. Theory Comput. 2011, 7, 931-948 is referred to as Gaus2011 below. Default: 2. DAMPXH = a flag to include the damping function for X-H atomic pair in DFTB3. See also DAMPEX, and eq 21 in Gaus2011. The damping function is used when at least one atom in a pair is "H". "HYDROGEN" and any other name will turn off the damping. Default: .FALSE. DAMPEX = an exponent used in the damping function for X-H atomic pairs. The default value is 4.0 (taken from the 3OB parameter set). SRSCC = a flag to perform shell-resolved SCC calculation. If set to .FALSE., the code uses the Hubbard value for an s orbital for p and d orbitals, ignoring their Hubbard values defined in Slater- Koster tables. Using .TRUE. enables the use of proper Hubbard values for p and d orbitals, implemented only for DFTB1 and DFTB2. Default: .FALSE. ITYPMX Convergence method of SCC calculations. = -1 Use standard GAMESS convergence methods. SOSCF and DIIS are supported, but DEM is not. = 0 Broyden's method. Interpolation is applied for atomic (or shell-resolved when SRSCC=.TRUE.) charges, but not Hamiltonian matrix. = 1 (reserved) = 2 DIIS for charges. Default: 0. ETEMP = electronic temperature in Kelvin. Non-zero values of ETEMP help SCF convergence of nearly-degenerate systems by smearing occupation numbers around the Fermi-level. Only the Fermi-Dirac distribution function is available as a smearing function. The default value is 0 Kelvin, meaning the smearing function is not used. ETEMP is implemented only for SCFTYP=RHF and when FMO is not used. DISP dispersion model for DFTB. = NONE no Dispersion correction. = UFF UFF-type dispersion correction. Parameters for atomic numbers up to 54 are available internally or can be supplied in DISPPR for any atom. Built-in parameters are taken from Rappe et al. J. Am. Chem. Soc. 1992, 114, 10024. = SK The Slater-Kirkwood type dispersion correction omitting the change polarizability depending on the number of bonds. No default values of DISPPR are available. Some are listed in the manual of the DFTB+ program. = SKHP The Slater--Kirkwood type dispersion with the dependence of polarizabilities on the number of bonds. = DFT Use so-called DFT-D. See $DFT for further details. DISP=GRIMME is a synonym. DISPPR an array of parameters used for dispersion correction, listed in sets for each species. For DISP=UFF, DISPPR(1) and DISPPR(2) define the non-bonded distance (Angs.) and energy (kcal/mol) for the first species, respectively, and so on. For DISP=SK, a set for a species has 3 parameters, the polarizability (Angstrom^3), cutoff length (Angstrom), and atomic charge. For DISP=SKHP, a set for a species has 14 parameters. The first six are the polarizabilities depending on the number of bonds, and the next six are cutoff length, and the last is atomic charge. Default: see DISP. HUBDER an array of Hubbard derivatives for each species (1 per species) used only for DFTB3 calculations. Default values are set for the elements included in the 3OB parameters (Br, C, Ca, Cl, F, H, I, K, Mg, N, Na, O, P, S, Zn). MAXANG array of maximum angular momentum of each species, which determines the number of basis functions. DFTB uses only valence orbitals and electrons! Most elements have proper default values, but for some atomic types (i.e., species) you need to manually define the values. QREF array of the number of reference electrons of each species. QREF is usually automatically taken from Slater-Koster parameters, so this option is seldom used. SPNCST an array of spin constants used in unrestricted (UHF) DFTB calculation. Provide 6 spin constants, W_{ss}, W_{sp}, W_{pp}, W_{sd}, W_{pd}, & W_{dd}, for each species in a continuous array. Constants for some elements can be found in the manual of the DFTB+ program. PARAM specifies the directory from which DFTB parameters are taken. If you wish to mix parameters from different directories, this option cannot be used. Specifying PARAM means no $DFTBSK; otherwise, $DFTBSK is read. Nota bene-bene: the actual path for parameters includes $DFTBPAR, defined in rungms. All directory names used in PARAM should be ** UPPER CASE **, as 3OB-3-1 in ~/gamess/dftb/param/3OB-3-1 where $DFTBPAR=~/gamess/auxdata/DFTB PARAM=3OB-3-1 The length of PARAM is maximum 8 characters! Each parameter file name has a limit of 150 characters. GAMESS includes 3OB-3-1 and MATSCI03 (properly called matsci-0-3), which you may specify in PARAM. 3ob-3-1 should be used with DFTB3 (biochemistry+water). matsci-0-3 should be used with DFTB2 (iorganics). You can find more parameter sets at dftb.org. Before using DFTB parameters, ~/gamess/dftb/README.dftb should be consulted regarding lisense and citations. Default: "", meaning that $DFTBSK is read. ISPDMP An array of integer specifying species X to which the X--H damping function (DAMPXH) is applied. By default, with DAMPXH=.TRUE., ISPDMP for all elements is 1 (apply). Setting 1 for H does not do anything. * * * The following options are FMO-DFTB specific (Nishimoto, Y. et al. J. Chem. Theory Comput. 2014, 10, 4801-4812.). FMO-DFTB has many limitations and some FMO options are not supported (for instance, multilayer FMO etc). Only single layer, restricted closed-shell FMO2/3-DFTB1/2/3 are implemented at present. SRSCC, ETEMP etc are not available. The analytic gradient is available for FMO-DFTB, requiring solving SCZV as in other FMO methods. MODESD = controls the behavior of ES-DIM (electrostatic dimer) approximation (bit additive). 1 Calculate interfragment repulsive energy for ES dimers (almost never used). 2 Add up all ES-DIM energies. This means that individual ES dimer energies are not calculated, but only their total lump sum, computed with the dynamic load balancing. 4 Lump ES-DIM routine with static load balancing. The bits of 2 or 4 are mutually exclusive. Default: 0 (i.e., individual ES dimer energies). MODGAM = controls the calculation of gamma values (interatomic 1/R-like function) in FMO-DFTB2 and FMO-DFTB3 (bit additive). 0 Calculate gamma values on the fly. 1 Calculate once and prestore gamma values in triangular matrix. 2 Calculate once and prestore gamma values in square matrix. 4 With the bits of 1 or 2, the calculation of gamma values is parallelized with GDDI. The bits of 1 or 2 are mutually exclusive. These options are faster but takes more memory. 8 Using this option omits computing ESP in dimer and trimer calculations by accumulating contributions of each fragment and subtracting double-counting contributions. Default: 8 * * * The following options are relevant to second- and third-order derivative calculations (RUNTYP=HESSIAN and RAMAN). CPCONV = Convergence criterion during coupled-pertrubed DFTB iterations, similar to CONV in $SCF. In DFTB, the program uses Mulliken charges for testing the convergence, but not density matrix itself. By default, CPCONV=1.0D-06. MXCPIT = Maximum number of coupled-perturbed DFTB iterations. By default, MXCPIT=50. DEGTHR = An array of two degeneracy thresholds. If the difference of two eigenvalues are less than the threshold, two orbitals are seen as degenerated. The first threshold is employed in solving coupled-perturbed equations, while the second threshold is in computing third-order derivatives analytically. By default, these are set to 1.0D-12 and 1.0D-08, whieh are usually reasonable. ARAMAN = A flag to compute third-order derivatives (static hyperpolarizability and polarizability derivative) analytically, in addition to Hessian. If this option is activated, users do not have to give $HESSIAN and $DIPDR in the input, and non-resonance Raman spectra can be simulated by a single run. This option requires that RUNTYP must be RAMAN. By default, ARAMAN=.FALSE. * * * The following optinos are relevant to long-ranged corrected DFTB. The formulation is based on Lutsker, V. et al. 2015, 143, 184107. With LC-DFTB, using ITYPMX=-1 options is highly recommended. LCDFTB = A flag to activate long-range correctios EMU = A parameter for long-range corrections. The meaning is very similar to MU in $DFT. By default, EMU=0.0, and this corresponds to regular DFTB. ICUT = A parameter applied in the screening using the Schwarz inequiality. The meaning is similar to ICUT in $CONTRL. By default, the screening is not employed, but this is usually as fast as ICUT=9, depending on the performance of the math library.

## Evaluating The Friction Of Rotary Joints In Molecular Machines

Howard Lovy's name came across my calendar (happy birthday!) and Josh Hall published a fresh post to "Where is my Flying Car?" – seems like a perfect time to post an everything-new-is-still-new-for-the-first-time-type update.

Published way back in 2017 in Molecular Systems Design & Engineering; reportedly one of the most-read Q3 2017 articles to boot. Also my first foray into arxiv territory for those limited in their journal access (arxiv.org/abs/1701.08202).

Tad Hogg^{*a} Matthew S. Moses^{*b} and Damian G. Allis^{*c}

a. Institute for Molecular Manufacturing, Palo Alto, USA

b. Independent Consultant, Lafayette, USA

c. Department of Chemistry, Syracuse University, Syracuse, USA

Abstract:A computationally-efficient method for evaluating friction in molecular rotary bearings is presented. This method estimates drag from fluctuations in molecular dynamics simulations via the fluctuation–dissipation theorem. This is effective even for simulation times short compared to a bearing's energy damping time and for rotation speeds comparable to or below typical thermal values. We apply this method to two molecular rotary bearings of similar size at 300 K: previously studied nested (9,9)/(14,14) double-walled carbon nanotubes and a hypothetical rotary joint consisting of single acetylenic bonds in a rigid diamondoid housing. The acetylenic joint has a rotational frictional drag coefficient of 2 × 10^-35 kg m2 s^-1. The friction for the nested nanotubes is 120 times larger, comparable to values reported by previous studies. This fluctuation-based method could evaluate dissipation in a variety of molecular systems with similarly rigid and symmetric bearings.

Mol. Syst. Des. Eng., 2017, 2, 235-252 (10.1039/C7ME00021A, direct link)