print_detailed_breakdown Subroutine

public subroutine print_detailed_breakdown(polymers, fragment_count, max_level, energies, delta_energies)

Print detailed energy breakdown for each fragment Shows full energy and deltaE correction for all monomers, dimers, trimers, etc. Uses int64 for fragment_count to handle large fragment counts that overflow int32.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: polymers(:,:)
integer(kind=int64), intent(in) :: fragment_count
integer, intent(in) :: max_level
real(kind=dp), intent(in) :: energies(:)
real(kind=dp), intent(in) :: delta_energies(:)

Calls

proc~~print_detailed_breakdown~~CallsGraph proc~print_detailed_breakdown print_detailed_breakdown header header proc~print_detailed_breakdown->header level_name level_name proc~print_detailed_breakdown->level_name proc~get_frag_level_name get_frag_level_name proc~print_detailed_breakdown->proc~get_frag_level_name verbose verbose proc~print_detailed_breakdown->verbose warning warning proc~print_detailed_breakdown->warning

Called by

proc~~print_detailed_breakdown~~CalledByGraph proc~print_detailed_breakdown print_detailed_breakdown proc~compute_mbe compute_mbe proc~compute_mbe->proc~print_detailed_breakdown 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
integer(kind=int64), private :: count_by_level
character(len=512), private :: energy_line
integer, private :: frag_level
integer, private :: fragment_size
character(len=512), private :: fragment_str
integer(kind=int64), private :: i
integer, private :: j

Source Code

   subroutine print_detailed_breakdown(polymers, fragment_count, max_level, energies, delta_energies)
      !! Print detailed energy breakdown for each fragment
      !! Shows full energy and deltaE correction for all monomers, dimers, trimers, etc.
      !! Uses int64 for fragment_count to handle large fragment counts that overflow int32.
      integer, intent(in) :: polymers(:, :), max_level
      integer(int64), intent(in) :: fragment_count
      real(dp), intent(in) :: energies(:), delta_energies(:)

      integer(int64) :: i
      integer :: fragment_size, j, frag_level
      character(len=512) :: fragment_str, energy_line
      integer(int64) :: count_by_level

      call logger%verbose(" ")
      call logger%verbose("============================================")
      call logger%verbose("Detailed Energy Breakdown by Fragment")
      call logger%verbose("============================================")

      ! Warn if we have very high fragmentation levels
      if (max_level > 10) then
         call logger%warning("Fragment levels exceed decamers (10-mers). Using generic N-mers notation.")
      end if

      do frag_level = 1, max_level
         count_by_level = 0_int64

         do i = 1_int64, fragment_count
            fragment_size = count(polymers(i, :) > 0)
            if (fragment_size == frag_level) count_by_level = count_by_level + 1_int64
         end do

         if (count_by_level > 0_int64) then
            call logger%verbose(" ")
            block
               character(len=256) :: header
               character(len=32) :: level_name
               level_name = get_frag_level_name(frag_level)
               write (header, '(a,a,i0,a)') trim(level_name), " (", count_by_level, " fragments):"
               ! Capitalize first letter
               if (len_trim(level_name) > 0) then
                  if (level_name(1:1) >= 'a' .and. level_name(1:1) <= 'z') then
                     header(1:1) = achar(iachar(header(1:1)) - 32)
                  end if
               end if
               call logger%verbose(trim(header))
            end block
            call logger%verbose("--------------------------------------------")

            do i = 1_int64, fragment_count
               fragment_size = count(polymers(i, :) > 0)

               if (fragment_size == frag_level) then
                  fragment_str = "["
                  do j = 1, fragment_size
                     if (j > 1) then
                        write (fragment_str, '(a,a,i0)') trim(fragment_str), ",", polymers(i, j)
                     else
                        write (fragment_str, '(a,i0)') trim(fragment_str), polymers(i, j)
                     end if
                  end do
                  write (fragment_str, '(a,a)') trim(fragment_str), "]"

                  if (frag_level == 1) then
                     write (energy_line, '(a,a,f20.10)') &
                        "  Fragment ", trim(adjustl(fragment_str)), energies(i)
                  else
                     write (energy_line, '(a,a,f20.10,a,f20.10)') &
                        "  Fragment ", trim(adjustl(fragment_str)), energies(i), &
                        "   deltaE: ", delta_energies(i)
                  end if
                  call logger%verbose(trim(energy_line))
               end if
            end do
         end if
      end do

      call logger%verbose(" ")
      call logger%verbose("============================================")

   end subroutine print_detailed_breakdown