Algorithms#

Algorithms in QC Lab define the sequence of operations that evolve the system defined by the model object (see Models) in time. They are composed of three recipes which define initialization steps, update steps, and collect steps that together define the desired algorithm. Each recipe is a list of “tasks” (see Tasks) which are executed in the order specified by the recipe list. Algorithms define the transient quantities of an algorithm in the state object, which is an instance of a dictionary.

Algorithm Objects#

Algorithm objects in QC Lab are instances of the qclab.Algorithm class. Each algorithm object is composed of three recipes: an initialization recipe algorithm.initialization_recipe, an update recipe algorithm.update_recipe, and a collect recipe algorithm.collect_recipe. Like a model object, an algorithm object has an instance of the Constants class algorithm.settings which contains the settings specific to the algorithm. Unlike the model object, algorithm objects do not have internal constants and so there is no initialization method as there is for model objects (see Models). Instead, the settings of the algorithm object are set directly by the user during or after instantiation of the algorithm object.

The empty Algorithm class is:

class Algorithm:
    """
    Algorithm class for defining and executing algorithm recipes.
    """

    def __init__(self, default_settings=None, settings=None):
        if settings is None:
            settings = {}
        if default_settings is None:
            default_settings = {}
        # Merge default settings with user-provided settings.
        settings = {**default_settings, **settings}
        # Construct a Constants object to hold settings.
        self.settings = Constants()
        # Put settings from the dictionary into the Constants object.
        for key, val in settings.items():
            setattr(self.settings, key, val)
        # Copy the recipes and output variables to ensure they are not shared
        # across instances.
        self.initialization_recipe = copy.deepcopy(self.initialization_recipe)
        self.update_recipe = copy.deepcopy(self.update_recipe)
        self.collect_recipe = copy.deepcopy(self.collect_recipe)

    initialization_recipe = []
    update_recipe = []
    collect_recipe = []

    def execute_recipe(self, sim, state, parameters, recipe):
        """
        Carry out the given recipe for the simulation by running
        each task in the recipe.
        """
        for func in recipe:
            state, parameters = func(sim, state, parameters)
        return state, parameters

After instantiating an algorithm object, the user can populate its recipes by assigning tasks to each recipe. For example, the mean-field algorithm can be defined from an empty Algorithm object as:

from qclab import Algorithm
import qclab.tasks as tasks
from functools import partial

# Create an empty algorithm object.
algorithm = Algorithm()
# Populate the initialization recipe.
algorithm.initialization_recipe = [
        tasks.initialize_variable_objects,
        tasks.initialize_norm_factor,
        tasks.initialize_z,
        tasks.update_h_q_tot,
]
# Populate the update recipe.
algorithm.update_recipe = [
    # Begin RK4 integration steps.
    # RK4 steps excluded for brevity.
    # End RK4 integration steps.
    tasks.update_wf_db_rk4,
    tasks.update_h_q
]
# Populate the collect recipe.
algorithm.collect_recipe = [
    tasks.update_t,
    tasks.update_dm_db_mf,
    tasks.update_quantum_energy,
    tasks.update_classical_energy,
    tasks.collect_t,
    tasks.collect_dm_db,
    tasks.collect_classical_energy,
    tasks.collect_quantum_energy,
]

Each recipe is executed by the method algorithm.execute_recipe. The initialization recipe is executed once at the beginning of the simulation, the update recipe is executed at each time step of the simulation, and the collect recipe is executed once at the end of the simulation to gather and process results.

Mean Field Example#

As an example of a complete algorithm we include the source code for the mean-field algorithm below. This algorithm is defined in the qclab.algorithms.MeanField module and uses tasks from the qclab.tasks module to populate its recipes.

FSSH Collected Observables#

Key

Description

quantum_energy

The quantum energy of the system.

classical_energy

The classical energy of the system.

dm_db

The diabatic density matrix of the quantum subsystem.

t

The time points of the simulation.
View full source
 1"""
 2This module contains the MeanField algorithm class.
 3"""
 4
 5from functools import partial
 6from qclab.algorithm import Algorithm
 7from qclab import tasks
 8
 9
10class MeanField(Algorithm):
11    """
12    Mean-field dynamics algorithm class.
13    """
14
15    def __init__(self, settings=None):
16        if settings is None:
17            settings = {}
18        self.default_settings = {}
19        super().__init__(self.default_settings, settings)
20
21    initialization_recipe = [
22        tasks.initialize_variable_objects,
23        tasks.initialize_norm_factor,
24        tasks.initialize_z,
25        tasks.update_h_q_tot,
26    ]
27
28    update_recipe = [
29        # Begin RK4 integration steps.
30        partial(tasks.update_classical_force, z_name="z"),
31        tasks.update_quantum_classical_force,
32        tasks.update_z_rk4_k123,
33        partial(tasks.update_classical_force, z_name="z_1"),
34        partial(
35            tasks.update_quantum_classical_force,
36            z_name="z_1",
37            wf_changed=False,
38        ),
39        partial(tasks.update_z_rk4_k123, z_name="z", z_k_name="z_2", k_name="z_rk4_k2"),
40        partial(tasks.update_classical_force, z_name="z_2"),
41        partial(
42            tasks.update_quantum_classical_force,
43            z_name="z_2",
44            wf_changed=False,
45        ),
46        partial(
47            tasks.update_z_rk4_k123,
48            z_name="z",
49            z_k_name="z_3",
50            k_name="z_rk4_k3",
51            dt_factor=1.0,
52        ),
53        partial(tasks.update_classical_force, z_name="z_3"),
54        partial(
55            tasks.update_quantum_classical_force,
56            z_name="z_3",
57            wf_changed=False,
58        ),
59        tasks.update_z_rk4_k4,
60        # End RK4 integration steps.
61        tasks.update_wf_db_rk4,
62        tasks.update_h_q_tot,
63    ]
64
65    collect_recipe = [
66        tasks.update_t,
67        tasks.update_dm_db_wf,
68        tasks.update_quantum_energy_wf,
69        tasks.update_classical_energy,
70        tasks.collect_t,
71        tasks.collect_dm_db,
72        tasks.collect_classical_energy,
73        tasks.collect_quantum_energy,
74    ]

