psmile_mg_prev_levels_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_Prev_levels_2d_dble
!
! !INTERFACE:

      subroutine psmile_mg_prev_levels_2d_dble (grid_id, lev, nlev, &
                      lstijk, xyz, found, newijk, nc)
!
! !USES:
!
      use PRISM_constants
!
      use PSMILe, dummy_interface => PSMILe_mg_prev_levels_2d_dble

      Implicit none
!
! !INPUT PARAMETERS:
!
      Integer, Intent (In)            :: grid_id
!
!     Grid Id for which the donor cells should be searched.
!
      Integer, Intent (In)            :: lev
!
!     Level number of the next (fine) grid
!
      Integer, Intent (In)            :: nlev
!
!     Total number of levels
!
      Integer, Intent (In)            :: nc
!
!     Number of points to be searched.
!
      Integer, Intent (In)            :: lstijk (ndim_2d, nc)
!
!     Last indices of cell on coarse level LEV+1.
!     !!! This routine assumes that all parts on the fine grid "LEV"
!     of these cells, were already controlled and point "xyz"
!     was not found.
!
      Double Precision, Intent (In)   :: xyz (ndim_2d, nc)
!
!     Coordinates of the points to be searched
!
! !OUTPUT PARAMETERS:
!
      Integer, Intent (Out)           :: found (nc)
!
!     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, nc)
!
!     Indices on level LEV if found > 0 is returned.
!
! !DEFINED PARAMETERS:
!
      Logical, Parameter              :: wide = .true.
      Logical, Parameter              :: small = .false.
!
! !LOCAL VARIABLES
!
      Integer                         :: i, n

!     ... for Grid

      Type (Grid), Pointer            :: grid_info
      Type (Enddef_mg), Pointer       :: mg_infos (:)

!     ... for levels

      Integer                         :: levc, levd, levdbg, levu, levubg
      Integer                         :: nlev1
      Integer                         :: nold (nc, lev:nlev-1)
      Integer                         :: ijk (ndim_2d, nc, lev:nlev-1)
      Integer                         :: ignore (ndim_2d, lev:nlev-1)
      Integer                         :: coarsf (ndim_2d)
!
#ifdef DEBUGX
      logical                         :: control
      common /huhu1/ control
#endif
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_mg_prev_levels_2d_dble" searches for the donor cell
! of points "xyz (*, 1:nc)". The points was last found on coarse level
! "lev+1" in donor cell "lstijk (*, 1:nc)". 
! If the i-th point is found in another cell of the coarser grid,
! "found(i)" > 0 is returned and
! "newijk (*, i)" returns the donor cell index on level "lev".
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.07.05    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
!  $Id: psmile_mg_prev_levels_2d_dble.F90 2325 2010-04-21 15:00:07Z valcke $
!  $Author: valcke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_mg_prev_levels_2d_dble.F90 2325 2010-04-21 15:00:07Z valcke $'
!
!----------------------------------------------------------------------
!
!  Initialization
!
#ifdef VERBOSE
      print 9990, trim(ch_id), lev, nc

      call psmile_flushstd
#endif /* VERBOSE */
!
      grid_info => Grids(grid_id)
      mg_infos => grid_info%mg_infos

      levc = lev + 1
      nlev1 = nlev - 1

#ifdef PRISM_ASSERTION
      if (levc >= nlev) then
         call psmile_assert (__FILE__, __LINE__, 'incorrect level')
      endif
#endif
!
!     make sure that ijk lies on coarse grid line
!     TODO: Store coarsening factor in grid info
!
         do i = 1, ndim_2d
         if (mg_infos(lev )%levdim(i) /= &
             mg_infos(levc)%levdim(i) ) then
            ijk (i, :, levc) = (lstijk (i, :) / 2) * 2
            coarsf(i) = 2
         else
            ijk (i, :, levc) =  lstijk (i, :)
            coarsf(i) = 1
         endif
         end do
!
         do levd = levc, nlev1
         nold (:, levd) = 0 ! not clear to which part "lstijk" belongs
         end do
!
!===> ... Compute indices on the coarser levels
!         for coarsening factor 2 !
!
         do levd = levc+1, nlev1
            do i = 1, ndim_2d
            if (mg_infos(levd  )%levdim(i) /= &
                mg_infos(levd-1)%levdim(i) ) then
!
!              make sure that ijk lies on coarse grid line
!
               ijk (i, :, levd) = (ijk (i, :, levd-1) / 4) * 2
            else
               ijk (i, :, levd) = ijk (i, :, levd-1)
            endif
            end do
         end do
!
!---------------------------------------------------------------------
!     Go down in the level hierarchie
!
!     TODO: transfer nc loop into PSMILe_mg_control_cell_2d_dble
!           Use mask for points to be controlled and/or
!           gather/scatter of points to be controlled.
!
!     Note: The "wide" control may cause some additional work
!           (repeated control of fine grid cells) 
!           at 2**i boundaries.
!           Otherwise, the "wide" control may find the correct
!           cell faster instead of looping till to the coarsest level.
!           It's not totally clear, whether "small" control would
!           be better should.
!           Should performance control done with example
!           "simple-orca-atmland" since it generates many controls.
!           Possibly, adaptive setting of "wide" or "small" control.
!---------------------------------------------------------------------
!
loop_n: do n = 1, nc

      levdbg = levc

