# -------------------------------------------------------------
# Analysis from File
# -------------------------------------------------------------
configuration = nlread('Au-C6H4S2-Au.nc', object_id='gID000')[0]
# -------------------------------------------------------------
# Total Energy
# -------------------------------------------------------------
total_energy = TotalEnergy(configuration)
nlsave('analysis.nc', total_energy)
nlprint(total_energy)
Next, check that the Object id corresponds to the self-consistent calculation you want to base the analysis on. If your calculation was a "standard run", then most likely the default gID001 is the one you want. However, if you did a geometry optimization, gID000 corresponds to the initial geometry, and the optimized geometry will be gID001. You can inspect the contents of the NetCDF file, and the object ids, in the Result Browser in the main VNL job.
Note: nlread() can in principle read many objects at once from a NetCDF file. Object ids must however be unique, so by using the object_id keyword, you are guaranteed to get only one object back. However, nlread() always returns a list of objects, even if there is only one match, so we need the [0] to pick out the first/only configuration matching the object_id.
configuration = nlread('Au-C6H4S2-Au.nc', object_id='gID000')[0]
# -------------------------------------------------------------
# Two-probe Configuration
# -------------------------------------------------------------
# -------------------------------------------------------------
# Left Electrode
# -------------------------------------------------------------
# Set up lattice
vector_a = [8.65127, 0.0, 0.0]*Angstrom
vector_b = [-4.32564, 7.49222, 0.0]*Angstrom
vector_c = [0.0, 0.0, 7.06373620597]*Angstrom
left_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
left_electrode_elements = [...
...
...]
# Define coordinates
left_electrode_coordinates = [...
...
...]*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 = [8.65127, 0.0, 0.0]*Angstrom
vector_b = [-4.32564, 7.49222, 0.0]*Angstrom
vector_c = [0.0, 0.0, 7.06373620597]*Angstrom
right_electrode_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
right_electrode_elements = [...
...
...]
# Define coordinates
right_electrode_coordinates = [...
...
...]*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 = [8.65127, 0.0, 0.0]*Angstrom
vector_b = [-4.32564, 7.49222, 0.0]*Angstrom
vector_c = [0.0, 0.0, 26.1041296975]*Angstrom
central_region_lattice = UnitCell(vector_a, vector_b, vector_c)
# Define elements
central_region_elements = [...
...
...]
# Define coordinates
central_region_coordinates = [...
...
...]
# 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]
)
nlsave('Au-C6H4S2-Au.nc', device_configuration) # the 1st one gID000
nlprint(device_configuration)
# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
#----------------------------------------
# Basis Set
#----------------------------------------
basis_set = [
GGABasis.Hydrogen_DoubleZetaPolarized(element=Hydrogen,filter_mesh_cutoff=200*Rydberg),
GGABasis.Carbon_DoubleZetaPolarized(element=Carbon,filter_mesh_cutoff=200*Rydberg),
GGABasis.Sulfur_DoubleZetaPolarized(element=Sulfur,filter_mesh_cutoff=200*Rydberg),
GGABasis.Gold_SingleZetaPolarized(element=Gold,filter_mesh_cutoff=200*Rydberg),
]
#----------------------------------------
# Exchange-Correlation
#----------------------------------------
exchange_correlation = GGA.PBE
#----------------------------------------
# Numerical Accuracy Settings
#----------------------------------------
left_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(6, 6, 20),
density_mesh_cutoff=200.0*Rydberg,
)
right_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(6, 6, 20),
density_mesh_cutoff=200.0*Rydberg,
)
device_numerical_accuracy_parameters = NumericalAccuracyParameters(
k_point_sampling=(6, 6, 20),
density_mesh_cutoff=200.0*Rydberg,
)
#----------------------------------------
# Iteration Control Settings
#----------------------------------------
left_electrode_iteration_control_parameters = IterationControlParameters(
max_steps=500,
)
right_electrode_iteration_control_parameters = IterationControlParameters(
max_steps=500,
)
device_iteration_control_parameters = IterationControlParameters(
max_steps=500,
)
#----------------------------------------
# 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]]
)
#----------------------------------------
# Contour Integral Settings
#----------------------------------------
contour_parameters = DoubleContourIntegralParameters(
integral_lower_bound=50.0*eV,
)
#----------------------------------------
# Electrode Calculators
#----------------------------------------
left_electrode_calculator = LCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=left_electrode_numerical_accuracy_parameters,
iteration_control_parameters=left_electrode_iteration_control_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,
iteration_control_parameters=right_electrode_iteration_control_parameters,
poisson_solver=right_electrode_poisson_solver,
)
#----------------------------------------
# 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,
contour_parameters=contour_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
device_configuration.setCalculator(calculator)
nlprint(device_configuration)
device_configuration.update()
nlsave('Au-C6H4S2-Au.nc', device_configuration) # the 2nd one gID001
nlprint( )
# -------------------------------------------------------------
# Transmission Spectrum
# -------------------------------------------------------------
transmission_spectrum = TransmissionSpectrum(
configuration=device_configuration,
energies=numpy.linspace(-2,2,400)*eV,
kpoints=MonkhorstPackGrid(6,6),
energy_zero_parameter=AverageFermiLevel,
infinitesimal=1.36057e-05*eV,
self_energy_calculator=RecursionSelfEnergy(),
)
nlsave('Au-C6H4S2-Au.nc', transmission_spectrum) # the 3rd one gID002
nlprint(transmission_spectrum)
# -------------------------------------------------------------
# Molecular Energy Spectrum
# -------------------------------------------------------------
molecular_energy_spectrum = MolecularEnergySpectrum(
configuration=device_configuration,
energy_zero_parameter=FermiLevel,
projection_list=ProjectionList(elements=[Carbon, Hydrogen, Sulfur])
)
nlsave('Au-C6H4S2-Au.nc', molecular_energy_spectrum) # the 4th one gID003
nlprint(molecular_energy_spectrum)
# -------------------------------------------------------------
# Eigenstate
# -------------------------------------------------------------
eigenstate = Eigenstate(
configuration=device_configuration,
projection_list=ProjectionList(elements=[Carbon, Hydrogen, Sulfur]),
quantum_number=19,
)
nlsave('Au-C6H4S2-Au.nc', eigenstate) # the 5th one gID004
nlprint( )
# -------------------------------------------------------------
# Eigenstate
# -------------------------------------------------------------
eigenstate = Eigenstate(
configuration=device_configuration,
projection_list=ProjectionList(elements=[Carbon, Hydrogen, Sulfur]),
quantum_number=20,
)
nlsave('Au-C6H4S2-Au.nc', eigenstate) # the 6th one gID005
nlprint( )
# -------------------------------------------------------------
# Analysis from File
# -------------------------------------------------------------
configuration = nlread('Au-C6H4S2-Au.nc', object_id='gID001')[0]
# -------------------------------------------------------------
# Total Energy
# -------------------------------------------------------------
total_energy = TotalEnergy(configuration)
nlsave('analysis.nc', total_energy)
nlprint(total_energy)
It simply depends on whether you read the whole file, or parts of it. The following points should be noted:
#----------------------------------------
# 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,
contour_parameters=contour_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
device_configuration.setCalculator(calculator)
nlprint(device_configuration)
device_configuration.update()
nlsave('Au-C6H4S2-Au-0.0.nc', device_configuration)
#----------------------------------------
# Device Calculator
#----------------------------------------
calculator = DeviceLCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=device_numerical_accuracy_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
electrode_voltages=( 0.0*Volt, 1.0*Volt)
)
device_configuration.setCalculator(calculator)
# -------------------------------------------------------------
# Initial State
# -------------------------------------------------------------
old_calculation = nlread('Au-C6H4S2-Au-0.0.nc', DeviceConfiguration)[0]
device_configuration.setCalculator(
calculator,
initial_state=old_calculation,
)
device_configuration.update()
nlsave('Au-C6H4S2-Au-1.0.nc', device_configuration)
nlprint(device_configuration)
# -------------------------------------------------------------
# Initial State
# -------------------------------------------------------------
old_calculation = nlread('Au-C6H4S2-Au-0.0.nc', DeviceConfiguration)[0]
device_configuration.setCalculator(
calculator(electrode_voltages=( 0.0*Volt, 1.0*Volt)),
initial_state=old_calculation,
)
device_configuration.update()
nlsave('Au-C6H4S2-Au-1.0.nc', device_configuration)
nlprint(device_configuration)
With respect to electrode_voltages=( 0.0*Volt, 1.0*Volt), how would ATK do for the shift of chemicalpotentials of the electrodes relative to the equilibrium ones? Will ATK only shift the right electrode down by 1.0 eV? Or instead, ATK would shift down the right electrode by 0.5 eV and shift up the left electrode by 0.5 eV. I think the latter technique would be favorable to the convergence.
old_calculation = nlread('Au-C6H4S2-Au-0.0.nc', DeviceConfiguration)[0]
device_configuration.setCalculator(
calculator,
initial_state=old_calculation,
)
device_configuration.update()
nlsave('Au-C6H4S2-Au-1.0.nc', device_configuration)
nlprint(device_configuration)
device_configuration = nlread('Au-C6H4S2-Au-0.0-checkpoint.nc')[0]
device_configuration.update(force_restart=True)
nlsave('Au-C6H4S2-Au-1.0.nc', device_configuration)
nlprint(device_configuration)
old_calculation = nlread('Au-C6H4S2-Au-0.0.nc', DeviceConfiguration)[0]
device_configuration.setCalculator(
calculator,
initial_state=old_calculation,
)
device_configuration.update()
nlsave('Au-C6H4S2-Au-1.0.nc', device_configuration)
nlprint(device_configuration)
#----------------------------------------
# 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,
contour_parameters=contour_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
device_configuration.setCalculator(calculator)
device_configuration = nlread("Au-C6H4S2-Au.nc")[1]
calculator=calculator(
electrode_voltages=(0.0*Volt, 0.0*Volt))
device_configuration.setCalculator(
calculator(),
initial_state=device_configuration)
device_configuration.update()
nlsave("Au-C6H4S2-Au-0.0.nc", device_configuration)
device_configuration.update()
nlsave('Au-C6H4S2-Au.nc', device_configuration)
#----------------------------------------
# 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,
contour_parameters=contour_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
device_configuration.setCalculator(calculator)
device_configuration = nlread("checkpoint-0.4.nc")[0]
device_configuration.update(force_restart=True)
nlsave("Au-C6H4S2-Au-restart.nc", device_configuration)
transmission_spectrum = TransmissionSpectrum(
configuration=device_configuration,
energies=numpy.linspace(-2,2,100)*eV,
kpoints=MonkhorstPackGrid(6,6),
energy_zero_parameter=AverageFermiLevel,
infinitesimal=1.36057e-05*eV,
self_energy_calculator=RecursionSelfEnergy(),
)
nlsave("Au-C6H4S2-Au-restart.nc", transmission_spectrum)
nlprint(transmission_spectrum)
transmission_spectrum.current()
#----------------------------------------
# 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,
contour_parameters=contour_parameters,
electrode_calculators=
[left_electrode_calculator, right_electrode_calculator],
)
device_configuration.setCalculator(calculator)
old_calculation = nlread("Au-C6H4S2-Au-0.4.nc", DeviceConfiguration)[0]
device_configuration.setCalculator(
calculator(electrode_voltages=(0.2*Volt, -0.2*Volt)),
initial_state=old_calculation,
)
device_configuration.update()
nlsave("Au-C6H4S2-Au-initial2.nc", device_configuration)
transmission_spectrum = TransmissionSpectrum(
configuration=device_configuration,
energies=numpy.linspace(-2,2,400)*eV,
kpoints=MonkhorstPackGrid(6,6),
energy_zero_parameter=AverageFermiLevel,
infinitesimal=1.36057e-05*eV,
self_energy_calculator=RecursionSelfEnergy(),
)
nlsave("Au-C6H4S2-Au-initial2.nc", transmission_spectrum)
nlprint(transmission_spectrum)
transmission_spectrum.current()
old_calculation = nlread("Au-C6H4S2-Au-0.4.nc", DeviceConfiguration)[0]
old_calculation = nlread("Au-C6H4S2-Au-0.4.nc")[0]
device_configuration = nlread("file.nc")[-1]
transmission_spectrum = TransmissionSpectrum(
configuration=device_configuration,
....,
)
nlsave("file.nc", transmission_spectrum)
nlprint(transmission_spectrum)
print transmission_spectrum.current()
device_configuration = nlread("checkpoint-0.4.nc")[0]
device_configuration.setCalculator(device_configuration.calculator(), initial_state=device_configuration)
device_configuration.update(force_restart=True)
nlsave("Au-C6H4S2-Au-restart.nc", device_configuration)
Just to clarify: force_restart=True does force the calculation to restart from scratch. But since you use the state from the checkpoint file as initial state, "from scratch" here means from that initial guess rather than from a superposition of atomic densities :)
device_configuration = nlread("checkpoint-0.4.nc")[0]
device_configuration.setCalculator(device_configuration.calculator(), initial_state=device_configuration)
device_configuration.update()
nlsave("Au-C6H4S2-Au-restart.nc", device_configuration)
(2) And I am puzzling about the force_restart=True in device_configuration.update(), that is what is the difference between device_configuration.update() and device_configuration.update(force_restart=True) since they both trigger a calculation?
device_configuration = nlread("checkpoint-0.4.nc")[0]
device_configuration.setCalculator(device_configuration.calculator(), initial_state=device_configuration)
device_configuration.update()
nlsave("Au-C6H4S2-Au-restart.nc", device_configuration)
device_configuration = nlread("checkpoint-0.4.nc")[0]
device_configuration.setCalculator(device_configuration.calculator(), initial_state=device_configuration)
device_configuration.update(force_restart=True)
nlsave("Au-C6H4S2-Au-restart.nc", device_configuration)
configuration = nlread('GaAs.nc', BulkConfiguration)[-1]
configuration.update()
configuration = nlread('GaAs.nc', BulkConfiguration)[-1]
configuration.setCalculator(configuration.calculator(), initial_state=configuration)
configuration.update()