psmile_get_intersect.F90

Go to the documentation of this file.

!-----------------------------------------------------------------------
! Copyright 2006-2010, NEC Europe Ltd., London, UK.
! All rights reserved. Use is subject to OASIS4 license terms.
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: PSMILe_Get_intersect
!
! !INTERFACE:

      subroutine psmile_get_intersect (ninter, nmyint, nnull, &
                                       num_intersect_per_grid, &
                                       num_dummy_intersect_per_grid, &
                                       lastag, ierror)
!
! !USES:
!
      use PRISM_constants
      use PSMILe_common, only : PSMILe_nullify_Enddef_search, num_req_types
      use PSMILe, dummy_interface => PSMILe_get_intersect
      implicit none
!
! !INPUT PARAMETERS:
!
      Integer, Intent (In)            :: ninter

!     Contains total number of intersections found for the
!     actual process.

      Integer, Intent (In)            :: nmyint

!     Contains total number of intersections found which are located
!     also on the actual process.

      Integer, Intent (In)            :: nnull

!     Contains total number of dummy (empty) intersections found
!     which were sent to a destination process

      Integer, Intent (In)            :: num_intersect_per_grid(:), 
                                         num_dummy_intersect_per_grid(:)

      Integer, Intent (In)            :: lastag

!     Message tag used by PSMILe_find_intersect
!
!
! !OUTPUT PARAMETERS:
!
      Integer, Intent (Out)           :: ierror

!     Returns the error code of PSMILe_get_intersect;
!             ierror = 0 : No error
!             ierror > 0 : Severe error
!
! !DEFINED PARAMETERS:
!
      Logical, Parameter           :: new_search = .true.
!
!  TOL must become input parameter
!
      Double Precision, parameter  :: tol = 1d-6
!
! !LOCAL VARIABLES
!
!  ... for intersections sent
!
      Integer                      :: ipart
!
!  ... for grid data sent by a partner process
!
      Integer                      :: i
      Type (Enddef_search)         :: search
      Logical                      :: new_intersection

!
!  ... for communication
!
      Integer                      :: index
      Integer                      :: lstatus (MPI_STATUS_SIZE)
#ifndef PRISM_with_MPI2
      Integer, Allocatable         :: statuses (:, :)
#endif
!
!  ... for finalize messages of global search
!
      Integer                      :: dest, icomp, n
!
      Integer                      :: fin_extra (2)
      Integer, Allocatable         :: fin_requests (:)