Surface Hopping Example#

As an additional example of a complete algorithm we include the source code for the fewest-switches surface hopping algorithm below. This algorithm is defined in the qclab.algorithms.FewestSwitchesSurfaceHopping module and uses tasks from the qclab.tasks module to populate its recipes.

FSSH Collected Observables#

Key

Description

quantum_energy

The quantum energy of the system.

classical_energy

The classical energy of the system.

dm_db

The diabatic density matrix of the quantum subsystem.

t

The time points of the simulation.
View full source
  1"""
  2This module contains the FewestSwitchesSurfaceHopping algorithm class.
  3"""
  4
  5from functools import partial
  6from qclab.algorithm import Algorithm
  7from qclab import tasks
  8
  9
 10class FewestSwitchesSurfaceHopping(Algorithm):
 11    """
 12    Fewest switches surface hopping algorithm class.
 13    """
 14
 15    def __init__(self, settings=None):
 16        if settings is None:
 17            settings = {}
 18        self.default_settings = {
 19            "fssh_deterministic": False,
 20            "gauge_fixing": "sign_overlap",
 21            "use_gauge_field_force": False,
 22        }
 23        super().__init__(self.default_settings, settings)
 24
 25    initialization_recipe = [
 26        tasks.initialize_variable_objects,
 27        tasks.initialize_norm_factor,
 28        tasks.initialize_branch_seeds,
 29        tasks.initialize_z,
 30        tasks.update_h_q_tot,
 31        partial(
 32            tasks.diagonalize_matrix,
 33            matrix_name="h_q_tot",
 34            eigvals_name="eigvals",
 35            eigvecs_name="eigvecs",
 36        ),
 37        partial(
 38            tasks.update_eigvecs_gauge,
 39            gauge_fixing="phase_der_couple",
 40            eigvecs_previous_name="eigvecs",
 41        ),
 42        partial(tasks.copy_in_state, copy_name="eigvecs_previous", orig_name="eigvecs"),
 43        partial(
 44            tasks.update_vector_basis,
 45            input_vec_name="wf_db",
 46            basis_name="eigvecs",
 47            output_vec_name="wf_adb",
 48            adb_to_db=False,
 49        ),
 50        tasks.initialize_random_values_fssh,
 51        tasks.initialize_active_surface,
 52        tasks.initialize_dm_adb_0_fssh,
 53        tasks.update_act_surf_wf,
 54    ]
 55
 56    update_recipe = [
 57        partial(tasks.copy_in_state, copy_name="eigvecs_previous", orig_name="eigvecs"),
 58        # Begin RK4 integration steps.
 59        tasks.update_classical_force,
 60        partial(
 61            tasks.update_quantum_classical_force,
 62            wf_db_name="act_surf_wf",
 63            wf_changed=True,
 64        ),
 65        tasks.update_z_rk4_k123,
 66        partial(tasks.update_classical_force, z_name="z_1"),
 67        partial(
 68            tasks.update_quantum_classical_force,
 69            wf_db_name="act_surf_wf",
 70            z_name="z_1",
 71            wf_changed=False,
 72        ),
 73        partial(
 74            tasks.update_z_rk4_k123, z_name="z", z_k_name="z_2", k_name="z_rk4_k2"
 75        ),
 76        partial(tasks.update_classical_force, z_name="z_2"),
 77        partial(
 78            tasks.update_quantum_classical_force,
 79            wf_db_name="act_surf_wf",
 80            z_name="z_2",
 81            wf_changed=False,
 82        ),
 83        partial(
 84            tasks.update_z_rk4_k123,
 85            z_name="z",
 86            z_k_name="z_3",
 87            k_name="z_rk4_k3",
 88            dt_factor=1.0,
 89        ),
 90        partial(tasks.update_classical_force, z_name="z_3"),
 91        partial(
 92            tasks.update_quantum_classical_force,
 93            wf_db_name="act_surf_wf",
 94            z_name="z_3",
 95            wf_changed=False,
 96        ),
 97        tasks.update_z_rk4_k4,
 98        # End RK4 integration steps.
 99        tasks.update_wf_db_propagator,
100        tasks.update_h_q_tot,
101        partial(
102            tasks.diagonalize_matrix,
103            matrix_name="h_q_tot",
104            eigvals_name="eigvals",
105            eigvecs_name="eigvecs",
106        ),
107        tasks.update_eigvecs_gauge,
108        partial(
109            tasks.update_vector_basis,
110            input_vec_name="wf_db",
111            basis_name="eigvecs",
112            output_vec_name="wf_adb",
113            adb_to_db=False,
114        ),
115        tasks.update_hop_prob_fssh,
116        tasks.update_hop_inds_fssh,
117        tasks.update_hop_vals_fssh,
118        tasks.update_z_hop,
119        tasks.update_act_surf_hop,
120        tasks.update_act_surf_wf,
121    ]
122
123    collect_recipe = [
124        tasks.update_t,
125        tasks.update_dm_db_fssh,
126        tasks.update_quantum_energy_act_surf,
127        tasks.update_classical_energy_fssh,
128        tasks.collect_t,
129        tasks.collect_dm_db,
130        tasks.collect_quantum_energy,
131        tasks.collect_classical_energy,
132    ]