Changeset 7992 for branches/UKMO/dev_r5518_obs_oper_update/NEMOGCM/NEMO/OPA_SRC/OBS/obs_mpp.F90

Timestamp:

2017-05-02T13:21:57+02:00 (7 years ago)

Author:

jwhile

Message:

This version of the code seems to work correctly after some bug fixes

File:

• branches/UKMO/dev_r5518_obs_oper_update/NEMOGCM/NEMO/OPA_SRC/OBS/obs_mpp.F90 (modified) (6 diffs)

branches/UKMO/dev_r5518_obs_oper_update/NEMOGCM/NEMO/OPA_SRC/OBS/obs_mpp.F90

@@ -7 +7 @@
    !!             -   ! 2006-05  (K. Mogensen)  Reformatted
    !!             -   ! 2008-01  (K. Mogensen)  add mpp_global_max
+   !!            3.6  ! 2015-01  (J. Waters) obs_mpp_find_obs_proc 
+   !!                            rewritten to avoid global arrays
    !!----------------------------------------------------------------------
 #  define mpivar mpi_double_precision
@@ -12 +14 @@
    !! obs_mpp_bcast_integer : Broadcast an integer array from a processor to all processors
    !! obs_mpp_max_integer   : Find maximum on all processors of each value in an integer on all processors
-   !! obs_mpp_find_obs_proc : Find processors which should hold the observations
+   !! obs_mpp_find_obs_proc : Find processors which should hold the observations, avoiding global arrays
    !! obs_mpp_sum_integers  : Sum an integer array from all processors
    !! obs_mpp_sum_integer   : Sum an integer from all processors
@@ -96 +98 @@
       !
       INTEGER :: ierr 
-      INTEGER, DIMENSION(kno) ::   ivals
-      !
-INCLUDE 'mpif.h'
-      !!----------------------------------------------------------------------
+      INTEGER, DIMENSION(:), ALLOCATABLE ::   ivals
+      !
+INCLUDE 'mpif.h'
+      !!----------------------------------------------------------------------
+
+      ALLOCATE( ivals(kno) )
 
       ! Call the MPI library to find the maximum across processors
@@ -105 +109 @@
          &                mpi_max, mpi_comm_opa, ierr )
       kvals(:) = ivals(:)
+
+      DEALLOCATE( ivals )
 #else
       ! no MPI: empty routine
@@ -111 +117 @@
 
 
-   SUBROUTINE obs_mpp_find_obs_proc( kobsp, kobsi, kobsj, kno )
-      !!----------------------------------------------------------------------
-      !!               ***  ROUTINE obs_mpp_find_obs_proc ***
-      !!          
-      !! ** Purpose : From the array kobsp containing the results of the grid
+   SUBROUTINE obs_mpp_find_obs_proc( kobsp,kno )
+      !!----------------------------------------------------------------------
+      !!               ***  ROUTINE obs_mpp_find_obs_proc  ***
+      !!         
+      !! ** Purpose : From the array kobsp containing the results of the
       !!              grid search on each processor the processor return a
       !!              decision of which processors should hold the observation.
       !!
-      !! ** Method  : A temporary 2D array holding all the decisions is
-      !!              constructed using mpi_allgather on each processor.
-      !!              If more than one processor has found the observation
-      !!              with the observation in the inner domain gets it
-      !!
-      !! ** Action  : This does only work for MPI. 
+      !! ** Method  : Synchronize the processor number for each obs using
+      !!              obs_mpp_max_integer. If an observation exists on two 
+      !!              processors it will be allocated to the lower numbered
+      !!              processor.
+      !!
+      !! ** Action  : This does only work for MPI.
       !!              It does not work for SHMEM.
       !!
@@ -130 +136 @@
       !!----------------------------------------------------------------------
       INTEGER                , INTENT(in   ) ::   kno
-      INTEGER, DIMENSION(kno), INTENT(in   ) ::   kobsi, kobsj
       INTEGER, DIMENSION(kno), INTENT(inout) ::   kobsp
       !
 #if defined key_mpp_mpi
       !
-      INTEGER :: ji
-      INTEGER :: jj
-      INTEGER :: size
-      INTEGER :: ierr
-      INTEGER :: iobsip
-      INTEGER :: iobsjp
-      INTEGER :: num_sus_obs
-      INTEGER, DIMENSION(kno) ::   iobsig, iobsjg
-      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   iobsp, iobsi, iobsj
-      !!
-INCLUDE 'mpif.h'
-      !!----------------------------------------------------------------------
-
-      !-----------------------------------------------------------------------
-      ! Call the MPI library to find the maximum accross processors
-      !-----------------------------------------------------------------------
-      CALL mpi_comm_size( mpi_comm_opa, size, ierr )
-      !-----------------------------------------------------------------------
-      ! Convert local grids points to global grid points
-      !-----------------------------------------------------------------------
+      !
+      INTEGER :: ji, isum
+      INTEGER, DIMENSION(:), ALLOCATABLE ::   iobsp
+      !!
+      !!
+
+      ALLOCATE( iobsp(kno) )
+
+      iobsp(:)=kobsp(:)
+
+      WHERE( iobsp(:) == -1 )
+         iobsp(:) = 9999999
+      END WHERE
+
+      iobsp(:)=-1*iobsp(:)
+
+      CALL obs_mpp_max_integer( iobsp, kno )
+
+      kobsp(:)=-1*iobsp(:)
+
+      isum=0
       DO ji = 1, kno
