psmile_hash_extra.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_Hash_extra
!
! !INTERFACE:

      subroutine psmile_hash_extra (search, locs, hash, nlocs, &
                    mask_array, mask_shape, mask_available,    &
                    grid_valid_shape, ierror)
!
! !USES:
!
      use PRISM
!
      use PSMILe, dummy_interface => PSMILe_Hash_extra

      Implicit none
!
! !INPUT PARAMETERS:
!
      Integer,            Intent (in) :: nlocs
!
!     Number of points (interpolation bases) found
!
      Integer,            Intent (In) :: locs (ndim_3d, nlocs)
!
!     Info on the component in which the donor cells
!     should be searched.

      Integer, Intent (In)            :: mask_shape (2, ndim_3d)

!     Dimension of (method) mask array

      Logical, Intent (In)            :: mask_array (mask_shape (1,1): 
                                                     mask_shape (2,1), 
                                                     mask_shape (1,2): 
                                                     mask_shape (2,2), 
                                                     mask_shape (1,3): 
                                                     mask_shape (2,3))
!     Mask of the method

      Logical, Intent (In)            :: mask_available

!     Is mask specified in array "mask_array" ?
!     mask_available = .false. : Mask is not available
!                      .true.  : Mask is     available

      Integer,                     Intent (In)    :: grid_valid_shape (2, 
                                                             ndim_3d)
!
!     Specifies the valid block shape for the "ndim_3d"-dimensional block
!     (without halo/overlap regions) with
!
!     grid_valid_shape (1,i) = Lowest  "i"-th index of block
!     grid_valid_shape (2,i) = Highest "i"-th index of block
!
! !INPUT/OUTPUT PARAMETERS:
!
      Type (Enddef_global_search), Intent (InOut) :: search

!     Data on the points to be searched

! !OUTPUT PARAMETERS:
!
      Integer,                     Intent (Out)   :: hash (nlocs)
!
!     Hash value used (can be used as masked for further computations)
!     with
!     hash (i) > 0  : Point "locs(:, i)" was put into list of points
!                     to be sent.
!     hash (i) == 0 : Point "locs(:, i)" is masked out
!     hash (i) == -1: Point "locs(:, i)" is outside
!
      Integer,                     Intent (Out)   :: ierror

!     Returns the error code of PSMILe_Search_donor_extra_off;
!             ierror = 0 : No error
!             ierror > 0 : Severe error
!
! !LOCAL PARAMETERS:
!
! outside_ng   = Hash value for a point which is outside the block region
! masked_out   = Hash value for a point which is masked out
!
      Integer, Parameter           :: nd_hliste = 1024
!
      Integer, Parameter           :: outside_ng = -1
      Integer, Parameter           :: masked_out = 0
!
! !LOCAL VARIABLES
!
!     ... loop variables
!
      Integer                      :: i, k, l, n, nf
      Integer                      :: nbeg, nend
!
!     ... 
!
      Integer                      :: length (ndim_3d)
!
!     ... For searched points of control volume
!
      Integer, Allocatable         :: locs2liste (:)
      Integer, Allocatable         :: hash_liste (:)
      Integer, Allocatable         :: liste (:)
      Integer, Allocatable         :: next (:)
!
!     ... Stored in search
!
      Integer                      :: n_liste, n_found
!
      Integer, Pointer             :: neigh (:, :)
      Integer, Pointer             :: index_found (:)
!
!     ... for error parameters
!
      Integer, parameter           :: nerrp = 1
      Integer                      :: ierrp (nerrp)
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_Hash_extra" creates compact list of points
! which were found in the extra global search (i.e. which have to be sent
! to the requesting process) and stores the corresponding
! data in "search".
!
! Points which are outside the valid shape "grid_valid_shape" are
! not included to the compact list. Points which are masked out
! (if "mask_available" is true) are transferred and the corresponding
! hash value "hash(*) = masked_out = 0" is returned.
!
! The routine creates
! (*) a vector "index_found (1:n_found) which maps the sequence of 
!     "n_found" points to the compact list "1:n_liste" which contains
!     each found point only once.
! (*) the compact list "neigh_3d(1:n_liste)" which contains the
!     3d logical indices of the compact list.
! 
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.11.06    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
! $Id: psmile_hash_extra.F90 2082 2009-10-21 13:31:19Z hanke $
! $Author: hanke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_hash_extra.F90 2082 2009-10-21 13:31:19Z hanke $'
!
!----------------------------------------------------------------------
!
!  Initialization
!
#ifdef VERBOSE
      print 9990, trim(ch_id)

      call psmile_flushstd
#endif /* VERBOSE */
!
      ierror = 0
!
      length (1:ndim_3d) = grid_valid_shape(2,1:ndim_3d) - &
                           grid_valid_shape(1,1:ndim_3d) + 1
!
#ifdef PRISM_ASSERTION
!
!===> Internal control
!
      if (nlocs <= 0) then
         call psmile_assert ( __FILE__, __LINE__, &
                             'Number of points to be searched == 0 ?!?')
      endif
