psmile_mg_first_subgrid_3d_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_First_subgrid_3d_dble
!
! !INTERFACE:

      subroutine psmile_mg_first_subgrid_3d_dble (array, &
                    idlow, idhigh, jdlow, jdhigh, kdlow, kdhigh, &
                    nbr_corners,   &
                           range,                        &
                    chmin, chmax, midp,                  &
                    levdim1, levdim2, levdim3, ierror)
!
! !USES:
!
      use PRISM_constants
!
      use PSMILe, dummy_interface => PSMILe_MG_First_subgrid_3d_dble

      implicit none
!
! !INPUT PARAMETERS:
!
      Integer, Intent (In)             :: idlow, idhigh, jdlow, jdhigh
      Integer, Intent (In)             :: kdlow, kdhigh, nbr_corners

!     Dimensions of "array"

      Double Precision, Intent (In)    :: array (idlow:idhigh, jdlow:jdhigh, 
                                                 kdlow:kdhigh, nbr_corners)

!     Fully dimensioned arrays for which the first level
!     should be computed.

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

!     Definintion of the subgrid

      Integer, Intent (In)             :: levdim1, levdim2, levdim3

!     Dimension of the first MG level
!
! !OUTPUT PARAMETERS:
!
      Double Precision, Intent (Out)   :: chmin (0:levdim1, 0:levdim2, 0:levdim3)
!
!     Minimum extension of bounding box of each coordinate cell.
!
      Double Precision, Intent (Out)   :: chmax (0:levdim1, 0:levdim2, 0:levdim3)
!
!     Maximum extension of bounding box of each coordinate cell.
!
      Double Precision, Intent (Out)   :: midp  (0:levdim1, 0:levdim2, 0:levdim3)
!
!     Mid point of bounding box of each coordinate cell.
!
      integer, Intent (Out)            :: ierror

!     Returns the error code of PSMILe_MG_First_subgrid_3d_dble
!             ierror = 0 : No error
!             ierror > 0 : Severe error
!
! !DEFINED PARAMETERS:
!
      Double Precision, parameter      :: epsp1 = 1.0d0 + 5d-6
!
! !LOCAL VARIABLES
!
      integer                          :: i, ibeg, iend
      integer                          :: j, jbeg, jend
      integer                          :: k, kbeg, kend
!
      Double Precision                 :: r_nbr
!
! !DESCRIPTION:
!
! Subroutine "PSMILe_MG_First_subgrid_3d_dble" computes the
! bounding box for each control volume (cell) of the 3-dimensional grid
! defined by the corner array "array".
! The bounding boxes and mid points are used in the Multigrid search.
!
! !REVISION HISTORY:
!
!   Date      Programmer   Description
! ----------  ----------   -----------
! 03.06.25    H. Ritzdorf  created
!
!EOP
!----------------------------------------------------------------------
!
!  $Id: psmile_mg_first_subgrid_3d_dble.F90 2325 2010-04-21 15:00:07Z valcke $
!  $Author: valcke $
!
   Character(len=len_cvs_string), save :: mycvs = 
       '$Id: psmile_mg_first_subgrid_3d_dble.F90 2325 2010-04-21 15:00:07Z valcke $'
!
!----------------------------------------------------------------------
!
#ifdef VERBOSE
      print *, trim(ch_id), ': PSMILe_MG_First_subgrid_3d_dble'

      call psmile_flushstd
#endif /* VERBOSE */
!
!  Initialization
!
      ierror = 0
!
      ibeg = range (1, 1)
      iend = min (range(1, 1)+levdim1, range (2, 1))
!
      jbeg = range (1, 2)
      jend = min (range(1, 2)+levdim2, range (2, 2))
!
      kbeg = range (1, 3)
      kend = min (range(1, 3)+levdim3, range (2, 3))
!
      r_nbr = 1.0d0 / nbr_corners
!
!     ... Compute minimum extension of grid cell
!         ??? wirklich MINVAL, MAXVAL und SUM verwenden, oder
!             per hand von aussen den Werte berechnen ?
!
         do k = kbeg, kend
            do j = jbeg, jend
!cdir vector
               do i = ibeg, iend
               chmin (i-range(1,1), j-range(1,2), k-range(1,3)) = &
                  MINVAL (array (i, j, k, :))
               enddo
            enddo
         enddo
!
!     ... Compute maximum extension of grid cell
!
         do k = kbeg, kend
            do j = jbeg, jend
!cdir vector
               do i = ibeg, iend
               chmax (i-range(1,1), j-range(1,2), k-range(1,3)) = &
                  MAXVAL (array (i, j, k, :))
               enddo
            enddo
         enddo
!
!     ... Compute mid point of grid cell
!
         do k = kbeg, kend
            do j = jbeg, jend
!cdir vector
               do i = ibeg, iend
               midp (i-range(1,1), j-range(1,2), k-range(1,3))  = &
                  SUM (array (i, j, k, :)) * r_nbr
               enddo
            enddo
         enddo
!
!     ... Create dummy cells if necessary
!
      ibeg = ibeg - range(1,1)
      iend = iend - range(1,1)
      jbeg = jbeg - range(1,2)
      jend = jend - range(1,2)
      kbeg = kbeg - range(1,3)
      kend = kend - range(1,3)

      if (iend < levdim1) then
            do k = kbeg, kend
!cdir vector
               do j = jbeg, jend
               chmin (iend+1:levdim1, j, k) = chmax (iend, j, k) * epsp1
               chmax (iend+1:levdim1, j, k) = chmax (iend, j, k)
               midp  (iend+1:levdim1, j, k) = midp  (iend, j, k)
               enddo
            enddo
      endif
!
      if (jend < levdim2) then
            do k = kbeg, kend
!cdir vector
               do j = jend+1, levdim2
               chmin (ibeg:levdim1, j, k) = chmax (ibeg:levdim1, jend, k) * epsp1
               chmax (ibeg:levdim1, j, k) = chmax (ibeg:levdim1, jend, k)
               midp  (ibeg:levdim1, j, k) = midp  (ibeg:levdim1, jend, k)
               enddo
            enddo
      endif
!
      if (kend < levdim3) then
!cdir vector
            do k = kend+1, levdim3
            chmin (ibeg:levdim1, jbeg:levdim2, k) = &
            chmax (ibeg:levdim1, jbeg:levdim2, kend) * epsp1
            chmax (ibeg:levdim1, jbeg:levdim2, k) = &
            chmax (ibeg:levdim1, jbeg:levdim2, kend)
            midp  (ibeg:levdim1, jbeg:levdim2, k) = &
            midp  (ibeg:levdim1, jbeg:levdim2, kend)
            enddo
      endif
!
!===> All done
!
#ifdef VERBOSE
      print *, trim(ch_id), ': PSMILe_MG_First_subgrid_3d_dble eof', &
                            ': ierror =', ierror

      call psmile_flushstd
#endif /* VERBOSE */
!
      end subroutine PSMILe_MG_First_subgrid_3d_dble