the following is my script.
from ATK.TwoProbe import *
import numpy
import ATK
ATK.setVerbosityLevel(1)
# Read the atomic configuration from a VNL file
vnl_file = VNLFile("2-G-opt.vnl")
configurations = vnl_file.readAtomicConfigurations()
two_probe_conf = configurations["2-G-opt"]
###############################################
# setting global parameters
###############################################
kpoints = (50,1,30)
mesh_cutoff = 250.*Rydberg
xc = LDA.PZ
temperature = 300.*Kelvin
tolerance = 1e-6
diagonal_mixing_parameter = 0.02
history_steps = 15
max_steps = 1000
basis_set_parameters = [basisSetParameters(SingleZeta,element = Carbon),
basisSetParameters(SingleZeta,element = Nitrogen)]
###############################################
# setting two electrodes parameters
###############################################
electrode_brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters(
monkhorst_pack_parameters = (kpoints)
)
electrode_electron_density_parameters = electronDensityParameters(
mesh_cutoff = mesh_cutoff,
)
electrode_iteration_control_parameters = iterationControlParameters(
tolerance = tolerance,
criterion = IterationControl.Strict,
max_steps = max_steps
)
electrode_iteration_mixing_parameters = iterationMixingParameters(
algorithm = IterationMixing.Pulay,
diagonal_mixing_parameter = diagonal_mixing_parameter,
quantity = IterationMixing.Hamiltonian,
history_steps = history_steps
)
electrode_eigenstate_occupation_parameters = eigenstateOccupationParameters(
temperature = temperature
)
electrode_parameters = ElectrodeParameters(
brillouin_zone_integration_parameters = electrode_brillouin_zone_integration_parameters,
electron_density_parameters = electrode_electron_density_parameters,
eigenstate_occupation_parameters = electrode_eigenstate_occupation_parameters,
iteration_mixing_parameters = electrode_iteration_mixing_parameters,
iteration_control_parameters = electrode_iteration_control_parameters
)
###############################################
# setting scattering region parameters
###############################################
ite_mix_para = iterationMixingParameters(
diagonal_mixing_parameter = diagonal_mixing_parameter,
quantity = IterationMixing.Hamiltonian,
history_steps = history_steps
)
ite_con_para = iterationControlParameters(
tolerance = tolerance,
criterion = IterationControl.Strict,
max_steps = max_steps
)
ele_den_para = electronDensityParameters(
mesh_cutoff = mesh_cutoff,
)
two_probe_algorithm_parameters = twoProbeAlgorithmParameters(
electrode_constraint = ElectrodeConstraints.RealSpaceDensity,
initial_density_type = InitialDensityType.EquivalentBulk
)
energy_contour_integral_parameters = energyContourIntegralParameters(
circle_points = 300,
integral_lower_bound = 30.0*Rydberg,
fermi_line_points = 10,
fermi_function_poles = 4,
real_axis_infinitesimal = 0.01*electronVolt,
real_axis_point_density = 0.02*electronVolt
)
two_center_integral_parameters = twoCenterIntegralParameters(
cutoff = 2500.0*Rydberg,
points = 1024
)
use_old=False
#use_old=True
for voltage0 in numpy.arange(0.00, 6.00+0.50, 0.50): # gate_valtages
voltage = -voltage0
two_probe_method = TwoProbeMethod(
exchange_correlation_type = xc,
electrode_parameters = (electrode_parameters, electrode_parameters), ###############
electrode_voltages = (0.20, -0.20)*Volt, ## parameters##
iteration_mixing_parameters = ite_mix_para, ###############
iteration_control_parameters = ite_con_para,
energy_contour_integral_parameters = energy_contour_integral_parameters,
two_center_integral_parameters = two_center_integral_parameters,
electron_density_parameters = ele_den_para,
basis_set_parameters = basis_set_parameters,
algorithm_parameters = two_probe_algorithm_parameters
)
###############
#voltage = 0 ## parameters##
three_probe_method = GatedTwoProbeMethod( ###############
two_probe_method = two_probe_method, ###############
gate_voltage = voltage*Volt, ## parameters##
surface_atoms = (5, 5) ###############
)
if use_old:
scf = executeSelfConsistentCalculation(
atomic_configuration = two_probe_conf,
method = three_probe_method,
initial_calculation = Old_Scf,
runtime_parameters = runtimeParameters(verbosity_level = 1,
checkpoint_filename = 'GScf-gate-%.2f.nc' % voltage) ## parameters##
)
else:
scf = executeSelfConsistentCalculation(
atomic_configuration = two_probe_conf,
method = three_probe_method,
runtime_parameters = runtimeParameters(verbosity_level = 1,
checkpoint_filename = 'GScf-gate-%.2f.nc' % voltage)
)
use_old=True
Old_Scf=scf
current = calculateCurrent(scf)
print "%.1f\t\t%.2e" %(voltage, current.inUnitsOf(Ampere))