HybridCalculator class

The HybridCalculator provides a regular ASE calculator interface, but internally combines results from multiple different calculators and interactions between them. Currently you can use it for potential energy and force calculations (stress calculations not yet implemented). Like with any ASE calculator, you can also run molecular dynamics and local geometry optimizations with it, but global energy optimization schemes that require a stress implementation will naturally not work due to missing stress implementation.

The following figure demonstrates the simplified workflow comparison of a regular ASE calculation vs. Pysic’s hybrid QM/MM calculation.

HybridCalculator is meant to be used in conjunction with the classes SubSystem and Interaction. See also the run example.

Mechanical Embedding

The electrostatic interaction between subsystems is treated with the mechanical embedding scheme. The Coulomb interaction is thus calculated in the MM level by assigning point charges for the atoms in the primary and secondary system. The charges in the secondary system are static and can be obtained from experimental data or theoretical calculations. Rather than also parametrizing the charges in the primary system, in Pysic they can be dynamically calculated by using the electron density provided by the DFT calculator. The positive charge of the nuclei can be accurately modeled as an atom-centered point charge. The effect of the negative electron charge density is approximated also as a single atom-centered point charge. The value of this point charge is calculated by integrating the electron charge density inside a certain volume. Primarily Pysic uses the Bader partitioning scheme to determine these volumes, but a spherical partitioning scheme is also available.

Hydrogen links in Pysic

The covalent bonds between subsystems are treated with the link atom approach, i.e. hydrogen atoms are used to cap the bonds in the primary QM system. The hydrogen link atoms are placed on the line connecting the two atoms that form the bond. The exact position of the link atoms in this line can be controlled with the CHL parameter, given when defining the links in add_hydrogen_links(). The link atoms automatically keep their correct position when running dynamics or geometry optimization.

List of methods

The following methods form the basis for creating and performing hybrid simulations with HybridCalculator:

• add_subsystem()
• add_interaction()
• get_potential_energy()
• get_forces()
• print_energy_summary()
• print_force_summary()
• print_time_summary()
• print_interaction_charge_summary()
• view_subsystems()

Full documentation of the HybridCalculator class

class pysic.hybridcalculator.HybridCalculator(atoms=None)

Used to create and perform hybrid calculations.

This class is a fully compatible ASE calculator that provides the possibility of dividing the Atoms object into subsystems with different ASE calculators attached to them.

You can also define hydrogen link atoms in the interfaces of these subsystems or define any Pysic Potentials through which the subsystems interact with each other.

Attributes:

total_timer: :class:’~pysic.utility.timer.Timer’

Keeps track of the total time spent in the calculations.

subsystem_energies: dictionary

Name to energy.

subsystems: dictionary

Name to SubSystemInternal.

subsystem_interactions: dictionary

Pair of names to InteractionInternal.

subsystem_info: dictionary

Name to SubSystem.

interaction_info: dictionary

Pair of names to Interaction.

forces: numpy array

stress: None

potential_energy: float

system_initialized: bool

Indicates whether the SubSystemInternals and InteractionInternals have been constructed.

self.atoms: ASE Atoms

The full system.

add_interaction(interaction)

Used to add a interaction between two subsystems.

Parameters:

interaction: Interaction

The Interaction object containing the information.

add_subsystem(subsystem)

Used to define subsystems

You can define a subsystem with a oneliner, or then call setters on the SubSystemInfo object returned by this function.

Parameters:

subsystem: SubSystem object

calculate_forces()

Calculates the forces in the current structure.

This force includes the forces internal to the subsystems and the forces that bind the subsystems together.

calculate_potential_energy()

Calculates the potential energy in the current structure.

calculate_subsystem_interaction_charges(name)

Returns the calculated interaction charges for the given subsystem.

Before the charges can be calulated, one energy calculation has to be made so that the electron density is available

calculation_required(atoms, quantities=['forces', 'energy', 'stress'])

