Extending the framework

Extending the multisystem framework, in general, consists of several aspects:

For all of the above, there are certain routines which have to be implemented, and certain flags, which have to be defined. These will be outlined in this document.

Creating a new system

Create the Fortran class:

New systems must be created by extending the system_t class (or a class derived from system_t).

The class should own (by composition) the components and fields needed to describe the system, its (ground) state and its time propagation.

Every new system class must implement the deferred methods of interaction_partner_t and system_t class (or a class derived from system_t).

Write the constructor:

The ‘constructor’ must be a function which takes the namespace as argument and returns a pointer to the new object. As constructors are not overloading any deferred methods, they can take additional arguments, if needed.

The constructors must first allocate the new object, and assign the new objects namespace to the one given as argument. These are the minimum requirements to be operable with the system factory.

Here, we also need to define the supported interactions and declare exposed Quantities (see below).

In general, the constructor is also the right place to query related input variables and set the corresponding variables.

Here we can also define the supported interactions. These are stored in the arrays supported_interactions and supported_interactions_as_partner. These arrays need to be allocated first.

They can be allocated with size 0. See note below.

    allocate(this%supported_interactions(0))
    allocate(this%supported_interactions_as_partner(0))

We can then extend this list by adding new interactions, as shown below:

   this%supported_interactions = [this%supported_interactions, <NEW_INTERACTION>]
   this%supported_interactions_as_partner = [this%supported_interactions_as_partner, <NEW_INTERACTION>]

From the Fortran 2003 standard, this code automatically re-allocates the arrays with the new required size.

Define exposed quantities:

New systems will need to interact with other systems. This is done via interactions and so-called couplings.

A number of available Quantities, which can be used as couplings, are already defined in multisystem/core/quantity.F90:


  integer, public, parameter ::         &
    POSITION                     =  1,  &
    VELOCITY                     =  2,  &
    CURRENT                      =  3,  &
    DENSITY                      =  4,  &
    SCALAR_POTENTIAL             =  5,  &
    VECTOR_POTENTIAL             =  6,  &
    E_FIELD                      =  7,  &
    B_FIELD                      =  8,  &
    MASS                         =  9,  &
    CHARGE                       = 10,  &
    PERMITTIVITY                 = 11,  &
    PERMEABILITY                 = 12,  &
    E_CONDUCTIVITY               = 13,  &
    M_CONDUCTIVITY               = 14,  &
    DIPOLE                       = 15,  &
    MAX_QUANTITIES               = 15

If no definition matches your quantity, you can add it to this list.

The parent class interaction_partner_t provides an array of type quantity_t. For each quantity, the system exposes, we can define their properties (see quantity_t). For instance:

sys%quantities(DIPOLE)%updated_on_demand = .true.

indicates, that the dipole moment is used by an interaction, and needs to be updated in a routine system%update_quantity() (see below). Some quantities are needed for the propagation of the system itself (e.g. the position of a classical particle). For such quantities, we need to set update_on_demand = .false..

Add to system factory:

Octopus uses a so-called factory to create systems at run time. the system factory defines the allowed systems in


  integer, parameter, public ::             &
    SYSTEM_ELECTRONIC         = 1,  & !< electronic system (electrons_oct_m::electrons_t)
    SYSTEM_MAXWELL            = 2,  & !< maxwell system, (maxwell_oct_m::maxwell_t)
    SYSTEM_CLASSICAL_PARTICLE = 3,  & !< single classical particle (classical_particle_oct_m::classical_particle_t)
    SYSTEM_CHARGED_PARTICLE   = 4,  & !< single charged classical particle (charged_particle_oct_m::charged_particle_t)
    SYSTEM_DFTBPLUS           = 5,  & !< tight binding system (dftb_oct_m::dftb_t)
    SYSTEM_LINEAR_MEDIUM      = 6,  & !< linear medium for Maxwell calculations (linear_medium_oct_m::linear_medium_t)
    SYSTEM_MATTER             = 7,  & !< electrons including ions (matter_oct_m::matter_t)
    SYSTEM_DISPERSIVE_MEDIUM  = 8,  & !< dispersive medium for classical electrodynamics (dispersive_medium_oct_m::dispersive_medium_t)
    SYSTEM_MULTISYSTEM        = 9,  & !< container system. (multisystem_basic_oct_m::multisystem_basic_t)
    SYSTEM_IONS               = 10, & !< ensemble of charged classical particles (ions_oct_m::ions_t)
    SYSTEM_ENSEMBLE           = 11    !< ensemble container (ensemble_oct_m::ensemble_t)

