TURB_LesHij.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: Get the Leonard stress of the dynamic mixed model.
!
! Description: In this subroutine we determine the H_{ij} components in the
!              dynamic mixed formulation, containing the scale similarity 
!              contribution. We assume L_{ij} to be given already on the cell 
!              face and we generate and subtract Hij from it. On exit, Hij
!              is stored in Lij.
!
! Input: region = data of current region
!        ijk    = ijk-face is being treated
!
! Output: Hij term in dynamic mixed formulation.
!
! Notes: Hij is stored in turb%lij.
!
!******************************************************************************
!
! $Id: TURB_LesHij.F90,v 1.8 2008/12/06 08:44:41 mtcampbe Exp $
!
! Copyright: (c) 2001 by the University of Illinois
!
!******************************************************************************

SUBROUTINE turb_leshij( region,ijk )

  USE moddatatypes
  USE moddatastruct, ONLY : t_region
  USE modglobal, ONLY     : t_global
  USE turb_modinterfaces, ONLY : turb_leslij
#ifdef RFLO
  USE modinterfaces, ONLY : rflo_getdimensdummynodes, rflo_getnodeoffset
  USE turb_modinterfaces, ONLY : turb_flolesunifiltff, turb_flolesgenfiltff, & 
                                 turb_flolesgenc2f

#include "Indexing.h"
#endif
  USE moderror
  USE turb_modparameters
  IMPLICIT NONE

! ... parameters
#ifdef RFLO
  TYPE(t_region)          :: region
#endif
#ifdef RFLU
  TYPE(t_region), POINTER :: region
#endif
  INTEGER                 :: ijk

! ... loop variables
  INTEGER :: i, j, k, ijkn

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

  INTEGER :: ibn,ien,idbeg,idend
  INTEGER :: ndel(diri:dirk),tndel(diri:dirk)
  REAL(RFREAL), POINTER :: fvar(:,:),ffvar(:,:),lij(:,:),field(:,:)
  REAL(RFREAL), POINTER :: hij(:,:)
#ifdef RFLO
  INTEGER :: idnbeg,idnend,jdnbeg,jdnend,kdnbeg,kdnend
  INTEGER :: ilev,inoff,ijnoff
#endif

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

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

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

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

#ifdef RFLO
  ilev   = region%currLevel
  CALL rflo_getdimensdummynodes( region,ilev,idnbeg,idnend, &
                                 jdnbeg,jdnend,kdnbeg,kdnend )
  CALL rflo_getnodeoffset( region,ilev,inoff,ijnoff )
  fvar  => region%levels(ilev)%turb%fVar
  ffvar => region%levels(ilev)%turb%ffVar
  lij   => region%levels(ilev)%turb%lij

  IF (ijk==diri) THEN
    field => region%levels(ilev)%turb%fISij
  ELSEIF (ijk==dirj) THEN 
    field => region%levels(ilev)%turb%fJSij
  ELSEIF (ijk==dirk) THEN 
    field => region%levels(ilev)%turb%fKSij
  ENDIF
  ibn   = indijk(idnbeg,jdnbeg,kdnbeg,inoff,ijnoff)
  ien   = indijk(idnend,jdnend,kdnend,inoff,ijnoff)
#endif

! allocate temporary arrays

  ALLOCATE( hij(e11:e33,ibn:ien) )

! test filter width is twice bar filter width

  ndel(diri) = region%turbInput%filterWidth(diri)
#ifdef RFLO
  ndel(dirj) = region%turbInput%filterWidth(dirj)
  ndel(dirk) = region%turbInput%filterWidth(dirk)
#endif

  tndel(diri)=2*ndel(diri)
#ifdef RFLO
  tndel(dirj)=2*ndel(dirj)
  tndel(dirk)=2*ndel(dirk)
#endif

