psmile_trili_3d_extra_off.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_Trili_3d_extra_off
!
! !INTERFACE:

      subroutine psmile_trili_3d_extra_off (comp_info, search,         &
                      grid_id, mask_array, mask_shape, mask_available, &
                      ibuf, len_item, n_send, num_neigh,               &
                      ierror)
!
! !USES:
!
      use PRISM
!
      use PSMILe, dummy_interface => PSMILe_Trili_3d_extra_off

      implicit none
!
! !INPUT PARAMETERS:
!
      Type (Enddef_comp),           Intent (In)    :: comp_info
!
!     Info on the component in which the donor cells
!     should be searched.

      Integer,                      Intent (In)    :: grid_id

!     Grid Id of the source grid (grid on which the points are 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)     :: len_item
!
!     Number of pieces of information in integer buffer "ibuf".
!        1:3: (i,j,k) index of lower left corner
!        4  : Index in
!        5  : Code for points required
!
      Integer,                     Intent (In)    :: n_send
!
!     Number of extra points sent

      Integer,                     Intent (In)    :: num_neigh
!
!     (Maximal) Number of interpolations bases required.
!     Need to be the first "num_neigh" bases in "ijkstd" (see below).
!
!
! !INPUT/OUTPUT PARAMETERS:
!
      Type (Enddef_global_search), Intent (InOut) :: search

!     Data on the points to be searched

      Integer,                     Intent (InOut) :: ibuf (len_item, n_send)
!
!     Integer buffer containing the packed integer data.
!     On return, ibuf (5, *) contains the code of the points found.
!
! !OUTPUT PARAMETERS:
!
      Integer,                     Intent (Out)   :: ierror

!     Returns the error code of PSMILe_Search_donor_extra_off;
!             ierror = 0 : No error
!             ierror > 0 : Severe error
!
! !LOCAL PARAMETERS:
!
! n_corners_3d = Number of interpolation bases in 3d case
!
! masked_out   = Hash value for a point which is masked out
!
      Integer, Parameter           :: n_corners_3d = 2 * 2 * 2
!
      Integer, Parameter           :: masked_out = 0
!
! !LOCAL VARIABLES
!
!     ... loop variables
!
      Integer                      :: i, j, n
!
!     ... 
!
      Integer                      :: index0
      Integer                      :: length (ndim_3d)
!
!     ... For Grids
!
      Type (Grid), Pointer         :: grid_info
      Integer,     Pointer         :: grid_valid_shape (:, :)
      Logical,     Pointer         :: cyclic (:)
!
!     ... For searched points of control volume
!
      Integer                      :: code, n_corners, nlocs
!
      Integer                      :: ijkstd (ndim_3d, n_corners_3d)
      Integer                      :: ijkctl (ndim_3d, n_corners_3d)
      Integer, Allocatable         :: locs (:, :)
      Integer, Allocatable         :: hash (:)
      Logical, Allocatable         :: ijkmsk (:, :)
!
!     ... for error parameters
!
      Integer, parameter           :: nerrp = 1
      Integer                      :: ierrp (nerrp)
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_Trili_3d_extra_off" performs the additional (global)
! search for
!
! TODO: PSMILe_Trili_3d_extra_off und PSMILe_Tricu_3d_extra_off vereinigen
!       locs alloziieren und ausrechnen und danach gleiche Routine
!       verwenden
!       
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.07.05    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
! $Id: psmile_trili_3d_extra_off.F90 2082 2009-10-21 13:31:19Z hanke $
! $Author: hanke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_trili_3d_extra_off.F90 2082 2009-10-21 13:31:19Z hanke $'
!
!----------------------------------------------------------------------
!
!  Relative indices of corner
!
#ifdef FORTRAN_ORDER
      data ((ijkstd (i, n), i=1,ndim_3d), n = 1, n_corners_3d) &
                            /  0, 0, 0,   1, 0, 0, &
                               0, 1, 0,   1, 1, 0, &
                               0, 0, 1,   1, 0, 1, &
                               0, 1, 1,   1, 1, 1/

#else
!huhu Offenbar erwartet Damien's Routine eine andere Ordnung
      data ((ijkstd (i, n), i=1,ndim_3d), n = 1, n_corners_3d) &
                            /  0, 0, 0,   1, 0, 0, &
                               1, 1, 0,   0, 1, 0, &
                               0, 0, 1,   1, 0, 1, &
                               1, 1, 1,   0, 1, 1/
#endif
!
!----------------------------------------------------------------------
!
!  Initialization
!
#ifdef VERBOSE
      print 9990, trim(ch_id), comp_info%comp_id

      call psmile_flushstd
#endif /* VERBOSE */
!
      ierror = 0
      n_corners = num_neigh
!
!===> Get grid info
!
!  grid_valid_shape (2, ndim_3d) : Specifies the valid block shape
!                                  for the "ndim_3d"-dimensional block
!                                  (without halo/overlap regions)
!  cyclic (ndim_3d) : Does the block have cyclic coordinates
!                     in the corresponding direction ?
!
      grid_info        => Grids(grid_id)
!
      grid_valid_shape => grid_info%grid_shape
      cyclic           => grid_info%cyclic
!
      length (1:ndim_3d) = grid_valid_shape(2,1:ndim_3d) - &
                           grid_valid_shape(1,1:ndim_3d) + 1
!
#ifdef PRISM_ASSERTION
      if (grid_info%grid_type /= PRISM_Reglonlatvrt     .and. &
          grid_info%grid_type /= PRISM_Irrlonlat_Regvrt .and. &
          grid_info%grid_type /= PRISM_Irrlonlatvrt) then
