TURB_RansSAVisFluxPatch.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: compute SA viscous flux: nutot*d_j(tilde[nu])
!          through a patch
!
! Description: this routine works in the same way as TURB_RansSAVisFlux
!              but applied on region patches
!
! Input: region  = data of current region.
!        patch   = current patch.
!
! Output: region%levels%turb%diss = viscous flux added to SA dissipation.
!
! Notes: none.
!
!******************************************************************************
!
! $Id: TURB_RansSAVisFluxPatch.F90,v 1.10 2008/12/06 08:44:42 mtcampbe Exp $
!
! Copyright: (c) 2003 by the University of Illinois
!
!******************************************************************************

SUBROUTINE turb_ranssavisfluxpatch( region,patch )

  USE moddatatypes
  USE modbndpatch, ONLY   : t_patch
  USE moddatastruct, ONLY : t_region
#ifdef RFLO
  USE modinterfaces, ONLY : rflo_getpatchindices, rflo_getpatchdirection, &
                            rflo_getcelloffset, rflo_getnodeoffset

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

! ... parameters
  TYPE(t_region) :: region
  TYPE(t_patch)  :: patch

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

! ... local variables
  INTEGER :: bctype, ijkcb0, ijkcd, ijknb
  INTEGER :: ibeg, iend, jbeg, jend, kbeg, kend

  REAL(RFREAL)          :: cb2, opcb2, rsigma, beta, rhoa, mua, rnua, nuf, nuc
  REAL(RFREAL)          :: nutilx, nutily, nutilz, fd, sf(3), ac0, ac1
  REAL(RFREAL), POINTER :: cv(:,:), tv(:,:), tcv(:,:), tdiss(:,:)

#ifdef RFLO
  INTEGER :: inode, jnode, knode, idir, jdir, kdir
  INTEGER :: ilev, lbound, icoff, ijcoff, inoff, ijnoff, acid0, acid1
  REAL(RFREAL)          :: sgn
  REAL(RFREAL), POINTER :: avgco(:,:), sface(:,:), grad(:,:)
#endif
#ifdef RFLU
  INTEGER :: ifgbeg, ijknbg
  REAL(RFREAL), POINTER :: fn(:,:), grad(:,:,:) 
#endif

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

  CALL registerfunction( region%global,'TURB_RansSAVisFluxPatch',&
  'TURB_RansSAVisFluxPatch.F90' )

! get dimensions and pointers -------------------------------------------------

  bctype = patch%bcType

#ifdef RFLO
  ilev   = region%currLevel
  lbound = patch%lbound
  
  cv   => region%levels(ilev)%mixt%cv
  tv   => region%levels(ilev)%mixt%tv
  tcv  => region%levels(ilev)%turb%cv
  tdiss=> region%levels(ilev)%turb%diss  
#endif  
#ifdef RFLU
  cv   => region%mixt%cv
  tv   => region%mixt%tv
  tcv  => region%turb%cv
  tdiss=> region%turb%diss   
#endif

! get coefficients -----------------------------------------------------------

  cb2    = region%turbInput%const(mc_sa_cb2)
  rsigma = region%turbInput%const(mc_sa_rsig)
  beta   = region%mixtInput%betrk(region%irkStep)*rsigma
  opcb2  = 1._rfreal+cb2

#ifdef RFLO

! take the right face vector and make it point outwards ----------------------

  CALL rflo_getpatchindices( region,patch,ilev,ibeg,iend,jbeg,jend,kbeg,kend )
  CALL rflo_getpatchdirection( patch,idir,jdir,kdir )
  CALL rflo_getcelloffset( region,ilev,icoff,ijcoff )
  CALL rflo_getnodeoffset( region,ilev,inoff,ijnoff )

  sgn   = +1._rfreal
  inode = 0
  jnode = 0
  knode = 0
  acid0 = 2
  acid1 = 1
  IF (lbound==2 .OR. lbound==4 .OR. lbound==6) THEN
    sgn   = -1._rfreal
    inode = -idir
    jnode = -jdir
    knode = -kdir
    acid0 = 1
    acid1 = 2
  ENDIF

! get the appropriate face vector --------------------------------------------

  IF (lbound==1 .OR. lbound==2) THEN
    avgco => region%levels(ilev)%grid%c2fCoI
    sface => region%levels(ilev)%grid%si
    grad  => region%levels(ilev)%turb%gradi
  ELSE IF (lbound==3 .OR. lbound==4) THEN
    avgco => region%levels(ilev)%grid%c2fCoJ
    sface => region%levels(ilev)%grid%sj
    grad  => region%levels(ilev)%turb%gradj
  ELSE
    avgco => region%levels(ilev)%grid%c2fCoK
    sface => region%levels(ilev)%grid%sk
    grad  => region%levels(ilev)%turb%gradk
  ENDIF

! non-conforming region interface --------------------------------------------

  IF (bctype>=bc_regionint .AND. bctype<=bc_regionint+bc_range) THEN

  ELSE IF (bctype>=bc_regnonconf .AND. bctype<=bc_regnonconf+bc_range) THEN

