TURB_floLesAverageFace.F90

Go to the documentation of this file.

! *********************************************************************
! * Rocstar Simulation Suite                                          *
! * Copyright@2015, Illinois Rocstar LLC. All rights reserved.        *
! *                                                                   *
! * Illinois Rocstar LLC                                              *
! * Champaign, IL                                                     *
! * www.illinoisrocstar.com                                           *
! * sales@illinoisrocstar.com                                         *
! *                                                                   *
! * License: See LICENSE file in top level of distribution package or *
! * http://opensource.org/licenses/NCSA                               *
! *********************************************************************
! *********************************************************************
! * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,   *
! * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES   *
! * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND          *
! * NONINFRINGEMENT.  IN NO EVENT SHALL THE CONTRIBUTORS OR           *
! * COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER       *
! * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,   *
! * Arising FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE    *
! * USE OR OTHER DEALINGS WITH THE SOFTWARE.                          *
! *********************************************************************
!******************************************************************************
!
! Purpose: Average face variable to another face variable whose size depends
!          on the number of homogeneous directions. The field to be filtered 
!          is contained in faceVar and on exit the averaged data is stored in 
!          avgFaceVar
!
! Description: Only the inner field (including boundaries) is averaged.
!              If the averaging direction (denoted by homDir) is the same as 
!              the working face direction (denoted by ijk) the averaging 
!              interval ranges from the first face up to the last face, 
!              otherwise from the first center to the last center, for instance 
!              if homDir=(1,0,0) and ijk=DIRI the averaging interval is 
!              [ipnbeg:ipnend], whereas if homDir=(0,1,0) and ijk=DIRI it is 
!              [jpcbeg:jpcend]. If there is no homogeneous direction, a local
!              averaging by applying a test filter procedure and the averaged
!              data is contained in a volume array. If there is only one
!              averaging (homogeneous) direction, the averaged data is
!              contained in a surface array. For two averaging directions, it 
!              is a line array, whereas for three averaging directions 
!              (isotropic flows) the averaged data has one single value. On 
!              exit the averaged data, eiter volume, surface, line array or 
!              point (scalar) is stored in 1D array avgFaceVar. 
!              For each ijk case we perform line averaging-operation in three 
!              consecutive sweeps corresponding to the three computational 
!              directions, depending on whether a particular direction 
!              homogeneous (homDir(i)=1) or not (=0). Averaging is only 
!              performed in homogeneous direction. We use trapezoidal rule if 
!              the averaging direction is normal to the faces and ensemble 
!              rule otherwise. homDir(i)=0 for all three directions implies no 
!              homogeneous direction is present and local face-face filtering 
!              is used in this case
!
! Input: region  = current region data
!        ijk     = ijk-face is being treated
!        faceVar = face variables to be averaged 
!
! Output: avgFaceVar = the resulting averaged variables.
!
! Notes: Several generic averaging routines contained at the end of this file.
!
!******************************************************************************
!
! $Id: TURB_floLesAverageFace.F90,v 1.6 2008/12/06 08:44:43 mtcampbe Exp $
!
! Copyright: (c) 2001 by the University of Illinois
!
!******************************************************************************

SUBROUTINE turb_flolesaverageface( region,ijk,faceVar,avgFaceVar )

  USE moddatatypes
  USE moddatastruct, ONLY : t_region
  USE modglobal, ONLY     : t_global
  USE modinterfaces, ONLY : rflo_getdimensphysnodes,rflo_getdimensdummynodes, &
                            rflo_getnodeoffset
  USE turb_modinterfaces, ONLY : turb_flolesunifiltff, turb_flolesgenfiltff
  USE modturbulence
  USE moderror
  USE modparameters
  USE turb_modparameters
  IMPLICIT NONE

#include "Indexing.h"

! ... parameters
  TYPE(t_region)        :: region
  INTEGER               :: ijk
  REAL(RFREAL), POINTER :: facevar(:),avgfacevar(:)

