Dear all,
I have learning on how to restart/resore or even do the post analysis based on a converged *.nc file in ATK-2015.1.
I found here the tutorial very useful for the psot analysis based on a converged *.nc file
http://www.quantumwise.com/publications/tutorials/item/504-computing-quantities-from-converged-calculations
The trick to use a converged *.nc file doing a post analysis with out any futher SCFs is to use
nlread, eg., the following is to do a total energy calculation based on the already converged Au-C6H4S2-Au.nc file.
# -------------------------------------------------------------
# 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.
However, I think there maybe something worth to emphasize. That is,
configuration = nlread('Au-C6H4S2-Au.nc', object_id='gID000')[0]
we should specify the object_id properly and may not be the
gID000 one. It depends on your number of
nlsave used in the job generating this *.nc file and also the quantities you want to base for the post analysis, usually the object_id contains the converged electron density not others would be expected to used for the post analysis.
The number of object_id in the *.nc file depends on the
nlsave used in the job generating this *.nc. So, usually and most probably the first object_id that is
gID000 is not the one contains the converged electron density.
For example, the following *.py input file contains 6
nlsave sentances and the second object_id that is
gID001 is what we need to do the total energy post analysis.
# -------------------------------------------------------------
# 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( )
The first object_id, i.e.,
gID000 only contain only the configuration (lattice vectors and atom positions) of the two probe system.
The converged electron density is stored in the second object_id, i.e.,
gID001.
So, the post analysis input file should read
# -------------------------------------------------------------
# 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)