psmile_mg_control_cell_2d_dble.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_MG_Control_cell_2d_dble
!
! !INTERFACE:

      subroutine psmile_mg_control_cell_2d_dble (          &
                      chmin1, chmin2,                      &
                      chmax1, chmax2,                      &
                      midp1,  midp2,           levdim,     &
                      ijk, xyz, nold, ignore, wide, found, newijk)
!
! !USES:
!
      use PRISM_constants
!
      use PSMILe, dummy_interface => PSMILe_mg_control_cell_2d_dble

      implicit none
!
! !INPUT PARAMETERS:
!
      Integer, Intent (In)                        :: levdim (ndim_2d)

!     Dimensions of chmin, chmax and midp

      Double Precision, Intent (In)               :: chmin1 (0:levdim(1), 
                                                             0:levdim(2))
      Double Precision, Intent (In)               :: chmin2 (0:levdim(1), 
                                                             0:levdim(2))

!     Minimum of grid coordinates per grid cell

      Double Precision, Intent (In)               :: chmax1 (0:levdim(1), 
                                                             0:levdim(2))
      Double Precision, Intent (In)               :: chmax2 (0:levdim(1), 
                                                             0:levdim(2))

!     Maximum of grid coordinates per grid cell

      Double Precision, Intent (In)               ::  midp1 (0:levdim(1), 
                                                             0:levdim(2))
      Double Precision, Intent (In)               ::  midp2 (0:levdim(1), 
                                                             0:levdim(2))

!     Mid point of the cell

      Integer, Intent (In)                        :: ijk (ndim_2d)

!     Index of the cell for which the cell parts should be
!     controlled.

      Double Precision, Intent (In)               :: xyz (ndim_2d)

!     Corrdinates of the point for which is searched

      Integer, Intent (In)                        :: nold

!     Old number of the cell part where the point XYZ was found last.
!     And all other parts "i" with 
!                DIST(i) < DIST(NOLD) were already controlled.
!     NOLD = 0: Start first search in this cell.
!
      Integer, Intent (In)                        :: ignore (ndim_2d)
!
!     Cell indices to be ignored.
!
      logical, Intent (In)                        :: wide

!     Wide control ?
!
! !OUTPUT PARAMETERS:
!
      integer, Intent (Out)                       :: found

!     Was the point xyz found ?
!     found = 0: Not found
!     found > 0: Part of the cell where point XYZ was found.

      integer, Intent (Out)                       :: newijk (ndim_2d)

!     Indices on level LEV if found > 0 is returned.
!
! !DEFINED PARAMETERS:
!
      Double Precision, Parameter                 :: remax = 1.0e20
!
      Integer, Parameter                          :: ndtrys = 4
      Integer, Parameter                          :: ndtryw = 16
!
! !LOCAL VARIABLES
!
      Double Precision                            :: dist (ndtryw)
      Integer                                     :: i, n, ndtry
      Integer                                     :: nmin (1)
!
      Integer, Pointer                            :: ijkctl (:, :)
      Integer, Target                             :: ijkctls (ndim_2d, ndtrys)
      Integer, Target                             :: ijkctlw (ndim_2d, ndtryw)
      Integer                                     :: ijkdst (ndim_2d, ndtryw)
      Logical                                     :: within (ndtryw)
!
#ifdef DEBUGX
      logical control
      common /huhu1/ control
#endif
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_mg_control_cell_2d_dble" searches the point XYZ
! in the n the neihbourhood of a cell. The cell which is controlled
! starts in IJK. All cells with (i) in
!      (ijk(:)-1:ijk(:)+2) if wide is true.
!      (ijk(:)-1:ijk(:)+1) otherwise
! are controlled.
!
! ?!? Vielleicht koennte man alle Punkte zurueckgeben, die man
!     gefunden hat.
!
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.07.21    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
!  $Id: psmile_mg_control_cell_2d_dble.F90 2325 2010-04-21 15:00:07Z valcke $
!  $Author: valcke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_mg_control_cell_2d_dble.F90 2325 2010-04-21 15:00:07Z valcke $'
!
!----------------------------------------------------------------------
!
!  Relative control indices for sub-cells of coarse grid cells
!
      data ((ijkctls (i, n), i=1,ndim_2d), n = 1,ndtrys) &
                             / 0, 0,   1, 0, &
                               0, 1,   1, 1/
!
!  Relative control indices for sub-cells of coarse grid cells
!
!  27 (-1:1, -1:1, -1:1) Zellen sollte auch reichen
!  Die reihenfolge der ersten 8 muss so bleiben !
!
      data ((ijkctlw (i, n), i=1,ndim_2d), n = 1,ndtryw) &
                             / 0, 0,   1, 0,   0, 1,   1, 1, &
                              -1,-1,   0,-1,   1,-1,   2,-1, &
                              -1, 0,                   2, 0, &
                              -1, 1,                   2, 1, &
                              -1, 2,   0, 2,   1, 2,   2, 2/