! ... loop variables
  INTEGER :: i, j, k, ijkn, ijn, ikn, jkn

! ... local variables
  CHARACTER(CHRLEN)       :: rcsidentstring
  TYPE(t_global), POINTER :: global

  INTEGER           :: ibeg,iend,jbeg,jend,kbeg,kend
  INTEGER           :: ipnbeg,ipnend,jpnbeg,jpnend,kpnbeg,kpnend
  INTEGER           :: idnbeg,idnend,jdnbeg,jdnend,kdnbeg,kdnend
  INTEGER           :: ilev,inoff,jnoff,knoff,ijnoff,iknoff,jknoff
  INTEGER           :: ibn,ien,homsum,minl,maxl,actsize,refsize
  INTEGER           :: tndel(diri:dirk),homdir(diri:dirk)

  REAL(RFREAL)      :: avgpoint
  REAL(RFREAL), POINTER :: avgline(:),avgsurf(:)
  REAL(RFREAL), POINTER :: onedvec(:,:),avgonedvec(:,:)

!******************************************************************************

  rcsidentstring = '$RCSfile: TURB_floLesAverageFace.F90,v $'

  global => region%global
  CALL registerfunction( global,'TURB_FloLesAverageFace',&
  'TURB_floLesAverageFace.F90' )

! get indices and pointers --------------------------------------------------

  ilev   = region%currLevel
  CALL rflo_getdimensphysnodes( region,ilev,ipnbeg,ipnend, &
                                jpnbeg,jpnend,kpnbeg,kpnend )
  CALL rflo_getdimensdummynodes( region,ilev,idnbeg,idnend, &
                                jdnbeg,jdnend,kdnbeg,kdnend )
  CALL rflo_getnodeoffset( region,ilev,inoff,ijnoff )
  ibn    = indijk(idnbeg,jdnbeg,kdnbeg,inoff,ijnoff)
  ien    = indijk(idnend,jdnend,kdnend,inoff,ijnoff)
  jnoff  = jdnend-jdnbeg+1
  knoff  = kdnend-kdnbeg+1
  iknoff = inoff*(kdnend-kdnbeg+1)
  jknoff = jnoff*(kdnend-kdnbeg+1)
  minl   = min(idnbeg,jdnbeg,kdnbeg)
  maxl   = max(idnend,jdnend,kdnend)

  IF (ijk==diri) THEN
    ibeg = ipnbeg
    iend = ipnend
    jbeg = jpnbeg
    jend = jpnend-1
    kbeg = kpnbeg
    kend = kpnend-1
  ELSEIF (ijk==dirj) THEN
    ibeg = ipnbeg
    iend = ipnend-1
    jbeg = jpnbeg
    jend = jpnend
    kbeg = kpnbeg
    kend = kpnend-1
  ELSEIF (ijk==dirk) THEN
    ibeg = ipnbeg
    iend = ipnend-1
    jbeg = jpnbeg
    jend = jpnend-1
    kbeg = kpnbeg
    kend = kpnend
  ENDIF

  homdir(:) = region%turbInput%homDir(:)
  homsum    = homdir(diri) + homdir(dirj) + homdir(dirk)

! for local point-averaging we use test filter width twice grid-spacing

! - Volume filtering with possibly more terms filtered across regions
!  tNDel(DIRI)=2
!  tNDel(DIRJ)=2
!  tNDel(DIRK)=2

! - Surface filtering for least terms filtered across region boundaries
  IF (ijk==diri) THEN
    tndel(diri)=0
    tndel(dirj)=2
    tndel(dirk)=2
  ELSEIF (ijk==dirj) THEN
    tndel(diri)=2
    tndel(dirj)=0
    tndel(dirk)=2
  ELSEIF (ijk==dirk) THEN
    tndel(diri)=2
    tndel(dirj)=2
    tndel(dirk)=0
  ENDIF

  IF (region%turbInput%filterType == filtype_nonunif) THEN
    IF ((tndel(diri)/=2 .AND. tndel(diri)/=0) .OR. &
        (tndel(dirj)/=2 .AND. tndel(dirj)/=0) .OR. &
        (tndel(dirk)/=2 .AND. tndel(dirk)/=0)) THEN 
      CALL errorstop( global,err_turb_fixparam,__line__, &
      'filter coeff. possibly not allocated if tNDel =/ 2 for nonunif. filter.')
    ENDIF
  ENDIF