! for efficiency we keep 1/tbRho in tbRho, leftover from lesLij;
! with this information we can generate the first term in Hij
! hij-component; store tbRhoUi*tbRhoUj/tbRho in hij:

  DO ijkn = ibn,ien
    hij(e11,ijkn)=ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_xmom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
    hij(e12,ijkn)=ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_ymom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
    hij(e13,ijkn)=ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_zmom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
    hij(e22,ijkn)=ffvar(cv_turb_ymom,ijkn)*ffvar(cv_turb_ymom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
    hij(e23,ijkn)=ffvar(cv_turb_ymom,ijkn)*ffvar(cv_turb_zmom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
    hij(e33,ijkn)=ffvar(cv_turb_zmom,ijkn)*ffvar(cv_turb_zmom,ijkn)* &
                  ffvar(cv_turb_dens,ijkn)
  ENDDO

! and we double-filter this, first with bar filter

  idbeg = e11
  idend = e33
#ifdef RFLO
  IF (region%turbInput%filterType == filtype_uniform) then
    CALL turb_flolesunifiltff( region,ijk,ndel,idbeg,idend,hij,field )
  ELSE
    CALL turb_flolesgenfiltff( region,ijk,ndel,idbeg,idend,hij,field )
  ENDIF
#endif

! then with test-filter:

  idbeg = e11
  idend = e33
#ifdef RFLO
  IF (region%turbInput%filterType == filtype_uniform) then
    CALL turb_flolesunifiltff( region,ijk,tndel,idbeg,idend,field,hij )
  ELSE
    CALL turb_flolesgenfiltff( region,ijk,tndel,idbeg,idend,field,hij )
  ENDIF
#endif

! subtract this first term component of hij from lij
! and we put 1/tbRho back to tbRho for the next treatment

  DO ijkn = ibn,ien
    lij(e11,ijkn) = lij(e11,ijkn) - hij(e11,ijkn)
    lij(e12,ijkn) = lij(e12,ijkn) - hij(e12,ijkn)
    lij(e13,ijkn) = lij(e13,ijkn) - hij(e13,ijkn)
    lij(e22,ijkn) = lij(e22,ijkn) - hij(e22,ijkn)
    lij(e23,ijkn) = lij(e23,ijkn) - hij(e23,ijkn)
    lij(e33,ijkn) = lij(e33,ijkn) - hij(e33,ijkn)
    ffvar(cv_turb_dens,ijkn) = 1._rfreal/ffvar(cv_turb_dens,ijkn)
  ENDDO

! next, we build second term of the first component in hij for which 
! we apply bar-filter and test-filter to test-bar filtered fields;
! first we filter tbRho and tbRhoUi with bar-filter and store in fVar

  idbeg = cv_turb_dens
  idend = cv_turb_zmom
#ifdef RFLO
  IF (region%turbInput%filterType == filtype_uniform) then
    CALL turb_flolesunifiltff( region,ijk,ndel,idbeg,idend,ffvar,fvar )
  ELSE
    CALL turb_flolesgenfiltff( region,ijk,ndel,idbeg,idend,ffvar,fvar )
  ENDIF
#endif

! then apply test filter and store in ffVar again

#ifdef RFLO
  IF (region%turbInput%filterType == filtype_uniform) then
    CALL turb_flolesunifiltff( region,ijk,tndel,idbeg,idend,fvar,ffvar )
  ELSE
    CALL turb_flolesgenfiltff( region,ijk,tndel,idbeg,idend,fvar,ffvar )
  ENDIF
#endif

! we can now finish the first contribution in hij and subtract it from lij
! for efficiency we store 1/tb(tbRho) in tb(tbRho)

  DO ijkn = ibn,ien
    ffvar(cv_turb_dens,ijkn) = 1._rfreal/ffvar(cv_turb_dens,ijkn)
    lij(e11,ijkn) = lij(e11,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_xmom,ijkn)

    lij(e12,ijkn) = lij(e12,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_ymom,ijkn)

    lij(e13,ijkn) = lij(e13,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_xmom,ijkn)*ffvar(cv_turb_zmom,ijkn)

    lij(e22,ijkn) = lij(e22,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_ymom,ijkn)*ffvar(cv_turb_ymom,ijkn)

    lij(e23,ijkn) = lij(e23,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_ymom,ijkn)*ffvar(cv_turb_zmom,ijkn)

    lij(e33,ijkn) = lij(e33,ijkn) + ffvar(cv_turb_dens,ijkn)* &
                    ffvar(cv_turb_zmom,ijkn)*ffvar(cv_turb_zmom,ijkn)
  ENDDO

! we now come to the second component which involves the test-filtered
! Bardina term; first recompute density and velocities at faces and 
! generate the Bardina term and put the result in field

#ifdef RFLO
  CALL turb_flolesgenc2f( region,ijk )
#endif

  CALL turb_leslij( region,ijk,ndel,field )

! test-filter Bardina component stored in field and put in hij

  idbeg = e11
  idend = e33
#ifdef RFLO
  IF (region%turbInput%filterType == filtype_uniform) THEN
    CALL turb_flolesunifiltff( region,ijk,tndel,idbeg,idend,field,hij )
  ELSE
    CALL turb_flolesgenfiltff( region,ijk,tndel,idbeg,idend,field,hij )
  ENDIF
#endif

! complete hij and subtract it from lij

  DO ijkn = ibn,ien
    lij(e11,ijkn) = lij(e11,ijkn) + hij(e11,ijkn)
    lij(e12,ijkn) = lij(e12,ijkn) + hij(e12,ijkn)
    lij(e13,ijkn) = lij(e13,ijkn) + hij(e13,ijkn)
    lij(e22,ijkn) = lij(e22,ijkn) + hij(e22,ijkn)
    lij(e23,ijkn) = lij(e23,ijkn) + hij(e23,ijkn)
    lij(e33,ijkn) = lij(e33,ijkn) + hij(e33,ijkn)
  ENDDO

! deallocate temporary arrays

  DEALLOCATE( hij )

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

  CALL deregisterfunction( global )

END SUBROUTINE turb_leshij

!******************************************************************************
!
! RCS Revision history:
!
! $Log: TURB_LesHij.F90,v $
! Revision 1.8  2008/12/06 08:44:41  mtcampbe
! Updated license.
!
! Revision 1.7  2008/11/19 22:17:54  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.6  2004/07/30 22:35:11  wasistho
! replaced floLesUniC2F by floLesGenC2F
!
! Revision 1.5  2004/05/28 02:01:06  wasistho
! update unstructured grid LES
!
! Revision 1.4  2004/03/27 02:16:42  wasistho
! compiled with Rocflu
!
! Revision 1.3  2004/03/19 02:50:25  wasistho
! prepared for RFLU
!
! Revision 1.2  2004/03/12 02:55:35  wasistho
! changed rocturb routine names
!
! Revision 1.1  2004/03/05 04:37:00  wasistho
! changed nomenclature
!
! Revision 1.2  2003/05/15 02:57:06  jblazek
! Inlined index function.
!
! Revision 1.1  2002/10/14 23:55:30  wasistho
! Install Rocturb
!
!
!******************************************************************************