!
!----------------------------------------------------------------------
!
!  Initialization
!
#ifdef VERBOSE_extra
      print 9990, trim(ch_id)

      call psmile_flushstd
#endif /* VERBOSE */
!
!===> Compute indices to be controlled
!
      if (wide) then
         ndtry = ndtryw
         ijkctl => ijkctlw
      else
         ndtry = ndtrys
         ijkctl => ijkctls
      endif
!
#ifdef PRISM_ASSERTION
      if (nold < 0 .or. nold > ndtry) then
         call psmile_assert (__FILE__, __LINE__, &
                             'wrong nold')
      endif
!
      if (minval(ijk(:)) < 0 .or.  minval(levdim(:)-ijk(:)) < 0) then
         call psmile_assert (__FILE__, __LINE__, 'wrong ijk')
      endif
#endif
!
!===> Compute indices to be controlled
!
        do n = 1, ndtry
        ijkdst (:, n) = ijk(:) + ijkctl(:,n)
        end do
!
!===> Control index range
!
        do n = 1, ndtry
        within (n) = ijkdst (1,n) >= 0 .and. ijkdst(1,n) <= levdim(1) .and. &
                     ijkdst (2,n) >= 0 .and. ijkdst(2,n) <= levdim(2)
        end do
#ifdef DEBUGX
if (control) print *, 'within a:', within (1:ndtry)
#endif
!
!===> Control location of point XYZ in the boxes
!
!  ??? ist (sent-chmin)*(sent-chmax) <= 0.0 schneller ???
!
!cdir vector
        do n = 1, ndtry
        if (within(n)) then
           within (n) = xyz(1) >= chmin1(ijkdst(1,n), ijkdst(2,n)) .and. &
                        xyz(2) >= chmin2(ijkdst(1,n), ijkdst(2,n)) .and. &
                        xyz(1) <= chmax1(ijkdst(1,n), ijkdst(2,n)) .and. &
                        xyz(2) <= chmax2(ijkdst(1,n), ijkdst(2,n))
        endif
        end do
#ifdef DEBUGX
if (control) print *, 'within b:', within (1:ndtry)
#endif
!
!===> Compute distances to the cell midpoints
!
!cdir vector
        do n = 1, ndtry
        if (within(n)) then
            dist (n) = (xyz(1) - midp1(ijkdst(1,n), ijkdst(2,n)))**2 &
                     + (xyz(2) - midp2(ijkdst(1,n), ijkdst(2,n)))**2

        else

            dist (n) = remax

        endif
        end do
!
!===> Set old values to remax
!
      if (nold > 0) then

#ifdef PRISM_ASSERTIONX
         if (dist(nold) .eq. remax) then
            write (*, 9970) xyz
9970  format (1x, 1p, 3e18.9)
!
            call psmile_assert (__FILE__, __LINE__, &
                                'incorrect nold')
         endif
#endif
!
!===> Exclude all old distances
!
            do n = 1, nold-1
            if (dist(n) <= dist(nold)) dist (n) = remax
            end do
!
            do n = nold+1, ndtry
            if (dist(n) < dist(nold)) dist (n) = remax
            end do
!
         dist (nold) = remax
      endif
!
!===> Look for the minimum distance
!
100   continue
!
      nmin = MINLOC (dist(1:ndtry))
#ifdef MINLOCFIX
      if (nmin(1).eq.0) nmin=1
#endif /* MINLOCFIX */
!
!===> Was a valid cell found ?
!
      if (dist(nmin(1)) == remax) then
         found = 0
      else
!
         if (ijkdst(1, nmin(1)) == ignore (1) .and. &
             ijkdst(2, nmin(1)) == ignore (2)) then
!           Cell has to be  ignored

            dist (nmin(1)) = remax
            go to 100
         endif
!
         found = nmin (1)
!
         newijk (:) = ijkdst (:, nmin(1))

#ifdef PRISM_ASSERTION
         if ( xyz(1) < chmin1(newijk(1), newijk(2)) .or. &
              xyz(1) > chmax1(newijk(1), newijk(2)) .or. &
              xyz(2) < chmin2(newijk(1), newijk(2)) .or. &
              xyz(2) > chmax2(newijk(1), newijk(2))) then
            print *, 'xyz', xyz
            print *, 'chmin', chmin1(newijk(1), newijk(2)), &
                              chmin2(newijk(1), newijk(2))
            print *, 'chmax', chmax1(newijk(1), newijk(2)), &
                              chmax2(newijk(1), newijk(2))

            call psmile_assert (__FILE__, __LINE__, &
                      'not contained in computed cell')
        endif
#endif
      endif
!
!===> All done
!
#ifdef VERBOSE_extra
      print 9980, trim(ch_id)

      call psmile_flushstd
#endif /* VERBOSE */
!
      return
!
!  Formats:
!
9990 format (1x, a, ': PSMILe_mg_control_cell_2d_dble:')
9980 format (1x, a, ': PSMILe_mg_control_cell_2d_dble: eof')

      end subroutine PSMILe_mg_control_cell_2d_dble