psmile_recv_req_mask.F90

Go to the documentation of this file.
00001 !-----------------------------------------------------------------------
00002 ! Copyright 2006-2010, NEC Europe Ltd., London, UK.
00003 ! All rights reserved. Use is subject to OASIS4 license terms.
00004 !-----------------------------------------------------------------------
00005 !BOP
00006 !
00007 ! !ROUTINE: PSMILe_Recv_req_mask
00008 !
00009 ! !INTERFACE:
00010 
00011       subroutine psmile_recv_req_mask (sender, tag,      &
00012                                        search, recv_req, &
00013                                        receive_mask, new_search, ierror)
00014 !
00015 ! !USES:
00016 !
00017       use PRISM_constants
00018 !
00019       use PSMILe, dummy_interface => PSMILe_Recv_req_mask
00020 
00021       implicit none
00022 !
00023 ! !INPUT PARAMETERS:
00024 !
00025       Integer, Intent (In)                :: sender
00026 
00027 !     Specifies the sender of the subgrid.
00028 
00029       Integer, Intent (In)                :: tag
00030 
00031 !     Specifies the message tag used
00032 
00033       Logical, Intent (In)                :: receive_mask
00034 
00035 !     Should mask array be received ?
00036 !
00037       Logical, Intent (In)                :: new_search
00038 
00039 !     Option:
00040 !     new_search = .true.  : New search; fill search structure
00041 !     new_search = .fasle. : Get only data for a new field
00042 !
00043 ! !INPUT/OUTPUT PARAMETERS:
00044 !
00045       Type (Enddef_search), Intent (InOut) :: search
00046 
00047 !     Data on the points to be searched
00048 !
00049 ! !OUTPUT PARAMETERS:
00050 !
00051       integer, Intent (Out)               :: recv_req (ndim_3d+1, search%search_data%npart)
00052 
00053 !     Receive requests of non-blocking receives
00054 !     The receive requests are stored in recv_req (ndim_3d+1, :).
00055 
00056       integer, Intent (Out)               :: ierror
00057 
00058 !     Returns the error code of PSMILe_Recv_req_mask;
00059 !             ierror = 0 : No error
00060 !             ierror > 0 : Severe error
00061 !
00062 ! !LOCAL VARIABLES
00063 !
00064 !     ... for the intersections
00065 !
00066       integer                      :: ipart
00067       integer                      :: len_mask
00068 !
00069       logical, save                :: dummy_mask
00070 !
00071 !     ... for error parameters
00072 !
00073       integer, parameter           :: nerrp = 3
00074       integer                      :: ierrp (nerrp)
00075 !
00076 ! !DESCRIPTION:
00077 !
00078 ! Subroutine "PSMILe_Recv_req_mask" receives the mask requested
00079 ! from the sending process. The mask is sent after the corresponding
00080 ! request was received by routine "PSMILe_Get_intersect".
00081 !
00082 ! One receive is submitted at least in order to ensure progress
00083 ! in psmile_get_intersect.
00084 !
00085 ! !REVISION HISTORY:
00086 !
00087 !   Date      Programmer   Description
00088 ! ----------  ----------   -----------
00089 ! 03.07.10    H. Ritzdorf  created
00090 !
00091 !EOP
00092 !----------------------------------------------------------------------
00093 !
00094 ! $Id: psmile_recv_req_mask.F90 3248 2011-06-23 13:03:19Z coquart $
00095 ! $Author: coquart $
00096 !
00097    Character(len=len_cvs_string), save :: mycvs = 
00098        '$Id: psmile_recv_req_mask.F90 3248 2011-06-23 13:03:19Z coquart $'
00099 !
00100 !----------------------------------------------------------------------
00101 !
00102 !  Initialization
00103 !
00104 #ifdef VERBOSE
00105       print 9990, trim(ch_id), sender, search%search_data%npart, receive_mask
00106 
00107       call psmile_flushstd
00108 #endif /* VERBOSE */
00109 !
00110       ierror = 0
00111 !
00112       recv_req (:, :) = MPI_REQUEST_NULL
00113 !
00114       if (receive_mask) then
00115 !
00116 ! ====================================================================
00117 !     Get dimensions of arrays to be allocated and
00118 !     allocate the arrays.
00119 ! ====================================================================
00120 
00121          do ipart = 1, search%search_data%npart
00122 
00123 
00124 !        search_dims(:) = search%search_data%range (2,:, ipart) - search%search_data%range (1,:, ipart) + 1
00125 !
00126          len_mask = (search%search_data%range(2,1, ipart) - search%search_data%range(1,1, ipart) + 1) * &
00127                     (search%search_data%range(2,2, ipart) - search%search_data%range(1,2, ipart) + 1) * &
00128                     (search%search_data%range(2,3, ipart) - search%search_data%range(1,3, ipart) + 1)
00129 
00130 #ifdef DEBUG
00131            PRINT*, 'Ipart in recv_req_mask :',ipart
00132            PRINT*, 'search%search_data%range(1,1, ipart),search%search_data%range(2,1, ipart)', &
00133                    search%search_data%range(1,1, ipart),search%search_data%range(2,1, ipart)
00134            PRINT*, 'search%search_data%range(1,2, ipart),search%search_data%range(2,2, ipart)', &
00135                    search%search_data%range(1,2, ipart),search%search_data%range(2,2, ipart)
00136            PRINT*, 'search%search_data%range(1,3, ipart),search%search_data%range(2,3, ipart)', &
00137                    search%search_data%range(1,3, ipart),search%search_data%range(2,3, ipart)
00138            PRINT*, 'len_mask :',len_mask
00139            call psmile_flushstd
00140 #endif
00141 !
00142 !===> ... Post receive for the mask if required
00143 !
00144          if (new_search .or. .not. Associated (search%search_mask(ipart)%vector)) &
00145             then
00146 
00147             Allocate (search%search_mask(ipart)%vector(len_mask), STAT = ierror)
00148             if ( ierror > 0 ) then
00149                ierrp (1) = ierror
00150                ierrp (2) = len_mask
00151                ierror = PRISM_Error_Alloc
00152                call psmile_error ( ierror, 'search%search_mask(ipart)%vector', &
00153                                    ierrp, 2, __FILE__, __LINE__ )
00154                return
00155             endif
00156          endif
00157 
00158          call MPI_Irecv (search%search_mask(ipart)%vector, len_mask, &
00159                          MPI_LOGICAL, &
00160                          sender, tag, comm_psmile, recv_req (ndim_3d+1, ipart), ierror)
00161          if (ierror /= MPI_SUCCESS) then
00162             ierrp (1) = ierror
00163             ierrp (2) = sender
00164             ierrp (3) = tag
00165             ierror = PRISM_Error_Recv
00166 
00167             call psmile_error (ierror, 'MPI_Irecv of mask', &
00168                                ierrp, 3, __FILE__, __LINE__ )
00169             return
00170          endif
00171          end do ! ipart
00172 
00173       else
00174 
00175 !
00176 !        ... Receive is not necessary; free previously allocated masks
00177 !
00178 
00179          if ( Associated(search%search_mask) .and. .not. new_search) then
00180                do ipart = 1, search%search_data%npart
00181                if (Associated (search%search_mask(ipart)%vector)) then
00182                   Deallocate  (search%search_mask(ipart)%vector)
00183                endif
00184 
00185                Nullify (search%search_mask(ipart)%vector)
00186                end do
00187          endif
00188 
00189 !
00190 !        ... Setup dummy receive in order to get progress in psmile_get_intersect
00191 !
00192          call MPI_Irecv (dummy_mask, 0, MPI_LOGICAL, &
00193                          sender, tag, comm_psmile, recv_req (ndim_3d+1, 1), &
00194                          ierror)
00195          if (ierror /= MPI_SUCCESS) then
00196             ierrp (1) = ierror
00197             ierrp (2) = sender
00198             ierrp (3) = tag
00199             ierror = PRISM_Error_Recv
00200 
00201             call psmile_error (ierror, 'MPI_Irecv of mask', &
00202                                ierrp, 3, __FILE__, __LINE__ )
00203             return
00204          endif
00205 
00206       endif
00207 
00208 #ifdef DEBUGX
00209       PRINT*, 'In recv_req_mask : recv_req (ndim_3d+1, :):',recv_req (ndim_3d+1, 1)
00210       call psmile_flushstd
00211 #endif
00212 
00213 !
00214 !===> All done
00215 !
00216 #ifdef VERBOSE
00217       print 9980, trim(ch_id), ierror
00218 
00219       call psmile_flushstd
00220 #endif /* VERBOSE */
00221 !
00222 !  Formats:
00223 !
00224 9990 format (1x, a, ': psmile_recv_req_mask: ', &
00225                     'sender =', i4, ', npart =', i3, ', receive_mask', l2)
00226 9980 format (1x, a, ': psmile_recv_req_mask: eof ierror =', i4)
00227 
00228       end subroutine PSMILe_Recv_req_mask

Generated on 1 Dec 2011 for Oasis4 by  doxygen 1.6.1