build_mbe_lookup_table Subroutine

private subroutine build_mbe_lookup_table(polymers, fragment_count, max_level, lookup, error)

Uses

  • proc~~build_mbe_lookup_table~~UsesGraph proc~build_mbe_lookup_table build_mbe_lookup_table module~mqc_error mqc_error proc~build_mbe_lookup_table->module~mqc_error

Build hash table for fast fragment lookups

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: polymers(:,:)
integer(kind=int64), intent(in) :: fragment_count
integer, intent(in) :: max_level
type(fragment_lookup_t), intent(inout) :: lookup
type(error_t), intent(out), optional :: error

Calls

proc~~build_mbe_lookup_table~~CallsGraph proc~build_mbe_lookup_table build_mbe_lookup_table debug debug proc~build_mbe_lookup_table->debug get_elapsed_time get_elapsed_time proc~build_mbe_lookup_table->get_elapsed_time proc~error_add_context error_t%error_add_context proc~build_mbe_lookup_table->proc~error_add_context proc~error_has_error error_t%error_has_error proc~build_mbe_lookup_table->proc~error_has_error proc~fragment_lookup_init fragment_lookup_t%fragment_lookup_init proc~build_mbe_lookup_table->proc~fragment_lookup_init proc~fragment_lookup_insert fragment_lookup_t%fragment_lookup_insert proc~build_mbe_lookup_table->proc~fragment_lookup_insert start start proc~build_mbe_lookup_table->start to_char to_char proc~build_mbe_lookup_table->to_char proc~next_prime_internal next_prime_internal proc~fragment_lookup_init->proc~next_prime_internal fnv_1a_hash fnv_1a_hash proc~fragment_lookup_insert->fnv_1a_hash proc~error_set error_t%error_set proc~fragment_lookup_insert->proc~error_set sort sort proc~fragment_lookup_insert->sort

Called by

proc~~build_mbe_lookup_table~~CalledByGraph proc~build_mbe_lookup_table build_mbe_lookup_table proc~compute_mbe compute_mbe proc~compute_mbe->proc~build_mbe_lookup_table 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, private :: fragment_size
integer(kind=int64), private :: i
type(error_t), private :: insert_error
type(timer_type), private :: lookup_timer

Source Code

   subroutine build_mbe_lookup_table(polymers, fragment_count, max_level, lookup, error)
      !! Build hash table for fast fragment lookups
      use mqc_error, only: error_t
      integer, intent(in) :: polymers(:, :)
      integer(int64), intent(in) :: fragment_count
      integer, intent(in) :: max_level
      type(fragment_lookup_t), intent(inout) :: lookup
      type(error_t), intent(out), optional :: error

      integer(int64) :: i
      integer :: fragment_size
      type(timer_type) :: lookup_timer
      type(error_t) :: insert_error

      call lookup_timer%start()
      call lookup%init(fragment_count)
      do i = 1_int64, fragment_count
         fragment_size = count(polymers(i, :) > 0)
         call lookup%insert(polymers(i, :), fragment_size, i, insert_error)
         if (insert_error%has_error()) then
            if (present(error)) then
               error = insert_error
               call error%add_context("build_mbe_lookup_table")
            end if
            return
         end if
      end do
      call lookup_timer%stop()
      call logger%debug("Time to build lookup table: "//to_char(lookup_timer%get_elapsed_time())//" s")
      call logger%debug("Hash table size: "//to_char(lookup%table_size)// &
                        ", entries: "//to_char(lookup%n_entries))
   end subroutine build_mbe_lookup_table