!
!  ... for error parameters
!
      Integer, Parameter           :: nerrp = 3
      Integer                      :: ierrp (nerrp)
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_get_intersect" receives the messages sent from
! routine "PSMILe_find_intersect", performs the actions required from
! other processes and performs the search of donor cells.
!
! Doxygen documentation (requires mscgen tool)
!> \file psmile_get_intersect.F90
!! The following is a general overview on the messages exchanged between source\n
!! and target. This routine mainly handles the header messages while the real\n
!! messages are being dealt with in psmile_enddef_action.\n
!! \msc
!!   width = "750";
!!   # entities
!!   Target,Source,"other Sources";
!!   ||| ;
!!   --- [label="psmile_find_intersect:"] ;
!!   Target => Source [label = "intersection (lasttag(100) reqidx(2)"] ;
!!   --- [label="psmile_get_intersect:"] ;
!!   Source rbox Source [label="initialise recv of data\n(psmile_recv_req_subgrid)"] ;
!!   Target <= Source [label="return previouse message as\nacknowledgment (regtag(101) reqidx(1))"] ;
!!   ||| ;
!!   Target => Source [label="send grid data (grdtag(102) reqidx(3))"] ;
!!   Source rbox Source [label="match grids\n(psmile_search_donor_cells)"] ;
!!   ||| ;
!!   Target :: "other Sources" [label="in case of global search for conservative remapping"] ;
!!   Target <= Source [label="header for cell search\n(loctag(110) reqidx(>=num_req_types))"] ;
!!   ||| ;
!!   Target <= Source [label="send information on target cells\nonly partially covered by local data\n(celltag(107) received by psmile_enddef_action_cell)"] ;
!!   ||| ;
!!   Target <= "other Sources" [label="header for cell search (loctag(110) reqidx(>=num_req_types))"] ;
!!   Target <= "other Sources" [label="send information on target cells only partially covered by local data\n(celltag(107) received by psmile_enddef_action_cell)"] ;
!!   ||| ;
!!   Target rbox Target [label="collect all information\non partially covered cells"] ;
!!   ||| ;
!!   Target => Source [label="inform about other processes covering\nthe same cells (celltag+1(108)\nrecevied by psmile_global_search_cell)"] ;
!!   ||| ;
!!   Target => "other Sources" [label="inform about other processes covering the same cells\n(celltag+1(108) recevied by psmile_global_search_cell)"] ;
!!   Target :: "other Sources" [label="end of global for converv. remap."] ;
!!   Target :: "other Sources" [label="in case of global search for points"] ;
!!   Source => "other Sources" [label="request for points not found in\nlocal data (exttag(104) reqidx(4))"] ;
!!   ||| ;
!!   Source <= "other Sources" [label="return results of local search\n(rexttag(105) received by psmile_global_search_*)"] ;
!!   Target :: "other Sources" [label="end of global search for points"] ;
!!   Target <= Source [label="inform about points/cells being\nprovided by directly or through coupler\n(loctag(110) reqidx(>=num_req_types))"] ;
!!   ||| ;
!!   Target <= "other Sources" [label="inform about points/cells being provided by directly or through coupler\n(loctag(110) reqidx(>=num_req_types))"] ;
!!   ||| ;
!!   Target rbox "other Sources" [label="exchange of finish messages (exttag(104) reqidx(4))"] ;
!! \endmsc
!<
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.07.06    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
! $Id: psmile_get_intersect.F90 2804 2010-12-07 10:07:10Z hanke $
! $Author: hanke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_get_intersect.F90 2804 2010-12-07 10:07:10Z hanke $'
!
!----------------------------------------------------------------------

#ifdef VERBOSE
      print 9990, trim(ch_id), ninter, nmyint, nnull

      call psmile_flushstd
#endif /* VERBOSE */
!
!  Initialization
!
      ierror = 0
!
      paction%n      = 0
      paction%ninter = ninter
      paction%lastag = lastag
!     paction%nmyint = nmyint
      paction%nmyint = 0
!
      paction%grid2receive = .false.
      allocate(paction%n_answer2recv_per_grid(Number_of_Grids_allocated))
      paction%n_answer2recv_per_grid(:) = num_intersect_per_grid(:) - num_dummy_intersect_per_grid(:)
      new_intersection = .true.
!
!  n_answer = Number of answers containing requests for grid data received.
!             If no grid data was required, the receiving process doesn't
!             send an answer. This are "nnull" messages.
!
      paction%n_answer = 0
!     paction%n_answer2recv = ninter - nmyint - nnull
      paction%n_answer2recv = ninter - nnull
!
!  nloc_recv = Number of processes which have sent the locations for
!              all grid functions (fields)
!
      paction%nloc_recv = 0
!
!===> Number of outstanding receives for nearest neighbour points
!
      paction%n_selected = 0
!
!===> Get number of finalize messages for global search
!     ??? Feststellen, ob global search specifiziert wurde ?
!
      paction%n_fin = 0
      paction%n_fin2recv = 0
      do icomp = 1, n_act_comp
         paction%n_fin2recv = paction%n_fin2recv + comp_infos(icomp)%size
      end do
!
!  lrequest (1)            = reqtag : Request to send grid data
!  lrequest (2)            = lastag : Receive data on grid intersection
!  lrequest (3)            = grdtag : Receive grid data
!  lrequest (4)            = exttag : Receive request for extra search
!  lrequest (5)            = seltag : Request to receive info on seleted points
!                                     of nearest neighbour search
!
!  lrequest (num_req_types:nreq) = loctag+: Receive data on locations found
!
      paction%nreq = num_req_types - 1 + paction%n_answer2recv