!
      if (mask_available) then
         if (grid_valid_shape (1,1) < mask_shape (1,1) .or. &
             grid_valid_shape (2,1) > mask_shape (2,1) .or. &
             grid_valid_shape (1,2) < mask_shape (1,2) .or. &
             grid_valid_shape (2,2) > mask_shape (2,2) .or. &
             grid_valid_shape (1,3) < mask_shape (1,3) .or. &
             grid_valid_shape (2,3) > mask_shape (2,3)) then
            print *, "grid_valid_shape", grid_valid_shape
            print *, "mask_shape      ", mask_shape

            call psmile_assert ( __FILE__, __LINE__, &
                    'grid_valid_shape must be in mask_shape')

         endif
      endif
#endif
!
!===> Compute hash value (index in fortran order) for each location
!     to be searched.
!     The hash value is always greater than 0.
!
!cdir vector
         do n = 1, nlocs
         hash(n) = ((locs (3,n) - grid_valid_shape(1,3))  * length(2) + &
                    (locs (2,n) - grid_valid_shape(1,2))) * length(1) + &
                    (locs (1,n) - grid_valid_shape(1,1)) + 1
         end do
!
!===> Remove points which are outside of block.
!     Set hash value to negative value in this case.
!
!cdir vector
         do n = 1, nlocs
         if (locs(1,n) < grid_valid_shape(1,1) .or. &
             locs(1,n) > grid_valid_shape(2,1) .or. &
             locs(2,n) < grid_valid_shape(1,2) .or. &
             locs(2,n) > grid_valid_shape(2,2) .or. &
             locs(3,n) < grid_valid_shape(1,3) .or. &
             locs(3,n) > grid_valid_shape(2,3)) then
            hash (n) = outside_ng
         endif
         end do
!
!===> Control Mask values
!     If mask value if not true, set hash value to 0 = masked_out.
!
      if (mask_available) then
!cdir vector
            do n = 1, nlocs
            if (hash (n) > 0) then
               if (.not. mask_array (locs(1, n), &
                                     locs(2, n), &
                                     locs(3, n)) ) then
                  hash (n) = masked_out
               endif
            endif
            end do
      endif
!
#ifdef DEBUGX
      print *, 'hash values and locs ', grid_valid_shape
         do n = 1, nlocs
         print *, 'n, hash(n), locs(:,n) ', &
                   n, hash(n), locs(:,n)
         end do
#endif
!
!------------------------------------------------------------------------
!     Collapse Points
!
! hash_liste(i) contains the last index "n" of "hash" 
!               with "hash (n) = i"
! next (n)      contains the next index "k" with "hash(k) = i"
!
! liste(1:n)    is a compact list of all indices of "hash"
!               which contain a significant point (not duplicated)
! locs2liste    contains the transformation
!               from an index in "locs" (or "hash")
!               into an index in "liste".
!------------------------------------------------------------------------
!
      Allocate (hash_liste (nd_hliste), stat = ierror)
      if (ierror /= 0) then
         ierrp (1) = nd_hliste

         ierror = PRISM_Error_Alloc

         call psmile_error ( ierror, 'nd_hliste', ierrp, 1, &
                 __FILE__, __LINE__ )
         return
      endif
!
      Allocate (next (nlocs), liste (nlocs), locs2liste (nlocs), &
                stat = ierror)
      if (ierror /= 0) then
         ierrp (1) = nlocs * 3

         ierror = PRISM_Error_Alloc

         call psmile_error ( ierror, 'next, liste, locs2liste', ierrp, 1, &
                 __FILE__, __LINE__ )
         return
      endif
!
      next (:) = 0
      hash_liste (:) = 0
!
#ifdef PRISM_ASSERTION
      locs2liste (:) = -11
#endif
!
!===> ... Map locations in locs to hash range
!         hash_liste (i) = Last index "n" in "locs" with
!         "mod (hash (n), nd_hliste) + 1 = i"
!
!cdir vector
         do n = 1, nlocs
         if (hash (n) > 0) then
            i = mod (hash(n), nd_hliste) + 1
            if (hash_liste(i) > 0) then
               next (n) = hash_liste(i)
            endif

            hash_liste(i) = n
         endif
         end do
!
#ifdef DEBUGX
         print *, 'length', length (:)
         print *, 'hash', hash (1:nlocs)
         print *, 'hashliste', hash_liste (1:nd_hliste)
#endif
!
!===> ... Create compact list "liste" of indices (in "locs")
!         of different "hash" values.
!         Create mapping "locs2liste" from "locs" to "liste"
!
! n_liste = Number of points to be sent to requesting process
!
      n = 0
!cdir vector
         do i = 1, nd_hliste
         if (hash_liste (i) > 0) then
            n = n + 1
!
!===> ... Store first index for hash value "hash_liste(i)"
!
            liste (n) = hash_liste (i)
            locs2liste (hash_liste (i)) = n

            nbeg = n
            nend = n
!
            k = hash_liste (i)
              do while (next (k) > 0) 
!
!===> ... Was the location "hash(k)" already stored in "liste" ?
!         If not, add index of location to "liste".
!         Otherwise, store mapping from "locs" to index in "liste"
!
              k = next (k)
                 do l = nbeg, nend
                 if (hash (k) == hash (liste(l))) exit
                 end do

              if (l > nend) then
                 n = n + 1
