process_intersection_derivatives Subroutine

private subroutine process_intersection_derivatives(inter_idx, k, sign_factor, intersection_results, intersection_sets, sys_geom, total_gradient, total_hessian, bonds, compute_grad, compute_hess, hess_dim, world_comm)

Uses

  • proc~~process_intersection_derivatives~~UsesGraph proc~process_intersection_derivatives process_intersection_derivatives module~mqc_config_parser mqc_config_parser proc~process_intersection_derivatives->module~mqc_config_parser module~mqc_error mqc_error proc~process_intersection_derivatives->module~mqc_error module~mqc_gmbe_utils mqc_gmbe_utils proc~process_intersection_derivatives->module~mqc_gmbe_utils module~mqc_physical_fragment mqc_physical_fragment proc~process_intersection_derivatives->module~mqc_physical_fragment module~mqc_result_types mqc_result_types proc~process_intersection_derivatives->module~mqc_result_types module~mqc_config_parser->module~mqc_error module~mqc_config_parser->module~mqc_physical_fragment module~mqc_calc_types mqc_calc_types module~mqc_config_parser->module~mqc_calc_types module~mqc_calculation_defaults mqc_calculation_defaults module~mqc_config_parser->module~mqc_calculation_defaults module~mqc_geometry mqc_geometry module~mqc_config_parser->module~mqc_geometry module~mqc_method_types mqc_method_types module~mqc_config_parser->module~mqc_method_types pic_types pic_types module~mqc_config_parser->pic_types module~mqc_gmbe_utils->module~mqc_error module~mqc_gmbe_utils->module~mqc_physical_fragment module~mqc_gmbe_utils->module~mqc_result_types module~mqc_combinatorics mqc_combinatorics module~mqc_gmbe_utils->module~mqc_combinatorics pic_io pic_io module~mqc_gmbe_utils->pic_io pic_logger pic_logger module~mqc_gmbe_utils->pic_logger module~mqc_gmbe_utils->pic_types module~mqc_physical_fragment->module~mqc_error module~mqc_cgto mqc_cgto module~mqc_physical_fragment->module~mqc_cgto module~mqc_elements mqc_elements module~mqc_physical_fragment->module~mqc_elements module~mqc_physical_fragment->module~mqc_geometry module~mqc_physical_constants mqc_physical_constants module~mqc_physical_fragment->module~mqc_physical_constants module~mqc_xyz_reader mqc_xyz_reader module~mqc_physical_fragment->module~mqc_xyz_reader module~mqc_physical_fragment->pic_types module~mqc_result_types->module~mqc_error pic_mpi_lib pic_mpi_lib module~mqc_result_types->pic_mpi_lib module~mqc_result_types->pic_types module~mqc_calc_types->pic_types module~mqc_calculation_defaults->pic_types module~mqc_cgto->pic_types module~mqc_combinatorics->pic_types module~mqc_elements->pic_types pic_ascii pic_ascii module~mqc_elements->pic_ascii module~mqc_geometry->pic_types module~mqc_method_types->pic_types module~mqc_physical_constants->pic_types module~mqc_xyz_reader->module~mqc_error module~mqc_xyz_reader->module~mqc_geometry module~mqc_xyz_reader->pic_types

Process derivatives for a single intersection fragment

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: inter_idx
integer, intent(in) :: k
real(kind=dp), intent(in) :: sign_factor
type(calculation_result_t), intent(in) :: intersection_results(:)
integer, intent(in) :: intersection_sets(:,:)
type(system_geometry_t), intent(in) :: sys_geom
real(kind=dp), intent(inout) :: total_gradient(:,:)
real(kind=dp), intent(inout), optional :: total_hessian(:,:)
type(bond_t), intent(in), optional :: bonds(:)
logical, intent(in) :: compute_grad
logical, intent(in) :: compute_hess
integer, intent(in) :: hess_dim
type(comm_t), intent(in), optional :: world_comm

Calls

proc~~process_intersection_derivatives~~CallsGraph proc~process_intersection_derivatives process_intersection_derivatives abort_comm abort_comm proc~process_intersection_derivatives->abort_comm error error proc~process_intersection_derivatives->error proc~build_fragment_from_atom_list build_fragment_from_atom_list proc~process_intersection_derivatives->proc~build_fragment_from_atom_list proc~error_get_full_trace error_t%error_get_full_trace proc~process_intersection_derivatives->proc~error_get_full_trace proc~error_has_error error_t%error_has_error proc~process_intersection_derivatives->proc~error_has_error proc~find_fragment_intersection find_fragment_intersection proc~process_intersection_derivatives->proc~find_fragment_intersection proc~fragment_destroy physical_fragment_t%fragment_destroy proc~process_intersection_derivatives->proc~fragment_destroy proc~get_monomer_atom_list get_monomer_atom_list proc~process_intersection_derivatives->proc~get_monomer_atom_list proc~redistribute_cap_gradients redistribute_cap_gradients proc~process_intersection_derivatives->proc~redistribute_cap_gradients proc~redistribute_cap_hessian redistribute_cap_hessian proc~process_intersection_derivatives->proc~redistribute_cap_hessian warning warning proc~process_intersection_derivatives->warning proc~build_fragment_from_atom_list->proc~error_has_error proc~add_hydrogen_caps add_hydrogen_caps proc~build_fragment_from_atom_list->proc~add_hydrogen_caps proc~check_duplicate_atoms check_duplicate_atoms proc~build_fragment_from_atom_list->proc~check_duplicate_atoms proc~count_hydrogen_caps count_hydrogen_caps proc~build_fragment_from_atom_list->proc~count_hydrogen_caps proc~error_add_context error_t%error_add_context proc~build_fragment_from_atom_list->proc~error_add_context proc~fragment_compute_nelec physical_fragment_t%fragment_compute_nelec proc~build_fragment_from_atom_list->proc~fragment_compute_nelec proc~error_get_full_trace->proc~error_has_error proc~basis_set_destroy molecular_basis_type%basis_set_destroy proc~fragment_destroy->proc~basis_set_destroy proc~atomic_basis_destroy atomic_basis_type%atomic_basis_destroy proc~basis_set_destroy->proc~atomic_basis_destroy proc~check_duplicate_atoms->error proc~element_number_to_symbol element_number_to_symbol proc~check_duplicate_atoms->proc~element_number_to_symbol proc~error_set error_t%error_set proc~check_duplicate_atoms->proc~error_set to_char to_char proc~check_duplicate_atoms->to_char proc~cgto_destroy cgto_type%cgto_destroy proc~atomic_basis_destroy->proc~cgto_destroy