Thus list needs to be extended to allow for new types.

the systems are created using the factory%create(namespace, type) call

Expand for the source of system_factory_create()

Implement algorithmic steps:

The propagation (and also the SCF algorithm) of the system is performed using abstract algorithms. These algorithms define algorithmic operations, which the systems might support.

For each system of a calculation, the multisystem framework calls the function system%do_algorithmic_operation(), which needs to be implemented by each system. This function consists of a large select case construct, which distinguishes the various algorithmic steps, the system can perform.

The function needs to return a logical indicating whether the step was successfully performed, and an array updated_quantities, which lists the non-on-demand quantities, which have been updated in the performed operation.

Creating a new interaction

Create the interaction class:

New interactions are derived from the abstract interaction_t class. For specific types of interactions, it might be beneficial to extend one of the derived classes, such as force_interaction_t , density_interaction_t or potential_interaction_t .

Interaction classes, in general’ keep copies of the partner’s couplings, which allows the partners to be updated without invalidating the interaction.

The constructor

Similar to systems, interactions need a constructor or init function, which is called by the interaction factory. This function must return a pointer to the new object, and take the interaction partner as argument.

Inside the constructor we need to set the label, and assign the partner pointer to the actual interaction partner. For the gravity_t interaction, this reads like

  function gravity_constructor(partner) result(this)
    class(interaction_partner_t), target, intent(inout) :: partner
    class(gravity_t),                     pointer       :: this

    PUSH_SUB(gravity_constructor)

    allocate(this)

    this%label = "gravity"

    this%partner => partner

    ! Gravity interaction needs two quantities from each system: the position and the mass
    ! From the sytem:
    this%system_quantities = [POSITION, MASS]

    ! From the partner:
    this%couplings_from_partner = [POSITION, MASS]

    POP_SUB(gravity_constructor)
  end function gravity_constructor

Like for systems, new interaction classes also must implement the deferred methods of interaction_t.

Specify couplings:

As seen in the example above, the interactions must declare which couplings they require from the system (on which the interaction has an effect), and the interaction partner (from which the interaction originates). This is done by adding the respective quantity ID’s to the arrays:

this%system_quantities
this%couplings_from_partner

Write initialization routines:

Certain properties of the interaction can only be initialized from the system or the partner.

For this, interaction systems and partners provide these deferred methods:

In general, this routine should call the init() routine of the supported interactions, if necessary perform other related tasks, and throw an fatal error if called for an unsupported interaction.

In general, interactions need to have data components with dimensions, which depend on the partner, e.g. the spatial dimensions. These quantities should be allocated and initialized in this routine.

Add to the interaction factory

Similar to the systems, interactions are created by an interaction factory. In order to register a new interaction, we must add it to the interactions enumerator in interactions/interactions_enum.F90


  integer, parameter, public :: &
    GRAVITY          = 1,       &
    LORENTZ_FORCE    = 2,       &
    COULOMB_FORCE    = 3,       &
    LINEAR_MEDIUM_TO_EM_FIELD = 4, &
    CURRENT_TO_MXLL_FIELD  = 5, &
    MXLL_E_FIELD_TO_MATTER   = 6,  &
    MXLL_B_FIELD_TO_MATTER   = 7,  &
    MXLL_VEC_POT_TO_MATTER   = 8,  &
    LENNARD_JONES            = 9

and to the function interactions_factory_create():

Source of interactions_factory_create

Connect interactions to the systems

Finally, the partners need to copy their couplings to the interaction and systems need to act on the interaction.

Here, the partners should copy their couplings to the corresponding variables of the interaction.

See, e.g. for classical_particle_t

Systems need to copy the data (e.g. forces) from the interaction. Usually, this is done when needed, as part of an algorithmic operation.

Some notes on on_demand quantities

For a more detailled discussion on how quantities are updated, see Quantities