flow_global_two_dots_owned Subroutine

public subroutine flow_global_two_dots_owned(flow, a1, b1, a2, b2, results)

Uses

  • proc~~flow_global_two_dots_owned~~UsesGraph proc~flow_global_two_dots_owned mod_mpi_flow::flow_global_two_dots_owned module~mod_profiler mod_profiler proc~flow_global_two_dots_owned->module~mod_profiler iso_fortran_env iso_fortran_env module~mod_profiler->iso_fortran_env module~mod_kinds mod_kinds module~mod_profiler->module~mod_kinds mpi_f08 mpi_f08 module~mod_profiler->mpi_f08 module~mod_kinds->iso_fortran_env

Computes two global dot products without constructing temporary full-size batches.

Arguments

Type IntentOptional Attributes Name
type(flow_mpi_t), intent(in) :: flow
real(kind=rk), intent(in) :: a1(:)
real(kind=rk), intent(in) :: b1(:)
real(kind=rk), intent(in) :: a2(:)
real(kind=rk), intent(in) :: b2(:)
real(kind=rk), intent(out) :: results(2)

Calls

proc~~flow_global_two_dots_owned~~CallsGraph proc~flow_global_two_dots_owned mod_mpi_flow::flow_global_two_dots_owned mpi_allreduce mpi_allreduce proc~flow_global_two_dots_owned->mpi_allreduce proc~check_mpi~2 mod_mpi_flow::check_mpi proc~flow_global_two_dots_owned->proc~check_mpi~2 proc~profiler_start mod_profiler::profiler_start proc~flow_global_two_dots_owned->proc~profiler_start proc~profiler_stop mod_profiler::profiler_stop proc~flow_global_two_dots_owned->proc~profiler_stop proc~fatal_error mod_kinds::fatal_error proc~check_mpi~2->proc~fatal_error mpi_wtime mpi_wtime proc~profiler_start->mpi_wtime proc~find_or_create_timer mod_profiler::find_or_create_timer proc~profiler_start->proc~find_or_create_timer proc~profiler_stop->mpi_wtime proc~profiler_stop->proc~find_or_create_timer proc~record_edge mod_profiler::record_edge proc~profiler_stop->proc~record_edge

Called by

proc~~flow_global_two_dots_owned~~CalledByGraph proc~flow_global_two_dots_owned mod_mpi_flow::flow_global_two_dots_owned proc~solve_pressure_correction mod_flow_projection::solve_pressure_correction proc~solve_pressure_correction->proc~flow_global_two_dots_owned proc~advance_projection_step mod_flow_projection::advance_projection_step proc~advance_projection_step->proc~solve_pressure_correction program~lowmach_react_hex lowmach_react_hex program~lowmach_react_hex->proc~advance_projection_step

Source Code

   subroutine flow_global_two_dots_owned(flow, a1, b1, a2, b2, results)
      use mod_profiler, only : profiler_start, profiler_stop
      type(flow_mpi_t), intent(in) :: flow
      real(rk), intent(in) :: a1(:), b1(:), a2(:), b2(:)
      real(rk), intent(out) :: results(2)
      real(rk) :: local_dots(2)
      integer :: c, ierr

      local_dots = zero
      do c = flow%first_cell, flow%last_cell
         local_dots(1) = local_dots(1) + a1(c) * b1(c)
         local_dots(2) = local_dots(2) + a2(c) * b2(c)
      end do

      if (flow%nprocs == 1) then
         results = local_dots
         return
      end if

      call profiler_start('MPI_Communication')
      call MPI_Allreduce(local_dots, results, 2, MPI_DOUBLE_PRECISION, MPI_SUM, flow%comm, ierr)
      call check_mpi(ierr, 'MPI_Allreduce two dots')
      call profiler_stop('MPI_Communication')
   end subroutine flow_global_two_dots_owned