!
                 liste (n) = k
                 locs2liste (k) = n
                 nend = nend + 1
              else
                 locs2liste(k) = l
              endif
              end do
!
         endif
         end do
!
      n_liste = n
!
#ifdef PRISM_ASSERTION
!
!===> Check assertions
!
!cdir vector
         do n = 1, nlocs
         if (hash (n) > 0 .and. &
             (locs2liste (n) <= 0 .or. locs2liste(n) > n_liste)) exit
         end do
!
      if (n < nlocs) then
         print *, 'n, nlocs, locs2liste (n)', n, nlocs, locs2liste (n)
         call psmile_assert ( __FILE__, __LINE__, &
                             'locs2liste is not totally filled')
      endif
#endif
!
!===> Remove arrays which are not further required
!
      Deallocate (hash_liste, next)
!
!===> Store data on locations found
!
      if (n_liste > 0) then
#if 0
!
!===> ... Get number of points found
!
         n = 0
         n_found = 0
            do j = 1, n_corners
!cdir vector
               do i = 1, n_send
               if (ijkmsk (i,j)) then
                  n = n + 1
                  if (hash(n) >= 0) n_found = n_found + 1
               endif 
               end do
            end do
!
#ifdef PRISM_ASSERTION
         if (n /= nlocs) then
            print *, 'nlocs, n', nlocs, n
            call psmile_assert ( __FILE__, __LINE__, &
                                'n /= nlocs')
         endif
!
         if (n_found > nlocs) then
            print *, 'nlocs, n_found', nlocs, n_found
            call psmile_assert ( __FILE__, __LINE__, &
                                'n_found > nlocs')
         endif
#endif
#else
         n_found = Count (hash(:) >= 0)
!cdir vector
!        n_found = 0
!           do n = 1, nlocs
!           if (hash(n) >= 0) n_found = n_found + 1
!           end do
#endif
!
!===> ... Allocate vectors required for other routines ...
!
         Allocate (neigh (ndim_3d, n_liste), index_found (n_found), &
                   stat = ierror)
         if (ierror /= 0) then
            ierrp (1) = n_liste * ndim_3d + n_found

            ierror = PRISM_Error_Alloc

            call psmile_error ( ierror, 'n_liste', ierrp, 1, &
                    __FILE__, __LINE__ )
            return
         endif
!
!===> ... and store indices of points to be transferred to
!         the destination process
!
!cdir vector
            do n = 1, n_liste
            neigh (1, n) = locs (1, liste(n))
            neigh (2, n) = locs (2, liste(n))
            neigh (3, n) = locs (3, liste(n))
            end do
!
#ifdef DEBUGX
        print *, 'locations to be sent: nd_liste', n_liste
            do n = 1, n_liste
            print *, 'neigh', neigh (:, n)
            end do
#endif
!
!===> ... 
!         Store indices to "neigh" in corresponding sequence
!         of locations found.
!
         nf = 0
!cdir vector
            do n = 1, nlocs
            if (hash (n) >= 0) then
!
!===> ...... Required interpolation base was found
!            (may be masked out (hash(n) == masked_out).
!            Corresponding value is stored in position "index_found(n)".
!
               nf = nf + 1
!
               if (hash(n) > 0) then
                  index_found (nf) = locs2liste (n)
               else
                  index_found (nf) = masked_out
               endif
            endif
!
            end do ! n

#ifdef PRISM_ASSERTION
        if (nf /= n_found) then
           print *, 'nf, n_found, index', nf, n_found
           call psmile_assert ( __FILE__, __LINE__, &
                   'nf /= n_found')
        endif
!
           do i = 1, n_found
           if (index_found(i) /= masked_out .and. &
              (index_found(i) < 1 .or. index_found(i) > n_liste)) exit
           end do
!
        if (i < n_found) then
           print *, 'i, n_found, index', i, n_found, index_found (i)
           call psmile_assert ( __FILE__, __LINE__, &
                   'invalid found index generated')
        endif
#endif
!
#ifdef DEBUGX
        print *, 'index_found', n_found
            do i = 1, n_found
            print *, 'i', i, index_found (i)
            end do
#endif 

      else
!
!===> No point found
!
         NULLIFY (neigh)
         NULLIFY (index_found)

         n_found = 0

      endif
!
!===> Deallocate memory allocated
!
      Deallocate (locs2liste)
!
!===> Set return info
!
      search%n_liste  = n_liste
      search%neigh_3d => neigh
!
      search%n_found  = n_found
      search%index_found => index_found
!
!===> All done
!
#ifdef VERBOSE
      print 9980, trim(ch_id), ierror, n_found, nlocs

      call psmile_flushstd
#endif /* VERBOSE */
!
      return
!
!  Formats:
!
#ifdef VERBOSE
9990 format (1x, a, ': psmile_hash_extra')
9980 format (1x, a, ': psmile_hash_extra: eof', &
                    ', ierror =', i4, ', n_found =', i8, i8)
#endif /* VERBOSE */

      end subroutine PSMILe_Hash_extra