! allocate temporary arrays
    
  ALLOCATE( avgline(minl:maxl),avgsurf(ibn:ien) )
  ALLOCATE( onedvec(1,ibn:ien ),avgonedvec(1,ibn:ien) )

! perform averaging process depending on number of homogeneous directions

  IF (homsum == 3) THEN

! - all three directions are homogeneous

    actsize=1

    IF (homdir(diri)==active .AND. homdir(dirj)==active) THEN

      IF (ijk==diri) THEN
! ----- averaging in all direction; trapezoidal in i-direction and ensemble 
!       in other directions
        CALL trapezoid3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL ensemble2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )
        CALL ensemble1loop( kbeg,kend,avgline,avgpoint )

      ELSEIF (ijk==dirj) THEN
! ----- averaging in all direction; trapezoidal in j-direction and ensemble 
!       in other directions
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL trapezoid2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )
        CALL ensemble1loop( kbeg,kend,avgline,avgpoint )

      ELSEIF (ijk==dirk) THEN
! ----- averaging in all direction; trapezoidal in k-direction and ensemble 
!       in other directions
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL ensemble2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )
        CALL trapezoid1loop( kbeg,kend,avgline,avgpoint )

      ENDIF 
    ELSE
      CALL errorstop( global,err_turb_lesinput,__line__,'homDir is not consistent.')
    ENDIF       

    avgfacevar = avgpoint

  ELSEIF (homsum == 2) THEN

    IF (homdir(dirj)==active .AND. homdir(dirk)==active) THEN
      actsize = inoff

      IF (ijk==diri) THEN
! ----- averaging in j and k direction; ensemble in both direcions
        CALL ensemble3ljki( jbeg,jend,kbeg,kend,ibeg,iend,facevar,avgsurf )
        CALL ensemble2loops( kbeg,kend,ibeg,iend,knoff,avgsurf,avgline )

      ELSEIF (ijk==dirj) THEN
! ----- averaging in j and k direction; trapezoidal in j direction and
!       ensemble in k direcions
        CALL trapezoid3ljki( jbeg,jend,kbeg,kend,ibeg,iend,facevar,avgsurf )
        CALL ensemble2loops( kbeg,kend,ibeg,iend,knoff,avgsurf,avgline )

      ELSEIF (ijk==dirk) THEN
! ----- averaging in j and k direction; trapezoidal in k direction and
!       ensemble in j direcions
        CALL ensemble3ljki( jbeg,jend,kbeg,kend,ibeg,iend,facevar,avgsurf )
        CALL trapezoid2loops( kbeg,kend,ibeg,iend,knoff,avgsurf,avgline )

      ENDIF

! --- copy to avgFaceVar
      DO k=kbeg,kend
        DO j=jbeg,jend
          DO i=ibeg,iend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgline(i)
          ENDDO
        ENDDO
      ENDDO

    ELSEIF (homdir(diri)==active .AND. homdir(dirk)==active) THEN
      actsize = jnoff

      IF (ijk==diri) THEN
! ----- averaging in i and k direction; trapezoidal in i direction and
!       ensemble in k direcions:
        CALL trapezoid3likj( ibeg,iend,kbeg,kend,jbeg,jend,facevar,avgsurf )
        CALL ensemble2loops( kbeg,kend,jbeg,jend,knoff,avgsurf,avgline )

      ELSEIF (ijk==dirj) THEN
! ----- averaging in i and k direction; ensemble in both direcions
        CALL ensemble3likj( ibeg,iend,kbeg,kend,jbeg,jend,facevar,avgsurf )
        CALL ensemble2loops( kbeg,kend,jbeg,jend,knoff,avgsurf,avgline )

      ELSEIF (ijk==dirk) THEN
