RADI_FlimDiffFluxPatch.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 FLDTRAN diffusion flux: c*lamda(Er)/Kr*d_j(Er)
!          through a patch
!
! Description: this routine works in the same way as RADI_FlimDiffFlux
!              but applied on region patches
!
! Input: region  = data of current region.
!        patch   = current patch.
!
! Output: region%levels%radi%diss = diffusion flux added to FLD dissipation.
!
! Notes: none.
!
!******************************************************************************
!
! $Id: RADI_FlimDiffFluxPatch.F90,v 1.4 2008/12/06 08:44:37 mtcampbe Exp $
!
! Copyright: (c) 2003 by the University of Illinois
!
!******************************************************************************

SUBROUTINE radi_flimdifffluxpatch( 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 radi_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) :: beta, sounda, flima, coefa, diffcoef, modsf
  REAL(RFREAL) :: radenx, radeny, radenz, fd, sf(3), ac0, ac1
  REAL(RFREAL), POINTER :: dv(:,:), rcv(:,:), rdiss(:,:), flim(:), coef(:,:)
  REAL(RFREAL), POINTER :: qr(:)

#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,'RADI_FlimDiffFluxPatch',&
  'RADI_FlimDiffFluxPatch.F90' )

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

  bctype = patch%bcType

#ifdef RFLO
  ilev   = region%currLevel
  lbound = patch%lbound
  
  dv    => region%levels(ilev)%mixt%dv
  rcv   => region%levels(ilev)%radi%cv
  rdiss => region%levels(ilev)%radi%diss  
  flim  => region%levels(ilev)%radi%fluxLim
  coef  => region%levels(ilev)%radi%radCoef
#endif  
#ifdef RFLU
  dv    => region%mixt%dv
  rcv   => region%radi%cv
  rdiss => region%radi%diss   
  flim  => region%radi%fluxLim
  coef  => region%radi%radCoef
#endif

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

  beta   = region%mixtInput%betrk(region%irkStep)

#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)%radi%gradi
    qr    => region%levels(ilev)%radi%qri
  ELSE IF (lbound==3 .OR. lbound==4) THEN
    avgco => region%levels(ilev)%grid%c2fCoJ
    sface => region%levels(ilev)%grid%sj
    grad  => region%levels(ilev)%radi%gradj
    qr    => region%levels(ilev)%radi%qrj
  ELSE
    avgco => region%levels(ilev)%grid%c2fCoK
    sface => region%levels(ilev)%grid%sk
    grad  => region%levels(ilev)%radi%gradk
    qr    => region%levels(ilev)%radi%qrk
  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%radi%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
          sounda   = ac0*dv(dv_mixt_soun,ijkcb0)+ac1*dv(dv_mixt_soun,ijkcd)
          flima    = ac0*flim(ijkcb0)+ac1*flim(ijkcd)
          coefa    = ac0*coef(ijkcb0,radi_coeff_extinct)+ &
                     ac1*coef(ijkcd,radi_coeff_extinct)
          diffcoef = sounda*flima/coefa 

#ifdef RFLO
          radenx = grad(gr_radi_ex,ijknb)
          radeny = grad(gr_radi_ey,ijknb)
          radenz = grad(gr_radi_ez,ijknb)
#endif
#ifdef RFLU
          radenx = grad(xcoord,gr_radi_ex,ijknbg)
          radeny = grad(ycoord,gr_radi_ex,ijknbg)
          radenz = grad(zcoord,gr_radi_ex,ijknbg)
#endif

          fd = diffcoef* &
               (radenx*sf(1)+radeny*sf(2)+radenz*sf(3))

          rdiss(cv_radi_ener,ijkcb0) = rdiss(cv_radi_ener,ijkcb0)+fd*beta

! ------- store rad. flux normal to surf. in positive direction (hence sgn* )
          modsf     = sqrt( sf(1)*sf(1) + sf(2)*sf(2) + sf(3)*sf(3) )  
          qr(ijknb) = sgn*fd/modsf
#ifdef RFLO
        ENDDO   ! i
      ENDDO     ! j
    ENDDO       ! k
  ENDIF
#endif 
      
#ifdef RFLU
  ENDDO         ! iC
#endif

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

  CALL deregisterfunction( region%global )

END SUBROUTINE radi_flimdifffluxpatch

!******************************************************************************
!
! RCS Revision history:
! 
! $Log: RADI_FlimDiffFluxPatch.F90,v $
! Revision 1.4  2008/12/06 08:44:37  mtcampbe
! Updated license.
!
! Revision 1.3  2008/11/19 22:17:50  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.2  2006/08/19 15:40:15  mparmar
! Removed bf2bg
!
! Revision 1.1  2004/09/30 17:49:10  wasistho
! prepared for full FLD radiation model
!
!
!
!******************************************************************************