Titanium chain; lattice constant

[Hasan Sahin]

By using the examples for Li chain optimization in ATKmanual, I tried to calculate lattice constant of Titanium linear chain.
Lattice constant must be around 2.3 Angstrom for linear chains but my result is around 0.7 Angstrom. (see attachment)
What is wrong in this simple calculation ?



my script:


# Set the basis set to be SZ
basis_set = basisSetParameters(type=DoubleZetaDoublePolarized)

# Define the density grid cut-off
electron_density_params = electronDensityParameters(
    mesh_cutoff=150*Rydberg
    )

# Define the k-point sampling
bz_integration = brillouinZoneIntegrationParameters( (1,1,100) )

# Define the DFT method
dft_method = KohnShamMethod(
    basis_set_parameters=basis_set,
    brillouin_zone_integration_parameters=bz_integration,
    electron_density_parameters=electron_density_params
    )

# Open file for storing the results
f = open('lattice.dat', 'w')

# Loop over lattice constants
import numpy
for c in numpy.arange(0.5,5.5,0.1):

    # Create the unit cell for 1D chain
    unit_cell = [[10.0,  0.0, 0.0],
                 [ 0.0, 10.0, 0.0],
                 [ 0.0,  0.0,  c ]]*Angstrom

    # Create the atomic configuration
    ti_chain = PeriodicAtomConfiguration(unit_cell,[Titanium],
                                         [ (0.0, 0.0, 0.0)*Angstrom ] )
   
    # Execute the self-consistent calculation
    dft_calculation = dft_method.apply(ti_chain)
   
    # Calculate the system energy
    energy = calculateTotalEnergy(dft_calculation)
   
    # Write lattice constant and energy to file and screen
    s = '%.4f %.10f \n' % (c, energy.inUnitsOf(Rydberg))
    f.write(s)
    print s,

# Close file and exit
f.close()

[Nordland]

Hey Hasan Sahin!

There is nothing wrong with your script, but I have redone the calculation with GGA,
and I get pretty good results, so I suppose GGA is better for Titanium, but properly the
default pseudopotential for LDA Titanium is not the best for wire calculations.
I also made the calculation spin-polarized, but it is not needed to get good results,
I just like to include spin always

I have attached my results and a slightly modified version of your script - where I have improved the
performance slightly by restarting using the old density.

I hope you can use this

[Hasan Sahin]

many thanks

[Anders Blom]

The obvious suspicion in a case like this is indeed the pseudopotential. As a test, I ran bulk Ti (Mg-type hexagonal, Strukturbericht A3, hcp) for both LDA and GGA. In both cases there is a well-defined minimum in roughly the correct position, close to the experimental lattice constant (I only tested for the value of a). Often, there isn't a good minimum if the pseudopotential is bad, so at least the potentials are ok from this perspective.

The GGA mimimum is however clearly closer to the experimental minimum, so it seems that's your best choice for several reasons.

Regarding the linear chain script, I would myself (as Nordland also did), increase the XY cell size a bit; especially when you get up to 5 Å spacing, having 10 Å in the other two directions is a bit small, although of course around 3 Å, where it counts, it's not a big problem. But 15 Å is probably a better spacing, especially if you imagine using this cell later on in two-probe calculations. Linear chains are very bad at screening the electrostatic potential from the repeated images of the system in the XY directions.

Also, and especially if you stick to GGA, I would increase the mesh cut-off or, rather, make a series of calculations to figure out which mesh cut-off to use for Ti.

I attach the script I used for the bulk lattice constant test, in case someone finds it interesting. It also uses the trick mentioned by Nordland to take the converged density of the previous calculation as the initial guess for the next one, to save time; it often reduces the number of iterations down to 2-3 in a sweep like this.

Also attached is a plot (and the script to make the plot) of the energy vs lattice constant "a" keeping "c" fixed at 4.684 Å. The experiemental value of "a" is 2.951 Å.