!
      Allocate (paction%lrequest(paction%nreq), stat = ierror)
      if ( ierror /= 0 ) then
         ierrp (1) = paction%nreq
         ierror = PRISM_Error_Alloc
         call psmile_error ( ierror, 'paction%lrequest', &
                             ierrp, 1, __FILE__, __LINE__ )
         return
      endif
!
      if (paction%n_answer2recv > 0) then
         Allocate (paction%loc_messages(msgloc_size, paction%n_answer2recv), &
                   stat = ierror)
         if ( ierror /= 0 ) then
            ierrp (1) = msgloc_size * paction%n_answer2recv
            ierror = PRISM_Error_Alloc
            call psmile_error ( ierror, 'paction%loc_messages', &
                                ierrp, 1, __FILE__, __LINE__ )
            return
         endif
      end if
!
      paction%lrequest (:) = MPI_REQUEST_NULL
!
!===> Set up request for a grid transfer (answer to tag "lastag")
!
!     n_answer2recv - n_answer is the number of messages to be received
!
      if (paction%n_answer < paction%n_answer2recv) then
         call MPI_Irecv (paction%msgreq, nd_msgint, MPI_INTEGER,      &
                         MPI_ANY_SOURCE, reqtag, comm_psmile, &
                         paction%lrequest(1), ierror)
         if ( ierror /= MPI_SUCCESS ) then
            ierrp (1) = ierror
            ierror = PRISM_Error_MPI

            call psmile_error ( ierror, 'MPI_Irecv', &
                                ierrp, 1, __FILE__, __LINE__ )
            return
         endif
#ifdef DEBUG
         print *, ' Posting Irecv request(1) ', paction%lrequest(1), 'with tag ', reqtag
#endif /* DEBUG */
      endif
!
!===> ... Set up request for receive of an extra search request
!
      if (paction%n_fin2recv > 0) then
         call MPI_Irecv (paction%msg_extra, nd_msgextra, MPI_INTEGER, &
                         MPI_ANY_SOURCE, exttag, comm_psmile,         &
                         paction%lrequest(4), ierror)

         if ( ierror /= MPI_SUCCESS ) then
            ierrp (1) = ierror
            ierror = PRISM_Error_MPI

            call psmile_error ( ierror, 'MPI_Irecv', &
                                ierrp, 1, __FILE__, __LINE__ )
            return
         endif
#ifdef DEBUG
         print *, ' Posting Irecv request(4) ', paction%lrequest(4), 'with tag ', exttag
#endif /* DEBUG */
      endif

      call psmile_nullify_enddef_search(search)

!
!-----------------------------------------------------------------------
!     Wait for an action to be performed
!
!     Note: Either
!             lrequest (2)
!             wait for an intersection
!           or
!             lrequest (3) corresponds to grid2receive
!              wait for grid data to be sent for the intersection
!
!           is active
!-----------------------------------------------------------------------
!
      do while ( paction%n          < paction%ninter        .or. &
                 paction%n_answer   < paction%n_answer2recv .or. &
                 paction%nloc_recv  < paction%n_answer2recv .or. &
                 paction%n_selected > 0                     .or. &
                 paction%grid2receive )

! if psmile is ready to receive an intersection and if there still
! intersection messages left to be received
         if (paction%ninter > 0         .and. &
             paction%n < paction%ninter .and. &
             new_intersection) then
!
!===> ... Set up request for receive of an intersection
!
#ifdef PRISM_ASSERTION
            if (paction%lrequest(2) /= MPI_REQUEST_NULL) then
               call psmile_assert ( __FILE__, __LINE__, &
                                    'paction%lrequest(2) is not finished !')
            endif
#endif
#define CONS_REMAP_DEADLOCK_WORKAROUND
#ifndef CONS_REMAP_DEADLOCK_WORKAROUND
            call MPI_Irecv (paction%msgint, nd_msgint, MPI_INTEGER,      &
                           MPI_ANY_SOURCE, paction%lastag, comm_psmile, &
                           paction%lrequest(2), ierror)