! ----- averaging in i and k direction; trapezoidal in k direction and
!       ensemble in i direcions
        CALL ensemble3likj( ibeg,iend,kbeg,kend,jbeg,jend,facevar,avgsurf )
        CALL trapezoid2loops( kbeg,kend,jbeg,jend,knoff,avgsurf,avgline )

      ENDIF

! --- copy to avgFaceVar
      DO k=kbeg,kend
        DO j=jbeg,jend
          DO i=ibeg,iend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgline(j)
          ENDDO
        ENDDO
      ENDDO

    ELSEIF (homdir(diri)==active .AND. homdir(dirj)==active) THEN
      actsize = knoff

      IF (ijk==diri) THEN
! ----- averaging in i and j direction; trapezoidal in i direction and
!       ensemble in j direcions
        CALL trapezoid3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL ensemble2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )

      ELSEIF (ijk==dirj) THEN
! ----- averaging in i and j direction; trapezoidal in j direction and
!       ensemble in i direcions
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL trapezoid2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )

      ELSEIF (ijk==dirk) THEN
! ----- averaging in i and j direction; ensemble in both directions
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )
        CALL ensemble2loops( jbeg,jend,kbeg,kend,jnoff,avgsurf,avgline )

      ENDIF 

! --- copy to avgFaceVar
      DO k=kbeg,kend
        DO j=jbeg,jend
          DO i=ibeg,iend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgline(k)
          ENDDO
        ENDDO
      ENDDO

    ELSE
      CALL errorstop( global,err_turb_lesinput,__line__,'homDir is not consistent.')
    ENDIF       

  ELSEIF (homsum == 1) THEN

    IF (homdir(dirj)/=active .AND. homdir(dirk)/=active) THEN
      actsize = jknoff

      IF (ijk==diri) THEN
! ----- trapezoidal averaging in i direction
        CALL trapezoid3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )

      ELSEIF (ijk==dirj) THEN
! ----- ensemble averaging in i direction
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )

      ELSEIF (ijk==dirk) THEN
! ----- ensemble averaging in i direction.
        CALL ensemble3lijk( ibeg,iend,jbeg,jend,kbeg,kend,facevar,avgsurf )

      ENDIF 

! --- copy to avgFaceVar
      DO k=kbeg,kend
        DO j=jbeg,jend
          jkn = indij(j ,k ,jnoff)
          DO i=ibeg,iend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgsurf(jkn)
          ENDDO
        ENDDO
      ENDDO

    ELSEIF (homdir(diri)/=active .AND. homdir(dirk)/=active) THEN
      actsize = iknoff

      IF (ijk==diri) THEN
! ----- ensemble averaging in j direction.
        CALL ensemble3ljik( jbeg,jend,ibeg,iend,kbeg,kend,facevar,avgsurf )

      ELSEIF (ijk==dirj) THEN
! ----- trapezoidal averaging in j direction
        CALL trapezoid3ljik( jbeg,jend,ibeg,iend,kbeg,kend,facevar,avgsurf )

      ELSEIF (ijk==dirk) THEN
! ----- ensemble averaging in j direction
        CALL ensemble3ljik( jbeg,jend,ibeg,iend,kbeg,kend,facevar,avgsurf )

      ENDIF 

! --- copy to avgFaceVar
      DO k=kbeg,kend
        DO i=ibeg,iend
          ikn = indij(i ,k ,inoff)
          DO j=jbeg,jend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgsurf(ikn)
          ENDDO
        ENDDO
      ENDDO

    ELSEIF (homdir(diri)/=active .AND. homdir(dirj)/=active) THEN
      actsize = ijnoff

      IF (ijk==diri) THEN
! ----- ensemble averaging in k direction
        CALL ensemble3lkij( kbeg,kend,ibeg,iend,jbeg,jend,facevar,avgsurf )

      ELSEIF (ijk==dirj) THEN