Check if a calculation is required for any of the the given properties.

Check if the quantities in the quantities list have already been calculated for the atomic configuration atoms. The quantities can be one or more of: ‘energy’, ‘forces’, ‘stress’. This method is used to check if a quantity is available without further calculations. For this reason, calculators should react to unknown/unsupported quantities by returning True, indicating that the quantity is not available.

Two sets of atoms are considered identical if they have the same positions, atomic numbers, unit cell and periodic boundary conditions.

Parameters:

atoms: ASE Atoms

This structure is compared to the currently stored.

quantities: list of strings

list of keywords ‘energy’, ‘forces’, ‘stress’

Returns:

bool: True if the quantities need to be calculated, false otherwise.

check_subsystem_indices(atom_indices, name)

Check that the atomic indices of the subsystem are present.

check_subsystem_overlap(atom_indices, name)

Check that the subsystem doesn’t overlap with another one.

full_system_set()

Checks whether the full system has been defined.

generate_subsystem_indices(name)

Generates the indices for the given subsystem.

get_atoms()

Return a copy of the full system.

get_colors()

Returns a color set for AtomEyeViewer with different colours for the different subsystems.

When the subsystems have been defined with add_subsystem() and the atoms have been set, this function will return a list of colors for each atom. The different subsystems have different colours for easy identification. You can provide this list of colors to an AtomEyeViewer object for visualization.

get_forces(atoms=None)

Returns the forces acting on the atoms.

If the atoms parameter is given, it will be used for updating the structure assigned to the calculator prior to calculating the forces. Otherwise the structure already associated with the calculator is used.

The calculator checks if the forces have been calculated already via calculation_required(). If the structure has changed, the forces are calculated using calculate_forces()

Parameters:

atoms: ASE Atoms object

The structure to calculate the forces on.

get_potential_energy(atoms=None, force_consistent=False)

Returns the potential energy of the hybrid system.

get_stress(atoms=None, skip_charge_relaxation=False)

This function has not been implemented.

get_subsystem(name)

Returns a copy of the ASE Atoms object for a certain subsystem.

The returned subsystem can be used for e.g. visualization or debugging.

get_subsystem_indices(name)

Return the indices of the atoms in the subsystem in the full system.

You can ask the indices even if the subsystems have not been initialized, but the indices of different subsystems may overlap in this case. If the subsystems have been initialized this function will only return indices if they are valid.

get_subsystem_pseudo_density(name)

Returns the electron pseudo density for the given subsystem.

get_unsubsystemized_atoms()

Return a list of indices for the atoms not already in a subsystem.

identical_atoms(atoms)

Compares the given atoms to the stored atoms object. The Atoms are identical if the positions, atomic numbers, periodic boundary conditions and cell are the same.

initialize_interaction(info)

Initializes a InteractionInternal from the given Interaction.

Parameters:

info: Interaction

initialize_subsystem(info)

Initializes a SubsystemInternal from the given Subsystem.

Parameters:

info: SubSystem

initialize_system()

Initializes the subsystems and interactions.

Called once during the lifetime of the calculator. Typically when calling set_atoms for the first time, or when calculating any quantity for the first time. If the atoms in the simulation need to be updated, update_system() is used.

print_energy_summary()

Print a detailed summary of the different energies in the system. This includes the energies in the subsystems, interaction energies and possible energy corrections.

print_force_summary()

Print a detailed summary of forces in the system.

print_interaction_charge_summary()

Print a summary of the atomic charges that are used in the electrostatic interaction between subsystems.

print_time_summary()

Print a detailed summary of the time usage.

set_atoms(atoms)

Set the full system for hybrid calculations.

Use add_subsystem() for setting up the subsystems.

Parameters:

atoms: ASE Atoms

The full system.

subsystem_defined(name)

Checks that there is a subsystem with the given name.

update_system(atoms)

Update the subsystem atoms.

view_subsystems()

Views the subsystems with ASE:s built in viewer.