!
         print *, "grid_type", grid_info%grid_type
         call psmile_assert (__FILE__, __LINE__, &
                   "This routine is not designed for this grid type")
      endif
!
      if (num_neigh < 1 .or. num_neigh > n_corners_3d) then
         print *, 'num_neigh', num_neigh, n_corners_3d
         call psmile_assert ( __FILE__, __LINE__, &
                'Invalid value for num_neigh (number of interpolation bases')
      endif
#endif
!
!===> Allocate mask "ijkmsk" for points to be searched.
!     This corresponds to the decoded code "code".
!
      Allocate (ijkmsk (n_send, n_corners), stat = ierror)
      if (ierror /= 0) then
         ierrp (1) = n_send * n_corners

         ierror = PRISM_Error_Alloc

         call psmile_error ( ierror, 'ijkmsk', ierrp, 1, &
                 __FILE__, __LINE__ )
         return
      endif
!
!===> Decode code for searched points and
!     create list of locations "locs" of points required
!
      code = 1
!
         do j = 1, n_corners
!cdir vector
            do i = 1, n_send
            ijkmsk (i,j) = IAND (code, ibuf (5, i)) /= 0
            end do

         code = code * 2
         end do
!
      nlocs = Count (ijkmsk(:,:))
!
#ifdef PRISM_ASSERTION
      if (nlocs <= 0) then
         call psmile_assert ( __FILE__, __LINE__, &
                             'Number of points to be searched == 0 ?!?')
      endif
#endif
!
#ifdef DEBUGX
      print *, 'code for locations searched', n_send
         do i = 1, n_send
         print *, 'code', ibuf (5, i), ijkmsk(i, :)
         end do
#endif
!
      Allocate (locs (ndim_3d, nlocs), hash (nlocs), stat = ierror)
      if (ierror /= 0) then
         ierrp (1) = (ndim_3d+1) * nlocs

         ierror = PRISM_Error_Alloc

         call psmile_error ( ierror, 'locs, hash', ierrp, 1, &
                 __FILE__, __LINE__ )
         return
      endif
!
!===> For special cases (2d Hyperplanes)
!     ??? Braucht man hier nicht ein spezielles Flag,
!         dass das globale Gebiet eine Hyperplane ist ?
!
      ijkctl = ijkstd
!
         do j = 1, ndim_3d
         if (length(j) == 1) ijkctl (j, :) = 0
         end do
!
      n = 0
         do j = 1, n_corners
!cdir vector
            do i = 1, n_send
            if (ijkmsk (i,j)) then
               n = n + 1
               locs (1, n) = ibuf (1, i) + ijkctl (1, j)
               locs (2, n) = ibuf (2, i) + ijkctl (2, j)
               locs (3, n) = ibuf (3, i) + ijkctl (3, j)
            end if
            end do ! i
         end do ! j
!
#ifdef PRISM_ASSERTION
      if (n /= nlocs) then
         print *, 'n, nlocs', n, nlocs
         call psmile_assert ( __FILE__, __LINE__, &
                'Inconsistent values of n and nlocs; must be same value')
      endif
#endif
!
!===> Set indices for cyclic coordinate directions
!
         do i = 1, ndim_3d
         if (cyclic(i) .and. length(i) > 1) then
            index0 = grid_valid_shape(1,i)
!
!cdir vector
               do n = 1, nlocs
               locs (i, n) = index0 + mod (locs(i, n) - index0, length(i))
               end do
         endif
         end do
!
!===> Create compact list of points to be sent and 
!     return hash value "hash" (index in fortran order) for each location
!
      call psmile_hash_extra (search, locs, hash, nlocs,                &
              mask_array, mask_shape, mask_available, grid_valid_shape, &
              ierror)
      if (ierror > 0) return
!
#ifdef DEBUGX
      print *, 'hash values and locs ', grid_valid_shape
      n = 0
         do j = 1, n_corners
!cdir vector
            do i = 1, n_send
            if (ijkmsk (i,j)) then
               n = n + 1
               print *, 'i, j, n, hash(n), locs(:,n) ', &
                         i, j, n, hash(n), locs(:,n)
            end if
            end do ! i
         end do ! j
#endif
!
!===> Save coded mask of points found in "ibuf"
!     ??? Ob das so clever ist ?
!
      if (search%n_liste > 0) then
         ibuf (5, 1:n_send) = 0
!
         n = 0
         code = 1
            do j = 1, n_corners
!cdir vector
               do i = 1, n_send
               if (ijkmsk (i,j)) then

                  n = n + 1
!
                  if (hash (n) >= 0) then
!
!===> ...... Required "j-th" interpolation base was found
!            (may be masked out (hash(n) == masked_out == 0).
!            Corresponding value is stored in position "index_found(n)".
!
                     ibuf (5, i) = ibuf (5, i) + code
                  endif
!
               endif
               end do ! i

            code = code * 2
            end do ! j

#ifdef PRISM_ASSERTION
         if (n /= nlocs) then
            print *, 'nlocs, n', nlocs, n
            call psmile_assert ( __FILE__, __LINE__, &
                                'n /= nlocs')
         endif
#endif
!
#ifdef DEBUGX
        print *, 'code for locations found', n_send
            do i = 1, n_send
            print *, 'code', i, ibuf (5, i), ijkmsk(i, :)
            end do
#endif 

      endif
!
!===> Deallocate memory allocated
!
      Deallocate (locs, hash, ijkmsk)
!
!===> All done
!
#ifdef VERBOSE
      print 9980, trim(ch_id), ierror, search%n_found, n_send

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

      end subroutine PSMILe_Trili_3d_extra_off