! ----- ensemble averaging in k direction.
        CALL ensemble3lkij( kbeg,kend,ibeg,iend,jbeg,jend,facevar,avgsurf )

      ELSEIF (ijk==dirk) THEN
! ----- trapezoidal averaging in k direction.
        CALL trapezoid3lkij( kbeg,kend,ibeg,iend,jbeg,jend,facevar,avgsurf )

      ENDIF 

! --- copy to avgFaceVar
      DO j=jbeg,jend
        DO i=ibeg,iend
          ijn = indij(i ,j ,inoff)
          DO k=kbeg,kend
            ijkn = indijk(i ,j ,k ,inoff,ijnoff)
            avgfacevar(ijkn) = avgsurf(ijn)
          ENDDO
        ENDDO
      ENDDO

    ELSE
      CALL errorstop( global,err_turb_lesinput,__line__,'homDir is not consistent.')
    ENDIF       

  ELSEIF (homsum == 0) THEN

    IF (homdir(diri)/=active .AND. homdir(dirj)/=active) THEN

      actsize = ijnoff*(kdnend-kdnbeg+1)

      DO i = ibn,ien
        onedvec(1,i) = facevar(i) 
      ENDDO
      IF (region%turbInput%filterType == filtype_uniform) THEN
        CALL turb_flolesunifiltff( region,ijk,tndel,1,1,onedvec,avgonedvec )
      ELSE
        CALL turb_flolesgenfiltff( region,ijk,tndel,1,1,onedvec,avgonedvec )
      ENDIF
      DO i = ibn,ien
        avgfacevar(i) = avgonedvec(1,i)
      ENDDO
    ELSE
      CALL errorstop( global,err_turb_lesinput,__line__,'homDir is not consistent.')
    ENDIF       
  ELSE
    CALL errorstop( global,err_turb_lesinput,__line__,'homDir is not consistent.')
  ENDIF

! check size of averaged variable array

!  IF (refSize /= actSize) THEN
!    CALL ErrorStop( global,ERR_TURB_LESINPUT,__LINE__,'coefDim is not consistent.')
!  ENDIF 

! deallocate temporary arrays

  DEALLOCATE( avgline,avgsurf,onedvec,avgonedvec )

! finalize --------------------------------------------------------------------

  CALL deregisterfunction( global )

! ==============================================================================
! summation subroutines
! ==============================================================================

CONTAINS

  SUBROUTINE trapezoid1loop( ib,ie,lineVar,pointVar )

! ... parameters
    INTEGER               :: ib,ie
    REAL(RFREAL), POINTER :: linevar(:)
    REAL(RFREAL)          :: pointvar  

! ... local variables
    INTEGER               :: i
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    pointvar = 0._rfreal

    sum = 0._rfreal
    DO i=ib+1,ie-1
      sum = sum + linevar(i)
    ENDDO
    sum = sum + 0.5_rfreal*(linevar(ib)+linevar(ie))
    pointvar = sum/REAL(ie-ib)

  END SUBROUTINE trapezoid1loop

!#####################################################################

  SUBROUTINE ensemble1loop( ib,ie,lineVar,pointVar )

! ... parameters
    INTEGER               :: ib,ie
    REAL(RFREAL), POINTER :: linevar(:)
    REAL(RFREAL)          :: pointvar  

! ... local variables
    INTEGER               :: i
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    pointvar = 0._rfreal

    sum = 0._rfreal
    DO i=ib,ie
      sum = sum + linevar(i)
    ENDDO
    pointvar = sum/REAL(ie-ib+1)

  END SUBROUTINE ensemble1loop

