Time propagation

In the new multisystem framework, the time propagation is handled by the routine time_dependent_run_multisystem()

The main loop

  subroutine time_dependent_run_multisystem(systems, from_scratch)
    class(multisystem_basic_t), intent(inout) :: systems
    logical,                    intent(in)    :: from_scratch

    logical :: trigger_restart, stop_code, stop_loop
    logical :: restart_read

    PUSH_SUB(time_dependent_run_multisystem)

    call multisystem_debug_init("debug/multisystem_propagation.log", global_namespace, systems%grp)

    call messages_write('Info: Running Multi-system time evolution')
    call messages_info(namespace=systems%namespace)

    ! Read restart files or set initial conditions
    call systems%init_iteration_counters()
    if (.not. from_scratch) then
      restart_read = systems%restart_read()
    else
      restart_read = .false.
    end if
    if (restart_read) then
      message(1) = "Successfully read restart data for all system."
      call messages_info(1, namespace=systems%namespace)
    else
      call systems%initial_conditions()
    end if

    call systems%propagation_start()

    call multisystem_debug_write_marker(systems%namespace, event_marker_t("propagation_start"))

    ! The full TD loop
    stop_loop = .false.
    do while (.not. systems%algorithm_finished())
      ! Execute algorithm until next barrier
      call systems%execute_algorithm()

      ! determine cases in which to trigger writing restart files
      stop_code = clean_stop(systems%grp%comm) .or. walltimer_alarm(systems%grp%comm)
      trigger_restart = stop_code .or. restart_walltime_period_alarm(systems%grp%comm)

      if (trigger_restart .and. .not. stop_loop) then
        call systems%start_barrier(systems%next_time_on_largest_dt(), BARRIER_RESTART)
        stop_loop = stop_code
        trigger_restart = .false.
      end if
      if (systems%arrived_at_barrier(BARRIER_RESTART)) then
        call systems%restart_write()
        call systems%end_barrier(BARRIER_RESTART)
        if (stop_loop) exit
      end if
    end do

    if (systems%algorithm_finished()) then
      call systems%restart_write()
    end if

    call multisystem_debug_write_marker(systems%namespace, event_marker_t("propagation_finish"))

    call systems%propagation_finish()

    call multisystem_debug_end()

    POP_SUB(time_dependent_run_multisystem)
  end subroutine time_dependent_run_multisystem

The code is quite self-explanatory. The central part of the time propagation is the do ... while loop (see "! The full TD loop"). Here the function system%dt_operation() is called successively, triggering one algorithmic step of the multisystem propagator.

Updating the system

In the multisystem framework, the highest level system is of type multisystem_basic_t, which inherits the dt_operation() routine from the abstract multisystem_t type.

This routine first calls the general system_dt_operation() subroutine from the parent class system_t and then loops over all subsystems, which are part of the multisystem, and calls their specific system%dt_operations() routine.

The general system_dt_operation() subroutine is handling all general algorithmic operations, such as starting and ending the SCF loop, or updating interaction.

This routine fetches the next operation from the propagator and executes it. The following operations are implemented here:

• UPDATE_INTRERACTIONS
  • increment the propagator-clock by one tick (if all interactions can be updated)
  • try to update all interactions (this might not always be possible, if some interaction partners are still lagging behind)
  • if all interactions are updated succesfully, progress the propagator to the next step.
• FINISHED
  • increment the system-clock by one tick.
  • rewind the propagator
• START_SCF_LOOP and END_SCF_LOOP
  • handles self-consistency and convergence for predictor-corrector algorithms.
  • no clocks are updated.

All remaining operations are passed down to the system specific routine system%do_dt_operation(), which is the place where specific systems need to implement all algorithmic steps, required by the propagators, allowed for that system.

This routine is triggering the update of the interactions, via the call to system%update_interations(), which loops over the interactions associated with the system, and all required exposed quantities.

Updating the interactions

Each of the systems (i.e. the multisystem and its subsystems) are attempting to update their interactions with system%update_interactions().

The first part makes sure that the update_interactions_start() routine is only called when no interaction has been updated yet in this time step, iu.e. if their clocks are behind the system clock.

In the second part, we actually attempt to update the interactions, if needed. If an interaction is not up-to-date, it first needs to be ensured that the quantities, on which the interaction depends, are up-to-date. Once all quantities are updated, the update() routine of the interaction can be called with the current propagator time as argument.

Finally, if all interactions are updated successfully, there is a step update_interactions_finish(), which can be necessary for some systems.

This function is to be called with for a specific requested_time, which is the time at which the interaction is needed.

The function first tries to update the exposed quantities of the partner. If this is successfull, the function calls the calculate() routine of the interaction, sets the interaction clock to the requested time and returns with the result .true., indicating a successful update. In case the function was blocked by one of the exposed quantities, the interaction is not calculated and the function returns .false..

As can be seen, it is, possible that an interaction cannot be updated at a given time. This can be the case when the interaction also depends on quantities of the other interaction partner, and that partner is not yet at the requested time.

This routine demands some more comments:

Firstly, the routine is only allowed if the interaction is of type interaction_with_partner_t, as otherwise there is no partner to update exposed quantities. It is important to remember who called the function, and which quantities are supposed to be updated:

system_update_exposed_quantities() is called from Implementation of interaction_with_partner_update() for each interaction of a system for the interaction partner of that system with respect to a given interaction.

In the following situations, the routine is not allowd to update the exposed quantities:

• we are in an internal SCF loop of the propagation step
• one of the quantities cannot be updated becuase it is either
  • too far behind in time (where two far is more than one clock tick)
  • it is one clock tick behind, but the user requested OPTION__INTERACTIONTIMING__TIMING_EXACT In these cases, the routine will return a .false.. The propagator will continue, but the system in question will be held at this step. Other systems, however, can continue to progress, and might be able to get to the time, where then this quantity can be updated.

If the timing conditions are fulfilled (and we are not in an internal SCF loop) the system will call the specific routines to update the exposed quantities, and copies their values to the interaction.

Deadlocks and how to avoid them