compute_mbe_dipole Subroutine

private subroutine compute_mbe_dipole(fragment_idx, fragment, lookup, results, delta_dipoles, n, world_comm)

Uses

  • proc~~compute_mbe_dipole~~UsesGraph proc~compute_mbe_dipole compute_mbe_dipole module~mqc_result_types mqc_result_types proc~compute_mbe_dipole->module~mqc_result_types module~mqc_error mqc_error module~mqc_result_types->module~mqc_error pic_mpi_lib pic_mpi_lib module~mqc_result_types->pic_mpi_lib pic_types pic_types module~mqc_result_types->pic_types

Bottom-up computation of n-body dipole correction Exactly mirrors the energy MBE logic: deltaDipole = Dipole - sum(all subset deltaDipoles) Dipoles are additive vectors in the system frame, no coordinate mapping needed

Arguments

Type IntentOptional Attributes Name
integer(kind=int64), intent(in) :: fragment_idx
integer, intent(in) :: fragment(:)
type(fragment_lookup_t), intent(in) :: lookup
type(calculation_result_t), intent(in) :: results(:)
real(kind=dp), intent(inout) :: delta_dipoles(:,:)

(3, fragment_count)
integer, intent(in) :: n
type(comm_t), intent(in), optional :: world_comm

MPI communicator for abort

Calls

proc~~compute_mbe_dipole~~CallsGraph proc~compute_mbe_dipole compute_mbe_dipole abort_comm abort_comm proc~compute_mbe_dipole->abort_comm error error proc~compute_mbe_dipole->error proc~fragment_lookup_find fragment_lookup_t%fragment_lookup_find proc~compute_mbe_dipole->proc~fragment_lookup_find proc~get_next_combination get_next_combination proc~compute_mbe_dipole->proc~get_next_combination fnv_1a_hash fnv_1a_hash proc~fragment_lookup_find->fnv_1a_hash proc~arrays_equal_internal arrays_equal_internal proc~fragment_lookup_find->proc~arrays_equal_internal sort sort proc~fragment_lookup_find->sort

Called by

proc~~compute_mbe_dipole~~CalledByGraph proc~compute_mbe_dipole compute_mbe_dipole proc~compute_mbe compute_mbe proc~compute_mbe->proc~compute_mbe_dipole proc~global_coordinator global_coordinator proc~global_coordinator->proc~compute_mbe proc~serial_fragment_processor serial_fragment_processor proc~serial_fragment_processor->proc~compute_mbe interface~global_coordinator global_coordinator interface~global_coordinator->proc~global_coordinator interface~serial_fragment_processor serial_fragment_processor interface~serial_fragment_processor->proc~serial_fragment_processor proc~mbe_run_distributed mbe_context_t%mbe_run_distributed proc~mbe_run_distributed->interface~global_coordinator proc~mbe_run_serial mbe_context_t%mbe_run_serial proc~mbe_run_serial->interface~serial_fragment_processor

Variables

Type Visibility Attributes Name Initial
logical, private :: has_next
integer, private :: i
integer, private :: indices(MAX_MBE_LEVEL)
integer, private :: subset(MAX_MBE_LEVEL)
integer(kind=int64), private :: subset_idx
integer, private :: subset_size

Source Code

   subroutine compute_mbe_dipole(fragment_idx, fragment, lookup, results, delta_dipoles, n, world_comm)
      !! Bottom-up computation of n-body dipole correction
      !! Exactly mirrors the energy MBE logic: deltaDipole = Dipole - sum(all subset deltaDipoles)
      !! Dipoles are additive vectors in the system frame, no coordinate mapping needed
      use mqc_result_types, only: calculation_result_t
      integer(int64), intent(in) :: fragment_idx
      integer, intent(in) :: fragment(:), n
      type(fragment_lookup_t), intent(in) :: lookup
      type(calculation_result_t), intent(in) :: results(:)
      real(dp), intent(inout) :: delta_dipoles(:, :)  !! (3, fragment_count)
      type(comm_t), intent(in), optional :: world_comm  !! MPI communicator for abort

      integer :: subset_size, i
      integer :: indices(MAX_MBE_LEVEL), subset(MAX_MBE_LEVEL)  ! Stack arrays to avoid heap contention
      integer(int64) :: subset_idx
      logical :: has_next

      ! Start with the full n-mer dipole
      delta_dipoles(:, fragment_idx) = results(fragment_idx)%dipole

      ! Subtract all proper subsets (size 1 to n-1)
      do subset_size = 1, n - 1
         ! Initialize first combination
         do i = 1, subset_size
            indices(i) = i
         end do

         ! Loop through all combinations
         do
            ! Build current subset
            do i = 1, subset_size
               subset(i) = fragment(indices(i))
            end do

            ! Look up subset index
            subset_idx = lookup%find(subset(1:subset_size), subset_size)
            if (subset_idx < 0) then
               call logger%error("Subset not found in MBE dipole computation")
               if (present(world_comm)) then
                  call abort_comm(world_comm, 1)
               else
                  error stop "Subset not found in MBE dipole!"
               end if
            end if

            ! Subtract pre-computed delta dipole
            delta_dipoles(:, fragment_idx) = delta_dipoles(:, fragment_idx) - &
                                             delta_dipoles(:, subset_idx)

            ! Get next combination
            call get_next_combination(indices, subset_size, n, has_next)
            if (.not. has_next) exit
         end do
      end do

   end subroutine compute_mbe_dipole