2
« on: October 19, 2011, 13:06 »
I'm still deciding if I should purchase ATK or not. As one of the tests I'm trying to reproduce the experimental conductivity of simple metals.
I have been trying without success to reproduce the experimental IV curve for gold. Invariably, the calculates conductivity is more than 10 times smaller than the experimental one, ~4.5E7 S/m. I would expect it to be higher, since the system is defect free and single crystal.
Have anyone had success reproducing the experimental conductivity of simple metals? If so please tell me the recipe.
I used the script below:
# -------------------------------------------------------------
# TwoProbe configuration
# -------------------------------------------------------------
# -------------------------------------------------------------
# Left electrode
# -------------------------------------------------------------
# Set up lattice
vector_a = [5.02259, 0.0, 0.0]*Angstrom
vector_b = [1.79916535203e-16, 2.93826, 0.0]*Angstrom
vector_c = [0.0, 0.0, 14.77686]*Angstrom
left_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
left_electrode_elements = [Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold,
Gold, Gold]
# Define coordinates
left_electrode_coordinates = [[ 0.41676981, 0.7345657 , 1.23139213],
[ 2.92806367, 2.20369422, 1.23139396],
[ 1.25573221, 2.2036943 , 3.69415997],
[ 3.76683931, 0.7345657 , 3.69426439],
[ 2.09471311, 0.7345657 , 6.15693373],
[ 4.60582019, 2.20369429, 6.15703787],
[ 0.41676981, 0.7345657 , 8.61982213],
[ 2.92806367, 2.20369422, 8.61982396],
[ 1.25573221, 2.2036943 , 11.08258997],
[ 3.76683931, 0.7345657 , 11.08269439],
[ 2.09471311, 0.7345657 , 13.54536373],
[ 4.60582019, 2.20369429, 13.54546787]]*Angstrom
# Set up configuration
left_electrode = BulkConfiguration(
bravais_lattice=left_electrode_lattice,
elements=left_electrode_elements,
cartesian_coordinates=left_electrode_coordinates
)
# -------------------------------------------------------------
# Right electrode
# -------------------------------------------------------------
# Set up lattice
vector_a = [5.02259, 0.0, 0.0]*Angstrom
vector_b = [1.79916535203e-16, 2.93826, 0.0]*Angstrom
vector_c = [0.0, 0.0, 14.77686]*Angstrom
right_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
right_electrode_elements = [Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold,
Gold, Gold]
# Define coordinates
right_electrode_coordinates = [[ 0.41676981, 0.7345657 , 1.23139213],
[ 2.92806367, 2.20369422, 1.23139396],
[ 1.25573221, 2.2036943 , 3.69415997],
[ 3.76683931, 0.7345657 , 3.69426439],
[ 2.09471311, 0.7345657 , 6.15693373],
[ 4.60582019, 2.20369429, 6.15703787],
[ 0.41676981, 0.7345657 , 8.61982213],
[ 2.92806367, 2.20369422, 8.61982396],
[ 1.25573221, 2.2036943 , 11.08258997],
[ 3.76683931, 0.7345657 , 11.08269439],
[ 2.09471311, 0.7345657 , 13.54536373],
[ 4.60582019, 2.20369429, 13.54546787]]*Angstrom
# Set up configuration
right_electrode = BulkConfiguration(
bravais_lattice=right_electrode_lattice,
elements=right_electrode_elements,
cartesian_coordinates=right_electrode_coordinates
)
# -------------------------------------------------------------
# Central region
# -------------------------------------------------------------
# Set up lattice
vector_a = [5.02259, 0.0, 0.0]*Angstrom
vector_b = [1.79916535203e-16, 2.93826, 0.0]*Angstrom
vector_c = [2.26205428756e-15, 2.26205428756e-15, 36.94215]*Angstrom
central_region_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
central_region_elements = [Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold,
Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold,
Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold, Gold]
# Define coordinates
central_region_coordinates = [[ 0.41676981, 0.7345657 , 1.23139213],
[ 2.92806367, 2.20369422, 1.23139396],
[ 1.25573221, 2.2036943 , 3.69415997],
[ 3.76683931, 0.7345657 , 3.69426439],
[ 2.09471311, 0.7345657 , 6.15693373],
[ 4.60582019, 2.20369429, 6.15703787],
[ 0.41676981, 0.7345657 , 8.61982213],
[ 2.92806367, 2.20369422, 8.61982396],
[ 1.25573221, 2.2036943 , 11.08258997],
[ 3.76683931, 0.7345657 , 11.08269439],
[ 2.09471311, 0.7345657 , 13.54536373],
[ 4.60582019, 2.20369429, 13.54546787],
[ 0.41676981, 0.7345657 , 16.00825213],
[ 2.92806367, 2.20369422, 16.00825396],
[ 1.25573221, 2.2036943 , 18.47101997],
[ 3.76683931, 0.7345657 , 18.47112439],
[ 2.09471311, 0.7345657 , 20.93379373],
[ 4.60582019, 2.20369429, 20.93389787],
[ 0.41676981, 0.7345657 , 23.39668213],
[ 2.92806367, 2.20369422, 23.39668396],
[ 1.25573221, 2.2036943 , 25.85944997],
[ 3.76683931, 0.7345657 , 25.85955439],
[ 2.09471311, 0.7345657 , 28.32222373],
[ 4.60582019, 2.20369429, 28.32232787],
[ 0.41676981, 0.7345657 , 30.78511213],
[ 2.92806367, 2.20369422, 30.78511396],
[ 1.25573221, 2.2036943 , 33.24787997],
[ 3.76683931, 0.7345657 , 33.24798439],
[ 2.09471311, 0.7345657 , 35.71065373],
[ 4.60582019, 2.20369429, 35.71075787]]*Angstrom
# Set up configuration
central_region = BulkConfiguration(
bravais_lattice=central_region_lattice,
elements=central_region_elements,
cartesian_coordinates=central_region_coordinates
)
device_configuration = DeviceConfiguration(
central_region,
[left_electrode, right_electrode]
)
nlprint(device_configuration)
# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
basis_set = GGABasis.DoubleZetaPolarized
#----------------------------------------
# Exchange-Correlation
#----------------------------------------
exchange_correlation = GGA.PBE
#----------------------------------------
# Numerical Accuracy Settings
#----------------------------------------
left_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(3, 5, 150),
)
right_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(3, 5, 150),
)
device_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(3, 5, 1),
)
#----------------------------------------
# Iteration Control Settings
#----------------------------------------
device_iteration_control_parameters = IterationControlParameters(
tolerance=1e-07,
)
#----------------------------------------
# Poisson Solver Settings
#----------------------------------------
left_electrode_poisson_solver = FastFourier2DSolver(
boundary_conditions=[[PeriodicBoundaryCondition,PeriodicBoundaryCondition],
[PeriodicBoundaryCondition,PeriodicBoundaryCondition],
[PeriodicBoundaryCondition,PeriodicBoundaryCondition]]
)
right_electrode_poisson_solver = FastFourier2DSolver(
boundary_conditions=[[PeriodicBoundaryCondition,PeriodicBoundaryCondition],
[PeriodicBoundaryCondition,PeriodicBoundaryCondition],
[PeriodicBoundaryCondition,PeriodicBoundaryCondition]]
)
#----------------------------------------
# Electrode Calculators
#----------------------------------------
left_electrode_calculator = LCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=left_electrode_numerical_accuracy_parameters,
poisson_solver=left_electrode_poisson_solver,
)
right_electrode_calculator = LCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=right_electrode_numerical_accuracy_parameters,
poisson_solver=right_electrode_poisson_solver,
)
# Define output NetCDF file
scf_filename = 'Au-bulk_iv_scf_analysis-GGA2.nc'
#----------------------------------------
# Device Calculator
#----------------------------------------
calculator = DeviceLCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=device_numerical_accuracy_parameters,
iteration_control_parameters=device_iteration_control_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
nlprint(device_configuration)
for voltage in [0., 0.1, 0.2, 0.3, 0.4, 0.5]*Volt:
device_configuration.setCalculator(
calculator(electrode_voltages=(-0.5*voltage,0.5*voltage)),
initial_state=device_configuration
)
device_configuration.update()
nlsave(scf_filename, device_configuration)
#
# Calculate the transmission spectra
#
# Define output NetCDF file
analysis_filename = 'Au-bulk_iv_scf_analysis-GGA2.nc'
# -------------------------------------------------------------
# Transmission spectrum
# -------------------------------------------------------------
# Read all self-consistent calculations from NetCDF file
configurations = nlread(scf_filename, DeviceConfiguration)
for configuration in configurations:
# For each one, extract the bias,
calculator = configuration.calculator()
bias = calculator.electrodeVoltages()[1]-calculator.electrodeVoltages()[0]
# ... calculate and save the transmission spectrum,
transmission_spectrum = TransmissionSpectrum(
configuration=configuration,
energies=numpy.linspace(-2,2,101)*eV,
kpoints=MonkhorstPackGrid(9,15),
energy_zero_parameter=AverageFermiLevel,
infinitesimal=1e-06*eV,
self_energy_calculator=KrylovSelfEnergy(),
)
nlsave(analysis_filename, transmission_spectrum, object_id="Transmission %s" % bias)
potential = EffectivePotential(configuration)
# Calculate and save the voltage drop (except for zero bias, of course)
if float(bias)!=0.:
voltage_drop = potential - zero_bias_potential
nlsave(analysis_filename, voltage_drop, object_id="Voltage drop %s" % bias)
else:
zero_bias_potential = potential
# Copy geometry to analysis file, for plotting
zero_bias_calculation = nlread(scf_filename, DeviceConfiguration, read_state = False)[0]
nlsave(analysis_filename, zero_bias_calculation)