!     Define cells to be ignored

         do i = 1, ndim_2d
         if (mg_infos(lev )%levdim(i) > &
             mg_infos(levc)%levdim(i)) then
            ignore (i, lev)  = lstijk (i, n) * 2
         else
            ignore (i, lev)  = lstijk (i, n)
         endif
         enddo
!
      ignore (:, levc) = lstijk (:, n)
!
         do i = 1, ndim_2d
            do levd = levc+1, nlev1
            if (mg_infos(levd-1)%levdim(i) > &
                mg_infos(levd)%levdim(i)) then
               ignore (i, levd) = ignore (i, levd-1) / 2
            else
               ignore (i, levd) = ignore (i, levd-1)
            endif
            end do
         end do
!
#ifdef DEBUGX
      control = lev == 3 .and. n == 26
!
      if (control) then
         do levd = lev, nlev1
         print *, 'ignore lev', levd, 'ignore ', ignore (:, levd)
         enddo
      endif 
#endif
!
100   continue
!
         do levd = levdbg, nlev1
#ifdef DEBUGX
         if (control) then
            print *, 'levd', levd, ', ignore', ignore (:, levd), 'lstijk', lstijk (:, n)
         endif
#endif
!
         call psmile_mg_control_cell_2d_dble (                 &
             mg_infos(levd)%double_arrays%chmin(1)%vector,     &
             mg_infos(levd)%double_arrays%chmin(2)%vector,     &
             mg_infos(levd)%double_arrays%chmax(1)%vector,     &
             mg_infos(levd)%double_arrays%chmax(2)%vector,     &
             mg_infos(levd)%double_arrays%midp (1)%vector,     &
             mg_infos(levd)%double_arrays%midp (2)%vector,     &
             mg_infos(levd)%levdim,                            &
             ijk(1, n, levd), xyz(:,n), nold(n, levd),         &
             ignore (1, levd), wide, found (n), newijk (:, n))

#ifdef DEBUGX
             if (control) then
                print *, 'levd', levd, ', found', found(n), ', newijk', newijk (:, n)
             endif
#endif
!
         if (found (n) > 0) then
            go to 1990
         endif
!
         end do
!
!===> Not found
!
      cycle loop_n
!
!---------------------------------------------------------------------
!     Up again; startpoint is newijk on level "levd"
!---------------------------------------------------------------------
!
1990  nold (n, levd) = found (n)
!
      levdbg = levd
      levubg = levd - 1
!
      nold (n, levubg) = 0
!
         do i = 1, ndim_2d
         if (mg_infos(levubg  )%levdim(i) > &
             mg_infos(levubg+1)%levdim(i)) then
            ijk (i, n, levubg) = min (newijk(i,n)*2, &
                                      mg_infos(levubg)%levdim(i))
         else
            ijk (i, n, levubg) = newijk (i,n)
         endif
         enddo
!
1995  continue
!
         do levu = levubg, lev, -1
!
2000     continue
!
#ifdef DEBUGX
         if (control) then
            print *, 'levu', levu, ', ignore', ignore (:, levu)
         endif
#endif
!
         call psmile_mg_control_cell_2d_dble (                &
             mg_infos(levu)%double_arrays%chmin(1)%vector,    &
             mg_infos(levu)%double_arrays%chmin(2)%vector,    &
             mg_infos(levu)%double_arrays%chmax(1)%vector,    &
             mg_infos(levu)%double_arrays%chmax(2)%vector,    &
             mg_infos(levu)%double_arrays%midp (1)%vector,    &
             mg_infos(levu)%double_arrays%midp (2)%vector,    &
             mg_infos(levu)%levdim,                           &
             ijk(1, n, levu), xyz(:,n), nold(n, levu),        &
             ignore (1, levu), small, found(n), newijk (:,n))
!
#ifdef DEBUGX
         if (control) then
            print *, 'levu', levd, ', found', found(n)
         endif
#endif
!
         if (found(n) == 0) then
            if (levu < levd-1) then
!
!===> ... go back to the previous level; and try again
!
               levubg = levu + 1
               go to 1995
            else
!
!===> ... go down
!
               go to 100
            endif
!
         else if (levu > lev) then
!           Not the finest level
            nold (n, levu) = found (n)
!
!           Initialize search on next finer level "levu-1"
!
            nold (n, levu-1) = 0
!
               do i = 1, ndim_2d
               if (mg_infos(levu-1)%levdim(i) > &
                   mg_infos(levu  )%levdim(i)) then
                  ijk (i, n, levu-1) = min (newijk(i,n)*2, &
                                       mg_infos(levu-1)%levdim(i))
               else
                  ijk (i, n, levu-1) = newijk (i,n)
               endif
               enddo
!
         else ! levu == lev
!           Finest level
            if (lstijk(1, n) == newijk(1,n)/coarsf(1) .and. &
                lstijk(2, n) == newijk(2,n)/coarsf(2)) then
!                 Initial cell lstijk was found.
!                 Search again on same level.
!
               nold (n, levu) = found (n)
!
               go to 2000
            endif
         endif
         end do ! levu
!
         end do loop_n
!
!===> All done
!
#ifdef VERBOSE
      print 9980, trim(ch_id), lev

      call psmile_flushstd
#endif /* VERBOSE */
!
      return
!
!  Formats:
!
9990 format (1x, a, ': psmile_mg_prev_levels_2d_dble: level', i3, &
                    ', nc', i6)
9980 format (1x, a, ': psmile_mg_prev_levels_2d_dble: eof, level', i3)

      end subroutine PSMILe_mg_prev_levels_2d_dble