Parameters and parametrization

• Hotbit as such does not contain that many parametrizations
• Some DFTB parametrizations are available at [http://www.dftb.org]. These tables are compatible with Hotbit, but a seperate license agreement with the maintainers of these tables is required in order to use them.
• Some first-trial -parametrizations exist (in folder $HOTBIT_PARAMETERS\inofficial)
  • Do not use these for publication-quality data
• In general,
  • there are no unique parametrizations
  • never construct, use and publish results with parametrizations, unless you know what you are doing
  • if you construct new parametrizations, please deliver them to hotbit and general use

Fitting process

• Units:
  1. Bohr and Hartree in atom calculations (ASE interface not even needed)
  2. Bohr and Hartree in Slater-Koster table compilation (ASE interface not even needed)
  3. eV and Ångström in repulsion fitting (ASE interface is needed)
• Couple of examples of fitting come with the code (try shell guide; hotbit -e).
• Schematic of the fitting process:

1. Calculate free atom with DFT

• To get free-atom on-site energies for valence electrons
• LDA (Perder-Wang 92) calculation for atom without additional confinements
• Used to make the .elm file.

    from hotbit.parametrization.atom import KSAllElectron
    from hotbit.parametrization.util import IP_EA
    from time import asctime

    """ 
    Create the element file for C
    - contains the single-particle energy levels
    - Hubbard parameter U
    - FWHM-value (calculated from U)
    - possibly comments about the parametrization
    """

    # element to parametrize
    element    = 'C'
    # lowest orbital that can accept electrons
    add_orb    = '2p'
    # sometimes adding full electron causes the convergence problems,
    # so add a fraction of an electron and interpolate
    add        =  0.67
    # highest occupied orbital
    remove_orb = '2p'
    # how many electrons to remove
    remove     =  1

    # solve the atom with charges add, 0 and remove,
    # and calculate ionization potential and electron affinity
    IP, EA, atom = IP_EA(element, remove_orb, add_orb, remove, add)

    # the Hubbard parameter is then
    U = IP-EA

    # save the results into an element file
    f = open("C.elm",'w')
    print >> f, "symbol=%s" % element
    print >> f, "comment="
    print >> f, asctime()
    print >> f, "Energy levels from DFT"
    # an empty line after the comments
    print >> f, ""
    for va in atom.get_valence_energies():
        print >> f, "epsilon_%s= %0.6f" % va
    print >> f, "U= %0.6f" % U
    print >> f, "FWHM= %0.6f" % (1.32856/U)
    f.close()

    # save the radial parts of the wave functions
    atom.write_unl("C.elm")

2. Calculate confined atom with DFT

• To get wave functions for the atoms (in molecules)
• Compression mimics the fact, that while atoms are bound to molecules or solids, their wave functions are more compressed and not as diffuse as for a free atom.
• Use quadratic external potential V_ext(r)=-Z/r + (r/!r0)**2 (Frauenheim type confinement)
• r0 is the first parameter in the parametrization process; rule of thumb r0=1.85*(covalent radius)

    from hotbit.parametrization.atom import KSAllElectron
    from box.data import data
    Bohr = 0.52917725750691647

    # The element parameters, r0 must be in Bohrs
    elem1 = 'H'
    r0_cov_1 = data[elem1]['R_cov']/Bohr
    r0_factor_1 = 1.85
    r0_1 = r0_cov_1 * r0_factor_1

    elem2 = 'C'
    r0_cov_2 = data[elem2]['R_cov']/Bohr
    r0_factor_2 = 1.85
    r0_2 = r0_cov_2 * r0_factor_2

    # Calculate wave functions for the confined isolated atoms
    atom1 = KSAllElectron(elem1, confinement={'mode':'quadratic', 'r0':r0_1})
    atom1.run()

    atom2 = KSAllElectron(elem2, confinement={'mode':'quadratic', 'r0':r0_2})
    atom2.run()

    # Produce png-files from the radial parts of the wave functions
    atom1.plot_Rnl()
    atom2.plot_Rnl()

3. Compile Slater-Koster tables

• To get Hamiltonian and overlap matrix elements for elementary two-center integrals as a function of distance
• Use wave functions from previous calculation.
• Output is the {{{.par}}} file without short range repulsion.

    # Add this to the previous script
    from hotbit.parametrization.slako import SlaterKosterTable

    # Slater-Koster table limits and the number of points
    R1 = 1
    R2 = 10
    N = 50

    # The output file
    file = '%s_%s_no_repulsion.par' % (elem1, elem2)

    # Create the Slater-Koster table
    sk = SlaterKosterTable(atom1, atom2)
    sk.run(R1, R2, N)
    sk.write(file)
    sk.plot()

4. Fit the repulsive potential

• To get the simple short range repulsive potential (=repulsion of ionic cores + all possible corrections)
• Before repulsion fitting, the electronic structure (Slater-Koster tables, on-site energies, Hubbard-U, FWHM,...) should be completely fixed. (Otherwise, if electronic parameters change, you need to fit the repulsion all over again.)
• Use bunch of DFT-optimized structures; in principle E_DFT - E_(band structure + coulomb) = E_(repulsion). Use this information to construct the derivative of V_rep(R) (Collect points to (R,V_rep'(R)), fit a smooth polynomial to this data, and integrate to get V_rep(R)).
• This repulsion is added into the .par-file containing the Slater-Koster tables.

    from ase import *
    from hotbit import *
    from hotbit.parametrization.fitting import RepulsiveFitting

    # force the hotbit to use CC-table where the repulsion is missing
    tables = {'CC':'C_C_no_repulsion.par', 'others':'default'}

    # for parametrization calculations it might be good idea to increase
    # the convergence criteria
    mixer_values = {'name':'anderson', 'convergence':1e-9}

    # the fitting object (element1, element2, r_cut)
    rep=RepulsiveFitting('C','C', 2.0)
    
    # now use various methods if RepulsiveFitting to collect information about V_rep'(R)
    # For each fitting method one must provide a suitable tight-binding calculator,
    # i.e., one that has correct charge and other parameters.

    # for dimer append_energy_curve is suitable
    calc = Hotbit(mixer=mixer_values, tables=tables, txt='-')
    rep.append_energy_curve(1, calc, 'dimer.traj', label='dimer',color='yellow')
                        
    # append energy curve can be used if there are many bond lengths of the
    # same lengths, and the rest can be ignored (e.g. chain)
    rep.append_energy_curve(1, calc, 'chain.traj', label='chain')

    # if there are different bond lengths, then append_homogeneous_structure
    # can be used. If the reference (DFT) calculations was charged,
    # then another calculator must be created with correct charge.
    calc = Hotbit(mixer=mixer_values, tables=tables, txt='-', charge=-1)
    rep.append_homogeneous_cluster(0.5, calc, 'homog.traj',label='homogeneous', color='green')

    # The append_energy_curve method also works with periodic calculations.
    # ASE-trajectory contains the pbc-info, but for now one must disable
    # the SCC with bulk calculations.
    calc_bulk = Hotbit(mixer=mixer_values, tables=tables, txt='-', SCC=False)
    rep.append_energy_curve(1, calc_bulk, 'bulk.traj', label='bulk')

    # make the fit
    rep.fit()

    # write final parametrization file
    rep.write_par("C_C.par")

    # check if the fitting looks reasonable (produces a pdf)
    rep.plot()

About parametrization tuning

1. r0

• 1.85*covalent radius
• Things to look at:
  1. Bulk band structure
  2. Energy curves (e.g. dimer) beyond cutoff

2. Hubbard-U

• U=(ionization energy)-(electron affinity)=IE-EA
• Things to look at:
  1. Charge transfer for heteronuclear systems
  2. Ionization energies for clusters
  3. Linear response (LR-TD-DFTB)

3. FWHM

• U=1.32856/FWHM (from on-site interaction)
• Things to look at: (as with hubbard-U)

4. Repulsive fitting:

• Notes & things to look at:
  1. Repulsion fitting affects only geometry and energy. Electronic structure should already be tuned.
  2. smooting parameter s: larger s, smoother V_rep'(R) -curve
  3. Repulsion cutoff r_cut: affects also absolute energies
  4. If data points in V_rep'(R) are scattered, try improving electronic structure