#else
            call MPI_Irecv (paction%msgint, nd_msgint, MPI_INTEGER,      &
                           maxval(paction%intersect_ranks), paction%lastag, comm_psmile, &
                           paction%lrequest(2), ierror)
            paction%intersect_ranks(maxloc(paction%intersect_ranks)) = -1
#endif

            if ( ierror /= MPI_SUCCESS ) then
               ierrp (1) = ierror
               ierror = PRISM_Error_MPI

               call psmile_error ( ierror, 'MPI_Irecv', &
                                 ierrp, 1, __FILE__, __LINE__ )
               return
            endif
#ifdef DEBUG
         print *, ' Posting Irecv request(2) ', paction%lrequest(2), 'with tag ', paction%lastag
#endif /* DEBUG */
         endif

!
!===> ... Wait for a request
!
#ifdef DEBUG
         print 9960, trim(ch_id), paction%n, paction%ninter,                  &
                                  paction%n_answer, paction%n_answer2recv,    &
                                  paction%grid2receive,                       &
                                  paction%nloc_recv, paction%n_answer2recv,   &
                                  paction%n_selected

         call psmile_flushstd
         do i = 1, paction%nreq
            print *, '    list of requests ', i, paction%lrequest(i)
         enddo
         call psmile_flushstd
#endif /* DEBUG */
         call MPI_Waitany (paction%nreq, paction%lrequest, index, lstatus, ierror)

         if ( ierror /= MPI_SUCCESS ) then
            ierrp (1) = ierror
            ierrp (2) = lstatus (MPI_SOURCE)
            ierrp (3) = lstatus (MPI_TAG)

            ierror = PRISM_Error_MPI

            print *, 'lstatus', lstatus
            call psmile_error ( ierror, 'MPI_Waitany', &
                                ierrp, 3, __FILE__, __LINE__ )
            return
         endif

#ifdef DEBUG
         print *, trim(ch_id), ": n, ninter, sender, index", &
                  paction%n, paction%ninter, lstatus (MPI_SOURCE), index
         call psmile_flushstd
#endif /* DEBUG */

         if (index == MPI_UNDEFINED) then
#ifdef PRISM_ASSERTION
            call psmile_assert ( __FILE__, __LINE__, &
                                 'request list is empty')
#endif /* PRISM_ASSERTION */
            exit
         endif
!
!===> ... A request was completed
!
         if (index == 3) then
!
!-----------------------------------------------------------------------
!        index = 3: Grid data was received (tag ``grdtag'')
!                   Note: The first request of "recv_req" was already
!                         fulfilled.
!
!  Fehlt hier nicht noch ein counter fuer das while
!-----------------------------------------------------------------------
!
#ifdef PRISM_ASSERTION
            if (lstatus(MPI_TAG) /= grdtag) then
               call psmile_assert ( __FILE__, __LINE__, &
                                   'wrong message tag; must be grdtag')
            endif
#endif

#ifdef PRISM_with_MPI2
            call MPI_Waitall ((ndim_3d+2)*search%npart-1, paction%recv_req(2:), &
                              MPI_STATUSES_IGNORE, ierror)
#else
            Allocate (statuses (MPI_STATUS_SIZE,search%npart*(ndim_3d+2)-1), &
                      stat = ierror)
            if ( ierror /= 0 ) then
               ierrp (1) = MPI_STATUS_SIZE * (search%npart*(ndim_3d+2)-1)
               ierror = PRISM_Error_Alloc
               call psmile_error ( ierror, 'statuses', &
                                   ierrp, 1, __FILE__, __LINE__ )
               return
            endif

            call MPI_Waitall ((ndim_3d+2)*search%npart-1, paction%recv_req(2:), &
                              statuses, ierror)
#endif
            if ( ierror /= MPI_SUCCESS ) then
               ierrp (1) = ierror

               ierror = PRISM_Error_MPI

               call psmile_error ( ierror, 'MPI_Waitall', &
                                   ierrp, 1, __FILE__, __LINE__ )
               return
            endif
!
#ifndef PRISM_with_MPI2
            Deallocate (statuses)