-         IF ( ( kobsi(ji) >= 1 ) .AND. ( kobsi(ji) <= jpi ) .AND. &
-            & ( kobsj(ji) >= 1 ) .AND. ( kobsj(ji) <= jpj ) ) THEN
-            iobsig(ji) = mig( kobsi(ji) )
-            iobsjg(ji) = mjg( kobsj(ji) )
-         ELSE
-            iobsig(ji) = -1
-            iobsjg(ji) = -1
+         IF ( kobsp(ji) == 9999999 ) THEN
+            isum=isum+1
+            kobsp(ji)=-1
          ENDIF
-      END DO
-      !-----------------------------------------------------------------------
-      ! Get the decisions from all processors
-      !-----------------------------------------------------------------------
-      ALLOCATE( iobsp(kno,size) )
-      ALLOCATE( iobsi(kno,size) )
-      ALLOCATE( iobsj(kno,size) )
-      CALL mpi_allgather( kobsp, kno, mpi_integer, &
-         &                iobsp, kno, mpi_integer, &
-         &                mpi_comm_opa, ierr )
-      CALL mpi_allgather( iobsig, kno, mpi_integer, &
-         &                iobsi, kno, mpi_integer, &
-         &                mpi_comm_opa, ierr )
-      CALL mpi_allgather( iobsjg, kno, mpi_integer, &
-         &                iobsj, kno, mpi_integer, &
-         &                mpi_comm_opa, ierr )
-
-      !-----------------------------------------------------------------------
-      ! Find the processor with observations from the lowest processor 
-      ! number among processors holding the observation.
-      !-----------------------------------------------------------------------
-      kobsp(:) = -1
-      num_sus_obs = 0
-      DO ji = 1, kno
-         DO jj = 1, size
-            IF ( ( kobsp(ji) == -1 ) .AND. ( iobsp(ji,jj) /= -1 ) ) THEN
-               kobsp(ji) = iobsp(ji,jj)
-               iobsip = iobsi(ji,jj)
-               iobsjp = iobsj(ji,jj)
-            ENDIF
-            IF ( ( kobsp(ji) /= -1 ) .AND. ( iobsp(ji,jj) /= -1 ) ) THEN
-               IF ( ( iobsip /= iobsi(ji,jj) ) .OR. &
-                  & ( iobsjp /= iobsj(ji,jj) ) ) THEN
-                  IF ( ( kobsp(ji) < 1000000 ) .AND. &
-                     & ( iobsp(ji,jj) < 1000000 ) ) THEN
-                     num_sus_obs=num_sus_obs+1
-                  ENDIF
-               ENDIF
-               IF ( mppmap(iobsip,iobsjp) /= ( kobsp(ji)+1 ) ) THEN
-                  IF ( ( iobsi(ji,jj) /= -1 ) .AND. &
-                     & ( iobsj(ji,jj) /= -1 ) ) THEN
-                     IF ((mppmap(iobsi(ji,jj),iobsj(ji,jj)) == (iobsp(ji,jj)+1))&
-                        & .OR. ( iobsp(ji,jj) < kobsp(ji) ) ) THEN
-                        kobsp(ji) = iobsp(ji,jj)
-                        iobsip = iobsi(ji,jj)
-                        iobsjp = iobsj(ji,jj)
-                     ENDIF
-                  ENDIF
-               ENDIF
-            ENDIF
-         END DO
-      END DO
-      IF (lwp) WRITE(numout,*) 'Number of suspicious observations: ',num_sus_obs
-
-      DEALLOCATE( iobsj )
-      DEALLOCATE( iobsi )
+      ENDDO
+
+
+      IF ( isum > 0 ) THEN
+         IF (lwp) WRITE(numout,*) isum, ' observations failed the grid search.'
+         IF (lwp) WRITE(numout,*)'If ln_grid_search_lookup=.TRUE., try reducing grid_search_res'
+      ENDIF
+
       DEALLOCATE( iobsp )
+
 #else
       ! no MPI: empty routine
-#endif
-      !
+#endif     
+      
    END SUBROUTINE obs_mpp_find_obs_proc