! everything else

  ELSE

! flux in the direction normal to the patch ----------------------------------

    DO k=kbeg,kend
      DO j=jbeg,jend
        DO i=ibeg,iend

          ijkcb0 = indijk(i      ,j      ,k      ,icoff,ijcoff)  ! bnd cells
          ijkcd  = indijk(i-idir ,j-jdir ,k-kdir ,icoff,ijcoff)  ! dummy
          ijknb  = indijk(i+inode,j+jnode,k+knode,inoff,ijnoff)  ! bnd nodes
          ac0    = avgco(acid0,ijknb)
          ac1    = avgco(acid1,ijknb)
          sf(1)  = sgn*sface(xcoord,ijknb)
          sf(2)  = sgn*sface(ycoord,ijknb)
          sf(3)  = sgn*sface(zcoord,ijknb)
#endif
#ifdef RFLU
    ibeg   = 1
    iend   = patch%nBFaces
! TEMPORARY : removing usage of bf2bg from everywhere
!    ifgBeg = patch%bf2bg(BF2BG_BEG)
    ac0    = 0.5_rfreal
    ac1    = 0.5_rfreal

    grad => region%turb%bGradi
    fn   => patch%fn

    DO ic=ibeg,iend
          ijkcb0 = patch%bf2c(ic)
          ijkcd  = ijkcb0
          ijknb  = ic
          ijknbg = ic + ifgbeg-1

          sf(1)  = fn(xcoord,ijknb)*fn(xyzmag,ijknb)
          sf(2)  = fn(ycoord,ijknb)*fn(xyzmag,ijknb)
          sf(3)  = fn(zcoord,ijknb)*fn(xyzmag,ijknb)             
#endif
          rhoa = ac0*cv(cv_mixt_dens ,ijkcb0)+ac1*cv(cv_mixt_dens ,ijkcd)
          mua  = ac0*tv(tv_mixt_muel ,ijkcb0)+ac1*tv(tv_mixt_muel ,ijkcd)
          rnua = ac0*tcv(cv_sa_nutil ,ijkcb0)+ac1*tcv(cv_sa_nutil ,ijkcd)
          nuf  = (rnua + mua)/rhoa
          nuc  = (tcv(cv_sa_nutil ,ijkcb0) + tv(tv_mixt_muel ,ijkcb0))/ &
                  cv(cv_mixt_dens ,ijkcb0) 

#ifdef RFLO
          nutilx = grad(gr_sa_nutilx,ijknb)
          nutily = grad(gr_sa_nutily,ijknb)
          nutilz = grad(gr_sa_nutilz,ijknb)
#endif
#ifdef RFLU
          nutilx = grad(xcoord,gr_sa_nutilx,ijknbg)
          nutily = grad(ycoord,gr_sa_nutilx,ijknbg)
          nutilz = grad(zcoord,gr_sa_nutilx,ijknbg)
#endif

          fd = beta*(nuf*opcb2-nuc*cb2)* &
               (nutilx*sf(1)+nutily*sf(2)+nutilz*sf(3))

          tdiss(cv_sa_nutil,ijkcb0) = tdiss(cv_sa_nutil,ijkcb0)+fd

#ifdef RFLO
        ENDDO   ! i
      ENDDO     ! j
    ENDDO       ! k
  ENDIF
#endif 
      
#ifdef RFLU
  ENDDO         ! iC
#endif

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

  CALL deregisterfunction( region%global )

END SUBROUTINE turb_ranssavisfluxpatch

!******************************************************************************
!
! RCS Revision history:
! 
! $Log: TURB_RansSAVisFluxPatch.F90,v $
! Revision 1.10  2008/12/06 08:44:42  mtcampbe
! Updated license.
!
! Revision 1.9  2008/11/19 22:17:54  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.8  2006/08/19 15:41:01  mparmar
! Removed bf2bg
!
! Revision 1.7  2004/08/04 22:07:46  wasistho
! bugfixed: ac0 and ac1 are common to flo and flu
!
! Revision 1.6  2004/08/02 23:08:35  wasistho
! shift location of lines defining ac0 and ac1
!
! Revision 1.5  2004/08/02 21:55:46  wasistho
! replaced cell2face midpoint by linear averaging
!
! Revision 1.4  2004/03/25 04:40:41  wasistho
! prepared for RFLU
!
! Revision 1.3  2004/03/24 03:37:02  wasistho
! prepared for RFLU
!
! Revision 1.2  2004/03/20 03:28:29  wasistho
! prepared for RFLU
!
! Revision 1.1  2004/03/05 04:37:00  wasistho
! changed nomenclature
!
! Revision 1.4  2003/10/25 22:07:34  wasistho
! modified non-conservative diffusion term
!
! Revision 1.3  2003/10/20 20:28:28  wasistho
! made consistent with compressible SA formulation
!
! Revision 1.2  2003/10/10 20:35:13  wasistho
! multiplied SA viscous fluxes by 1/sigma through beta
!
! Revision 1.1  2003/10/07 02:17:03  wasistho
! initial installation of RaNS-SA and DES
!
!
!******************************************************************************