#endif
            Deallocate (paction%recv_req)

            paction%grid2receive = .false.
!
!           ... Search donor cells
!
            call psmile_search_donor_cells (search, tol, ierror )
            if (ierror > 0) return
!
! ... to do

!
!           ... Deallocate arrays allocated
!               The arrays and vectors were allocated in
!               PSMILe_Recv_req_subgrid and PSMILe_Recv_req_coords_????
!
            if (search%datatype == MPI_REAL) then
               do ipart = 1, search%npart
                  do i = 1, ndim_3d
                     Deallocate (search%search_real(i, ipart)%vector)
                  end do
               end do
!
               Deallocate (search%search_real)
!
            else if (search%datatype == MPI_DOUBLE_PRECISION) then
               do ipart = 1, search%npart
                  do i = 1, ndim_3d
                     Deallocate (search%search_dble(i, ipart)%vector)
                  end do
               end do
!
               Deallocate (search%search_dble)
!
#if defined ( PRISM_QUAD_TYPE )
            else if (search%datatype == MPI_REAL16) then
               do ipart = 1, search%npart
                  do i = 1, ndim_3d
                  Deallocate (search%search_quad(i, ipart)%vector)
                  end do
               end do
!
               Deallocate (search%search_quad)
#endif
            endif
!
            if (Associated(search%search_mask)) then
               do ipart = 1, search%npart
                  Deallocate(search%search_mask(ipart)%vector)
               end do
!
               Deallocate (search%search_mask)
            endif
!
            if (Associated(search%global_index)) then
               do ipart = 1, search%npart
                  Deallocate(search%global_index(ipart)%vector)
               end do
               Deallocate (search%global_index)
            endif

            Deallocate (search%msgint)
            Deallocate (search%dims)
            Deallocate (search%shape)
            Deallocate (search%range)
!
!===> ... Set up request for receive of an intersection
!
            new_intersection = .true.
!
         else

            new_intersection = .false.
!
!===> ... Perform default action
!
            call psmile_enddef_action (search, index, lstatus, ierror)
            if (ierror > 0) return
         endif
!
      end do ! do while (paction%n < paction%ninter .and. ...)
!
!===> All requests except for global search are fulfilled.
!     Send info to processes.
!     TODO: Send message over a tree in order to decrease the
!           number of sends.
!
      if (paction%n_fin2recv > 0) then

         Allocate (fin_requests (paction%n_fin2recv), stat = ierror)
         if ( ierror /= 0 ) then
            ierrp (1) = paction%n_fin2recv
            ierror = PRISM_Error_Alloc
            call psmile_error ( ierror, 'fin_requests', &
                                ierrp, 1, __FILE__, __LINE__ )
            return
         endif
!
         fin_extra (1) = PSMILe_Finalize_extra_search
         fin_extra (2) = 1 ! Number of processes finalized
!
         n = 0
         do icomp = 1, n_act_comp
            do i = 1, comp_infos(icomp)%size
               dest = comp_infos(icomp)%psmile_ranks(i)

               n = n + 1

               call MPI_Isend (fin_extra, 2, MPI_INTEGER, &
                               dest, exttag, comm_psmile, &
                               fin_requests(n), ierror)
!
               if ( ierror /= MPI_SUCCESS ) then
                  ierrp (1) = ierror
                  ierrp (2) = dest
                  ierrp (3) = exttag
                  ierror = PRISM_Error_Send

                  call psmile_error ( ierror, 'MPI_Isend', &
                                 ierrp, 3, __FILE__, __LINE__ )
                  return
               endif
#ifdef DEBUGX
               print *, ' Posting Isend request ', fin_requests(n), &
                        ' with tag ', exttag, ' to dest ', dest
               call psmile_flushstd
#endif /* DEBUGX */
            end do ! i = 1, comp_infos(icomp)%size
         end do ! icomp
!
!===> Wait that all finish messages of global search has been
!     received.
!     index 4: Message of extra_search (inclusive finish messages
!              with command "PSMILe_Finalize_extra_search")
!     index 5: Receive info on seleted points of nearest neighbour search
!
         do while (paction%n_fin < paction%n_fin2recv .or. &
                   paction%n_selected > 0)