Called by

proc~~process_intersection_derivatives~~CalledByGraph proc~process_intersection_derivatives process_intersection_derivatives proc~compute_gmbe compute_gmbe proc~compute_gmbe->proc~process_intersection_derivatives

Variables

Type Visibility Attributes Name Initial
integer, private, allocatable :: current_atoms(:)
integer, private :: current_n
logical, private :: has_intersection
type(error_t), private :: inter_error
type(physical_fragment_t), private :: inter_fragment
integer, private, allocatable :: intersect_atoms(:)
integer, private :: j
integer, private :: n_intersect
integer, private, allocatable :: next_atoms(:)
integer, private :: next_n
real(kind=dp), private, allocatable :: term_gradient(:,:)
real(kind=dp), private, allocatable :: term_hessian(:,:)

Source Code

   subroutine process_intersection_derivatives(inter_idx, k, sign_factor, intersection_results, &
                                               intersection_sets, sys_geom, total_gradient, &
                                               total_hessian, bonds, compute_grad, compute_hess, hess_dim, world_comm)
      !! Process derivatives for a single intersection fragment
      use mqc_result_types, only: calculation_result_t
      use mqc_physical_fragment, only: build_fragment_from_atom_list, &
                                       redistribute_cap_gradients, redistribute_cap_hessian
      use mqc_gmbe_utils, only: find_fragment_intersection
      use mqc_config_parser, only: bond_t
      use mqc_error, only: error_t

      integer, intent(in) :: inter_idx, k
      real(dp), intent(in) :: sign_factor
      type(calculation_result_t), intent(in) :: intersection_results(:)
      integer, intent(in) :: intersection_sets(:, :)
      type(system_geometry_t), intent(in) :: sys_geom
      real(dp), intent(inout) :: total_gradient(:, :)
      real(dp), intent(inout), optional :: total_hessian(:, :)
      type(bond_t), intent(in), optional :: bonds(:)
      logical, intent(in) :: compute_grad, compute_hess
      integer, intent(in) :: hess_dim
      type(comm_t), intent(in), optional :: world_comm

      integer :: j, current_n, next_n, n_intersect
      integer, allocatable :: current_atoms(:), next_atoms(:), intersect_atoms(:)
      real(dp), allocatable :: term_gradient(:, :), term_hessian(:, :)
      logical :: has_intersection
      type(physical_fragment_t) :: inter_fragment
      type(error_t) :: inter_error

      ! Build the intersection atom list
      call get_monomer_atom_list(sys_geom, intersection_sets(1, inter_idx), current_atoms, current_n)

      if (current_n > 0) then
         do j = 2, k
            call get_monomer_atom_list(sys_geom, intersection_sets(j, inter_idx), next_atoms, next_n)
            has_intersection = find_fragment_intersection(current_atoms, current_n, &
                                                          next_atoms, next_n, &
                                                          intersect_atoms, n_intersect)
            deallocate (current_atoms, next_atoms)

            if (.not. has_intersection) then
               current_n = 0
               exit
            end if

            current_n = n_intersect
            allocate (current_atoms(current_n))
            current_atoms = intersect_atoms
            deallocate (intersect_atoms)
         end do
      end if

      if (current_n > 0) then
         call build_fragment_from_atom_list(sys_geom, current_atoms, current_n, &
                                            inter_fragment, inter_error, bonds)
         if (inter_error%has_error()) then
            call logger%error(inter_error%get_full_trace())
            if (present(world_comm)) then
               call abort_comm(world_comm, 1)
            else
               error stop "Failed to build intersection fragment in GMBE"
            end if
         end if

         if (compute_grad .and. intersection_results(inter_idx)%has_gradient) then
            allocate (term_gradient(3, sys_geom%total_atoms))
            term_gradient = 0.0_dp
            call redistribute_cap_gradients(inter_fragment, intersection_results(inter_idx)%gradient, &
                                            term_gradient)
            total_gradient = total_gradient + sign_factor*term_gradient
            deallocate (term_gradient)
         end if

         if (compute_hess .and. intersection_results(inter_idx)%has_hessian) then
            allocate (term_hessian(hess_dim, hess_dim))
            term_hessian = 0.0_dp
            call redistribute_cap_hessian(inter_fragment, intersection_results(inter_idx)%hessian, &
                                          term_hessian)
            total_hessian = total_hessian + sign_factor*term_hessian
            deallocate (term_hessian)
         end if

         call inter_fragment%destroy()
      else
         call logger%warning("GMBE intersection has no atoms; skipping derivatives")
      end if

      if (allocated(current_atoms)) deallocate (current_atoms)
   end subroutine process_intersection_derivatives