!#####################################################################

  SUBROUTINE trapezoid2loops( ib,ie,jb,je,nOff,surfVar,lineVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,noff
    REAL(RFREAL), POINTER :: surfvar(:),linevar(:)

! ... local variables
    INTEGER               :: i,j,ijn,ijnb,ijne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    linevar = 0._rfreal

    DO j=jb,je
      sum = 0._rfreal
      DO i=ib+1,ie-1
        ijn = indij(i ,j ,noff)
        sum = sum + surfvar(ijn)
      ENDDO
      ijnb = indij(ib ,j ,noff) 
      ijne = indij(ie ,j ,noff)
      sum  = sum + 0.5_rfreal*(surfvar(ijnb)+surfvar(ijne))
      linevar(j) = sum/REAL(ie-ib)
    ENDDO

  END SUBROUTINE trapezoid2loops

!#####################################################################

  SUBROUTINE ensemble2loops( ib,ie,jb,je,nOff,surfVar,lineVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,noff
    REAL(RFREAL), POINTER :: surfvar(:),linevar(:)

! ... local variables
    INTEGER               :: i,j,ijn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    linevar = 0._rfreal

    DO j=jb,je
      sum = 0._rfreal
      DO i=ib,ie
        ijn = indij(i ,j ,noff) 
        sum = sum+surfvar(ijn)
      ENDDO
      linevar(j) = sum/REAL(ie-ib+1)
    ENDDO

  END SUBROUTINE ensemble2loops

!#####################################################################

  SUBROUTINE trapezoid3lijk( ib,ie,jb,je,kb,ke,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,jkn,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO k=kb,ke
      DO j=jb,je
        jkn = indij(j ,k ,jnoff) 
        sum = 0._rfreal
        DO i=ib+1,ie-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff)
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(ib ,j ,k ,inoff,ijnoff) 
        ijkne = ijknb + ie-ib
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(jkn) = sum/REAL(ie-ib)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3lijk
!===============================================================
  SUBROUTINE trapezoid3likj( ib,ie,kb,ke,jb,je,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,kjn,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO j=jb,je
      DO k=kb,ke
        kjn = indij(k ,j ,knoff) 
        sum = 0._rfreal
        DO i=ib+1,ie-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(ib ,j ,k ,inoff,ijnoff) 
        ijkne = ijknb + ie-ib
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(kjn) = sum/REAL(ie-ib)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3likj
!===============================================================
  SUBROUTINE trapezoid3ljik( jb,je,ib,ie,kb,ke,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,ikn,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO k=kb,ke
      DO i=ib,ie
        ikn = indij(i ,k ,inoff) 
        sum = 0._rfreal
        DO j=jb+1,je-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(i ,jb ,k ,inoff,ijnoff) 
        ijkne = ijknb + (je-jb)*inoff
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(ikn) = sum/REAL(je-jb)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3ljik
!===============================================================
  SUBROUTINE trapezoid3ljki( jb,je,kb,ke,ib,ie,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,kin,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO i=ib,ie
      DO k=kb,ke
        kin = indij(k ,i ,knoff) 
        sum = 0._rfreal
        DO j=jb+1,je-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(i ,jb ,k ,inoff,ijnoff) 
        ijkne = ijknb + (je-jb)*inoff
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(kin) = sum/REAL(je-jb)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3ljki
!===============================================================
  SUBROUTINE trapezoid3lkij( kb,ke,ib,ie,jb,je,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,ijn,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO j=jb,je
      DO i=ib,ie
        ijn = indij(i ,j ,inoff) 
        sum = 0._rfreal
        DO k=kb+1,ke-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(i ,j ,kb ,inoff,ijnoff) 
        ijkne = ijknb + (ke-kb)*ijnoff
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(ijn) = sum/REAL(ke-kb)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3lkij
!===============================================================
  SUBROUTINE trapezoid3lkji( kb,ke,jb,je,ib,ie,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,jin,ijkn,ijknb,ijkne
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO i=ib,ie
      DO j=jb,je
        jin = indij(j ,i ,jnoff) 
        sum = 0._rfreal
        DO k=kb+1,ke-1
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        ijknb = indijk(i ,j ,kb ,inoff,ijnoff) 
        ijkne = ijknb + (ke-kb)*ijnoff
        sum   = sum + 0.5_rfreal*(volvar(ijknb)+volvar(ijkne))
        surfvar(jin) = sum/REAL(ke-kb)
      ENDDO
    ENDDO          

  END SUBROUTINE trapezoid3lkji

!#####################################################################

  SUBROUTINE ensemble3lijk( ib,ie,jb,je,kb,ke,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,jkn,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO k=kb,ke
      DO j=jb,je
        jkn = indij(j ,k ,jnoff) 
        sum = 0._rfreal
        DO i=ib,ie
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(jkn) = sum/REAL(ie-ib+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3lijk
!===============================================================
  SUBROUTINE ensemble3likj( ib,ie,kb,ke,jb,je,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,kjn,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO j=jb,je
      DO k=kb,ke
        kjn = indij(k ,j ,knoff) 
        sum = 0._rfreal
        DO i=ib,ie
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(kjn) = sum/REAL(ie-ib+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3likj
!===============================================================
  SUBROUTINE ensemble3ljik( jb,je,ib,ie,kb,ke,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,ikn,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO k=kb,ke
      DO i=ib,ie
        ikn = indij(i ,k ,inoff) 
        sum = 0._rfreal
        DO j=jb,je
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(ikn) = sum/REAL(je-jb+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3ljik
!===============================================================
  SUBROUTINE ensemble3ljki( jb,je,kb,ke,ib,ie,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,kin,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO i=ib,ie
      DO k=kb,ke
        kin = indij(k ,i ,knoff) 
        sum = 0._rfreal
        DO j=jb,je
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(kin) = sum/REAL(je-jb+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3ljki
!===============================================================
  SUBROUTINE ensemble3lkij( kb,ke,ib,ie,jb,je,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,ijn,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO j=jb,je
      DO i=ib,ie
        ijn = indij(i ,j ,inoff) 
        sum = 0._rfreal
        DO k=kb,ke
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(ijn) = sum/REAL(ke-kb+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3lkij
!===============================================================
  SUBROUTINE ensemble3lkji( kb,ke,jb,je,ib,ie,volVar,surfVar )

! ... parameters
    INTEGER               :: ib,ie,jb,je,kb,ke
    REAL(RFREAL), POINTER :: volvar(:),surfvar(:)

! ... local variables
    INTEGER               :: i,j,k,jin,ijkn
    REAL(RFREAL)          :: sum

! - start procedure ---------------------------------------

    surfvar = 0._rfreal

    DO i=ib,ie
      DO j=jb,je
        jin = indij(j ,i ,jnoff) 
        sum = 0._rfreal
        DO k=kb,ke
          ijkn = indijk(i ,j ,k ,inoff,ijnoff) 
          sum  = sum + volvar(ijkn)
        ENDDO
        surfvar(jin) = sum/REAL(ke-kb+1)
      ENDDO
    ENDDO          

  END SUBROUTINE ensemble3lkji

END SUBROUTINE turb_flolesaverageface

!******************************************************************************
!
! RCS Revision history:
!
! $Log: TURB_floLesAverageFace.F90,v $
! Revision 1.6  2008/12/06 08:44:43  mtcampbe
! Updated license.
!
! Revision 1.5  2008/11/19 22:17:55  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.4  2004/06/03 02:14:54  wasistho
! provided surface and volume (commented) local filtering
!
! Revision 1.3  2004/05/29 03:39:55  wasistho
! used local surface filtering instead of volume filtering in case of no hom. direction
!
! Revision 1.2  2004/03/12 02:55:36  wasistho
! changed rocturb routine names
!
! Revision 1.1  2004/03/08 23:35:45  wasistho
! changed turb nomenclature
!
! Revision 1.4  2004/02/28 01:19:10  wasistho
! reset filterwidth ratio for local averaging to 2
!
! Revision 1.3  2004/02/18 21:27:51  wasistho
! used smaller stencil for local averaging by filtering
!
! Revision 1.2  2003/05/15 02:57:06  jblazek
! Inlined index function.
!
! Revision 1.1  2002/10/14 23:55:29  wasistho
! Install Rocturb
!
!******************************************************************************