#ifdef DEBUGX
            print 9950, trim(ch_id), paction%n_fin, paction%n_fin2recv
            do i = 4, 5
               print *, '    list of request ', i, paction%lrequest(i)
            enddo
            call psmile_flushstd
#endif /* DEBUGX */
            call MPI_Waitany (2, paction%lrequest(4:5), index, lstatus, ierror)

            if ( ierror /= MPI_SUCCESS ) then
               ierrp (1) = ierror
               ierrp (2) = lstatus (MPI_SOURCE)
               ierrp (3) = lstatus (MPI_TAG)

               ierror = PRISM_Error_MPI

               print *, 'lstatus', lstatus
               call psmile_error ( ierror, 'MPI_Waitany', &
                                   ierrp, 3, __FILE__, __LINE__ )
               return
            endif
#ifdef DEBUGX
            print *, ' returned index ', index
#endif /* DEBUGX */
            index = index + (4-1)
!
            call psmile_enddef_action (search, index, lstatus, ierror)
            if (ierror > 0) return
!
         enddo ! while (paction%n_fin < paction%n_fin2recv .or. paction%n_selected > 0)
!
!===> Wait that all finish messages have been sent
!
#ifdef PRISM_with_MPI2
         call MPI_Waitall (paction%n_fin2recv, fin_requests, &
                           MPI_STATUSES_IGNORE, ierror)
#else
         Allocate (statuses(MPI_STATUS_SIZE, paction%n_fin2recv), &
                   stat = ierror)
         if (ierror > 0) then
            ierrp (1) = ierror
            ierrp (2) = paction%n_fin2recv * MPI_STATUS_SIZE

            ierror = PRISM_Error_Alloc
            call psmile_error ( ierror, 'statuses', &
                                ierrp, 2, __FILE__, __LINE__ )
            return
         endif

         call MPI_Waitall (paction%n_fin2recv, fin_requests, statuses, ierror)
!
         Deallocate (statuses)
#endif
!
         if ( ierror /= MPI_SUCCESS ) then
            ierrp (1) = ierror

            ierror = PRISM_Error_MPI

            call psmile_error ( ierror, 'MPI_Waitall', &
                                ierrp, 1, __FILE__, __LINE__ )
            return
         endif

         Deallocate (fin_requests)
!
      endif ! paction%n_fin2recv > 0
!
!=======================================================================
!     Free
!=======================================================================
!
      if (paction%n_answer2recv > 0) then
         Deallocate (paction%loc_messages)
      endif
!
      Deallocate (paction%lrequest)
      Deallocate (paction%n_answer2recv_per_grid)
!
!===> All done
!
#ifdef DEBUG
!
!  Print send infos
!
      do i = 1, Number_of_Methods_allocated
         if (Methods(i)%status  == PSMILe_status_defined) then
            call psmile_print_send_info (i, 0, 'end of psmile_get_intersect')
         endif
      end do
#endif /* DEBUG */

#ifdef VERBOSE
      print 9980, trim(ch_id), ierror
      call psmile_flushstd
#endif /* VERBOSE */

!
!  Formats:
!

#ifdef VERBOSE

9990 format (1x, a, ': psmile_get_intersect: ninter', i4, &
                    ', nmyint', i4, ', nnull', i4)
9980 format (1x, a, ': psmile_get_intersect: eof ierror =', i3)

#endif /* VERBOSE */

#ifdef DEBUG
9960 format (1x, a, ': psmile_get_intersect: calling waitany: n, ninter', 2i4, ';',&
            /2x,    'n_answer, n_answer2recv', 2i4, &
                    '; grid2receive ', l1, &
                    '; nloc_recv, nreq-num_req_types+1', 2i4, &
                    '; n_selected', i4)
9950 format (1x, a, ': psmile_get_intersect: calling wait: n_fin, n_fin2recv', 2i4)
#endif /* DEBUG */
      end subroutine PSMILe_get_intersect