ExplicitMultistage.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: calculate solution at a new time level/iteration.
!
! Description: the governing equations are integrated in time using
!              an explicit, multistage (Runge-Kutta type) temporal scheme.
!              Spatial and temporal discretizations are independent of
!              each other based on the method of lines.
!
! Input: regions    = data of all regions
!        ftermNew   = compute new forcing term (true/false)
!        residFterm = add forcing term to residual (true/false).
!
! Output: regions%levels%mixt = new solution (and forcing term) after
!                               one time step/iteration.
!
! Notes: none.
!
! ******************************************************************************
!
! $Id: ExplicitMultistage.F90,v 1.15 2008/12/06 08:44:09 mtcampbe Exp $
!
! Copyright: (c) 2001-2005 by the University of Illinois
!
! ******************************************************************************

SUBROUTINE explicitmultistage( regions,ftermNew,residFterm )

  USE moddatatypes
  USE moddatastruct, ONLY : t_region
  USE modglobal, ONLY     : t_global
  USE moderror
  USE modmpi
  USE modparameters
#ifdef RFLO
  USE rflo_modboundaryconditions, ONLY : rflo_boundaryconditionsset, &
        rflo_boundaryconditionsrecv, rflo_boundaryconditionssend
  USE modinterfaces, ONLY : rflo_getdimensdummy, rflo_getcelloffset, &
        numericaldissipation, viscousfluxes, convectivefluxes, &
        mixtureproperties, rflo_residualsmoothing, rflo_checkvalidity, &
        rflo_zerodummycells, rflo_timestepinviscid, rflo_timestepviscous, &
        rflo_clearsendrequests, rflo_residualsmoothingcoeffs, &
        rflo_getboundaryvalues, rflo_sendboundaryvalues, sourceterms, &
        rflo_sendboundaryvaluesalpha, updatetbc, globalcommunicationmp
#include "Indexing.h"
#endif
#ifdef RFLU
  USE rflu_modconvertcv
  USE rflu_moddifferentiationcells
  USE rflu_modlimiters, ONLY: rflu_limitgradcellssimple
  USE rflu_modmpi
  USE modinterfaces, ONLY: convectivefluxes, &
                           mixtureproperties, &
                           numericaldissipation, &
                           rflu_checkpositivitywrapper, &
                           rflu_checkvaliditywrapper, &
                           rflu_setvars, & 
                           rflu_timestepinviscid, &
                           rflu_timestepviscous, &
                           rflu_zerovirtualcellvars, &
                           sourceterms, &
                           updatetbc, &
                           viscousfluxes
#endif
#ifdef PEUL
  USE modinterfaceseulerian, ONLY : peul_spectralradii, &
                                    peul_residualsmoothingcoeffs, &
                                    peul_residualsmoothing
#endif
#ifdef INRT
  USE modinterfacesinteract, ONLY : inrt_sourceterms
#endif
#ifdef RADI
  USE modinterfacesradiation, ONLY : radi_sourceterms
#endif
#ifdef PERI
  USE modinterfacesperiodic, ONLY : peri_sourceterms, peri_solutionupdate
!#ifdef RFLO 
!  USE PERI_ModHybridDES, ONLY : PERI_CoMeanCorrection 
!#endif 
#endif
#ifdef TURB
  USE modinterfacesturbulence, ONLY : turb_emsinit, turb_ransconvectivefluxes, &
                turb_ransnumericaldissipation,      turb_ranssourceterms, &
                turb_ranszerodummycells,            turb_solutionupdate
#ifdef RFLO
  USE modinterfacesturbulence, ONLY : turb_rflo_ransbndconditionssend, &    
          turb_rflo_ransbndconditionsset,     turb_rflo_ransbndconditionsrecv, &
          turb_rflo_ransclearsendrequests,    turb_rflo_ransressmoothing, &
          turb_rflo_ransressmoothingcoeff,    turb_rflo_ransspectralradii
#endif
#endif

  IMPLICIT NONE

#ifdef GENX
  include 'roccomf90.h'
#endif

! ... parameters
  TYPE(t_region), POINTER :: regions(:)

  LOGICAL :: ftermnew, residfterm

! ... loop variables
  INTEGER :: ireg, ireglocal, ic, istage

! ... local variables
#ifdef RFLO
  INTEGER :: idcbeg, idcend, jdcbeg, jdcend, kdcbeg, kdcend
  INTEGER :: ilev, icoff, ijcoff
#endif
  INTEGER :: ibc, iec, iectot, ldiss(5), flowmodel, gasmodel

  LOGICAL :: movegrid

#ifdef RFLO
  REAL(RFREAL) :: smoocf
#endif

  REAL(RFREAL)          :: cfl, ark(5), betrk(5), blend1, fac, adtv
  REAL(RFREAL)          :: alpha, time
  REAL(RFREAL), POINTER :: cv(:,:), cvold(:,:), dt(:), diss(:,:), rhs(:,:)
  REAL(RFREAL), POINTER :: vol(:), fterm(:,:)

  TYPE(t_global), POINTER :: global
#ifdef RFLU
  TYPE(t_region), POINTER :: pregion
#endif

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

  global => regions(1)%global

  CALL registerfunction( global,'ExplicitMultistage',&
  'ExplicitMultistage.F90' )

#ifdef RFLO
! set time for time-dependent BC data (for dual-time stepping)-----------------

  time  = global%currentTime + global%dtMin
  alpha = (time-global%timeStamp)/global%dTimeSystem
#endif

! loop over stages and regions ------------------------------------------------

  DO istage=1,regions(1)%global%nrkSteps
#ifdef RFLO
    DO ireg=1,global%nRegions
      IF (regions(ireg)%procid==global%myProcid .AND. & ! region active and
          regions(ireg)%active==active) THEN            ! on my processor
#endif
#ifdef RFLU
    DO ireglocal = 1,global%nRegionsLocal
      ireg = ireglocal
#endif
        regions(ireg)%irkStep = istage

! ----- get dimensions and pointers

        ldiss(:)  = regions(ireg)%mixtInput%ldiss(:)
        cfl       = regions(ireg)%mixtInput%cfl
        ark(:)    = regions(ireg)%mixtInput%ark(:)
        betrk(:)  = regions(ireg)%mixtInput%betrk(:)
        flowmodel = regions(ireg)%mixtInput%flowModel
        gasmodel  = regions(ireg)%mixtInput%gasModel
        movegrid  = regions(ireg)%mixtInput%moveGrid

#ifdef RFLO
        ilev = regions(ireg)%currLevel

        CALL rflo_getdimensdummy( regions(ireg),ilev,idcbeg,idcend, &
                                  jdcbeg,jdcend,kdcbeg,kdcend )
        CALL rflo_getcelloffset( regions(ireg),ilev,icoff,ijcoff )
        ibc    = indijk(idcbeg,jdcbeg,kdcbeg,icoff,ijcoff)
        iec    = indijk(idcend,jdcend,kdcend,icoff,ijcoff)
        iectot = iec

        smoocf = regions(ireg)%mixtInput%smoocf

        cv    => regions(ireg)%levels(ilev)%mixt%cv
        cvold => regions(ireg)%levels(ilev)%mixt%cvOld
        diss  => regions(ireg)%levels(ilev)%mixt%diss
        rhs   => regions(ireg)%levels(ilev)%mixt%rhs
        vol   => regions(ireg)%levels(ilev)%grid%vol
        dt    => regions(ireg)%levels(ilev)%dt
        IF (residfterm) fterm => regions(ireg)%levels(ilev)%mixt%fterm
#endif

#ifdef RFLU
        pregion => regions(ireg)

        ibc    = 1
        iec    = pregion%grid%nCells
        iectot = pregion%grid%nCellsTot

        cv    => pregion%mixt%cv
        cvold => pregion%mixt%cvOld
        diss  => pregion%mixt%diss
        rhs   => pregion%mixt%rhs
        vol   => pregion%grid%vol
        dt    => pregion%dt
        IF (residfterm) fterm => pregion%mixt%fterm
#endif

#ifdef RFLO
! ----- compute local time step (first stage only)

        IF (istage == 1) THEN
          IF (flowmodel == flow_euler) THEN
            CALL rflo_timestepinviscid( regions(ireg) )
          ELSE
            CALL rflo_timestepviscous( regions(ireg) )
          ENDIF
          IF (smoocf > 0._rfreal) THEN
            CALL rflo_residualsmoothingcoeffs( regions(ireg) )
          END IF ! smoocf
#ifdef PEUL
          IF (global%peulUsed) THEN
            CALL peul_spectralradii( regions(ireg) )

            IF (regions(ireg)%peulInput%smoocf > 0._rfreal) THEN
              CALL peul_residualsmoothingcoeffs( regions(ireg) )
            ENDIF ! smoocf
          ENDIF   ! peulUsed
#endif
#ifdef TURB
          IF (flowmodel == flow_navst .AND. &
              regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
            CALL turb_rflo_ransspectralradii( regions(ireg) )

            IF (regions(ireg)%turbInput%smoocf > 0._rfreal) THEN
              CALL turb_rflo_ransressmoothingcoeff( regions(ireg) )
            ENDIF ! smoocf
          ENDIF   ! flowModel
#endif
        ENDIF
#endif

! ----- store previous solution; set dissipation to zero (first stage)

        IF (istage == 1) THEN
          DO ic=ibc,iectot
            cvold(cv_mixt_dens,ic) = cv(cv_mixt_dens,ic)
            cvold(cv_mixt_xmom,ic) = cv(cv_mixt_xmom,ic)
            cvold(cv_mixt_ymom,ic) = cv(cv_mixt_ymom,ic)
            cvold(cv_mixt_zmom,ic) = cv(cv_mixt_zmom,ic)
            cvold(cv_mixt_ener,ic) = cv(cv_mixt_ener,ic)
            diss(cv_mixt_dens,ic)  = 0._rfreal
            diss(cv_mixt_xmom,ic)  = 0._rfreal
            diss(cv_mixt_ymom,ic)  = 0._rfreal
            diss(cv_mixt_zmom,ic)  = 0._rfreal
            diss(cv_mixt_ener,ic)  = 0._rfreal
          ENDDO
        ENDIF

! ----- initialize dissipation (later stages)

        IF (istage>1 .AND. ldiss(istage)/=0) THEN
          blend1 = 1._rfreal - betrk(istage)
          DO ic=ibc,iectot
            diss(cv_mixt_dens,ic) = blend1*diss(cv_mixt_dens,ic)
            diss(cv_mixt_xmom,ic) = blend1*diss(cv_mixt_xmom,ic)
            diss(cv_mixt_ymom,ic) = blend1*diss(cv_mixt_ymom,ic)
            diss(cv_mixt_zmom,ic) = blend1*diss(cv_mixt_zmom,ic)
            diss(cv_mixt_ener,ic) = blend1*diss(cv_mixt_ener,ic)
          ENDDO
        ENDIF

#ifdef RFLU
! ----- Compute cell gradients for higher-order scheme

        IF ( regions(ireg)%mixtInput%spaceOrder > 1 ) THEN
          pregion => regions(ireg)
          CALL rflu_convertcvcons2prim(pregion,cv_mixt_state_duvwp)
          CALL rflu_computegradcells(pregion,cv_mixt_dens,cv_mixt_pres, &
                                     grc_mixt_dens,grc_mixt_pres, &
                                     pregion%mixt%cv,pregion%mixt%gradCell)
          CALL rflu_computegradcellseno(pregion,grc_mixt_dens,grc_mixt_pres, &
                                        pregion%mixt%gradCell)
          CALL rflu_limitgradcellssimple(pregion,cv_mixt_dens,cv_mixt_pres, &
                                         grc_mixt_dens,grc_mixt_pres, &
                                         pregion%mixt%cv,pregion%mixt%cvInfo, &
                                         pregion%mixt%gradCell)
          CALL rflu_convertcvprim2cons(pregion,cv_mixt_state_cons)
        END IF ! regions
#endif

#ifdef TURB
          IF (flowmodel == flow_navst .AND. &
              regions(ireg)%mixtInput%turbModel /= turb_model_none) &
            CALL turb_emsinit( regions(ireg),istage )
#endif

! ----- compute numerical dissipation

        IF (ldiss(istage) /= 0) THEN
          CALL numericaldissipation( regions(ireg) )
        ENDIF

#ifdef TURB
        IF (ldiss(istage) /= 0 .AND. flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_ransnumericaldissipation( regions(ireg) )
        ENDIF
#endif

! ----- compute viscous fluxes

        IF (flowmodel==flow_navst .AND. ldiss(istage)/=0) THEN
          CALL viscousfluxes( regions(ireg) )
        ENDIF

! ----- compute convective fluxes; form residual

        CALL convectivefluxes( regions(ireg) )

#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_ransconvectivefluxes( regions(ireg) )
        ENDIF
#endif

! ----- add source terms

        CALL sourceterms( regions(ireg) )

#ifdef INRT
        IF (global%inrtUsed) THEN
          CALL inrt_sourceterms( regions(ireg) )
        ENDIF
#endif
#ifdef RADI
        CALL radi_sourceterms( regions(ireg) )
#endif
#ifdef PERI
        CALL peri_sourceterms( regions(ireg) )
#endif
#ifdef TURB
       IF (flowmodel == flow_navst .AND. &
           regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
         CALL turb_ranssourceterms( regions(ireg) )
       ENDIF
#endif

! ----- zero out residuals in dummy cells

#ifdef RFLO
        CALL rflo_zerodummycells( regions(ireg),rhs )
#endif
#ifdef RFLU
        pregion => regions(ireg)
        CALL rflu_zerovirtualcellvars(pregion,rhs)
#endif
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_ranszerodummycells( regions(ireg) )
        ENDIF
#endif

! ----- compute new forcing term

        IF (ftermnew) THEN
          DO ic=ibc,iectot
            fterm(cv_mixt_dens,ic) = fterm(cv_mixt_dens,ic) - &
                                       rhs(cv_mixt_dens,ic)
            fterm(cv_mixt_xmom,ic) = fterm(cv_mixt_xmom,ic) - &
                                       rhs(cv_mixt_xmom,ic)
            fterm(cv_mixt_ymom,ic) = fterm(cv_mixt_ymom,ic) - &
                                       rhs(cv_mixt_ymom,ic)
            fterm(cv_mixt_zmom,ic) = fterm(cv_mixt_zmom,ic) - &
                                       rhs(cv_mixt_zmom,ic)
            fterm(cv_mixt_ener,ic) = fterm(cv_mixt_ener,ic) - &
                                       rhs(cv_mixt_ener,ic)
          ENDDO
        ENDIF

        ftermnew = .false.

! ----- residual (+ forcing term) * time step / volume

        fac = ark(istage)*cfl

        IF (residfterm) THEN
          DO ic=ibc,iec
            adtv = fac*dt(ic)/vol(ic)
            rhs(cv_mixt_dens,ic) = adtv*(rhs(cv_mixt_dens,ic)+ &
                                         fterm(cv_mixt_dens,ic))
            rhs(cv_mixt_xmom,ic) = adtv*(rhs(cv_mixt_xmom,ic)+ &
                                         fterm(cv_mixt_xmom,ic))
            rhs(cv_mixt_ymom,ic) = adtv*(rhs(cv_mixt_ymom,ic)+ &
                                         fterm(cv_mixt_ymom,ic))
            rhs(cv_mixt_zmom,ic) = adtv*(rhs(cv_mixt_zmom,ic)+ &
                                         fterm(cv_mixt_zmom,ic))
            rhs(cv_mixt_ener,ic) = adtv*(rhs(cv_mixt_ener,ic)+ &
                                         fterm(cv_mixt_ener,ic))
          ENDDO
        ELSE
          DO ic=ibc,iec
            adtv = fac*dt(ic)/vol(ic)
            rhs(cv_mixt_dens,ic) = adtv*rhs(cv_mixt_dens,ic)
            rhs(cv_mixt_xmom,ic) = adtv*rhs(cv_mixt_xmom,ic)
            rhs(cv_mixt_ymom,ic) = adtv*rhs(cv_mixt_ymom,ic)
            rhs(cv_mixt_zmom,ic) = adtv*rhs(cv_mixt_zmom,ic)
            rhs(cv_mixt_ener,ic) = adtv*rhs(cv_mixt_ener,ic)
          ENDDO
        ENDIF

#ifdef RFLO
! ----- implicit residual smoothing

        IF (smoocf > 0._rfreal) THEN
          CALL rflo_zerodummycells( regions(ireg),rhs )
          CALL rflo_residualsmoothing( regions(ireg) )
        ENDIF
#ifdef PEUL
        IF (global%peulUsed .AND. &
            regions(ireg)%peulInput%smoocf > 0._rfreal) THEN
          CALL rflo_zerodummycells( regions(ireg), &
            regions(ireg)%levels(ilev)%peul%rhs )
          CALL peul_residualsmoothingcoeffs( regions(ireg) )
        ENDIF
#endif
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none .AND. &
            regions(ireg)%turbInput%smoocf > 0._rfreal) THEN
          CALL turb_ranszerodummycells( regions(ireg) )
          CALL turb_rflo_ransressmoothing( regions(ireg) )
        ENDIF
#endif
#endif

! ----- update solution

        IF (global%solverType == solv_implicit) THEN
          fac = 1.5_rfreal*ark(istage)*cfl/global%dtMin
          DO ic=ibc,iec
            adtv = 1._rfreal/(1._rfreal+fac*dt(ic))
            cv(cv_mixt_dens,ic) = cvold(cv_mixt_dens,ic) - &
                                  adtv*rhs(cv_mixt_dens,ic)
            cv(cv_mixt_xmom,ic) = cvold(cv_mixt_xmom,ic) - &
                                  adtv*rhs(cv_mixt_xmom,ic)
            cv(cv_mixt_ymom,ic) = cvold(cv_mixt_ymom,ic) - &
                                  adtv*rhs(cv_mixt_ymom,ic)
            cv(cv_mixt_zmom,ic) = cvold(cv_mixt_zmom,ic) - &
                                  adtv*rhs(cv_mixt_zmom,ic)
            cv(cv_mixt_ener,ic) = cvold(cv_mixt_ener,ic) - &
                                  adtv*rhs(cv_mixt_ener,ic)
          ENDDO
        ELSE
          DO ic=ibc,iec
            cv(cv_mixt_dens,ic) = cvold(cv_mixt_dens,ic) - rhs(cv_mixt_dens,ic)
            cv(cv_mixt_xmom,ic) = cvold(cv_mixt_xmom,ic) - rhs(cv_mixt_xmom,ic)
            cv(cv_mixt_ymom,ic) = cvold(cv_mixt_ymom,ic) - rhs(cv_mixt_ymom,ic)
            cv(cv_mixt_zmom,ic) = cvold(cv_mixt_zmom,ic) - rhs(cv_mixt_zmom,ic)
            cv(cv_mixt_ener,ic) = cvold(cv_mixt_ener,ic) - rhs(cv_mixt_ener,ic)
          ENDDO
        ENDIF

#ifdef PERI
        IF (regions(ireg)%periInput%flowKind /= off) THEN
           CALL peri_solutionupdate( regions(ireg) )
        ENDIF
#endif
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_solutionupdate( regions(ireg),istage,ibc,iec )
        ENDIF
#endif

! ----- check for positivity/validity -----------------------------------------

#ifdef RFLO
        CALL rflo_checkvalidity( regions(ireg) )
#endif

#ifdef RFLU
        pregion => regions(ireg)
        CALL rflu_checkvaliditywrapper(pregion)
        CALL rflu_checkpositivitywrapper(pregion)
#endif

#ifdef RFLO
! ----- send conservative variables to other processors

        CALL rflo_boundaryconditionssend( regions,ireg )

#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_rflo_ransbndconditionssend( regions,ireg )
        ENDIF
#endif
#endif
#ifdef RFLU
! ----- exchange conservative variables with other processors

        pregion => regions(ireg)
        CALL rflu_mpi_isendwrapper(pregion)
        CALL rflu_setvars(pregion,1,pregion%grid%nCells)
#endif

! ----- update dependent variables

#ifdef RFLO
        IF (gasmodel == gas_model_tcperf) THEN   ! cp, Mol=const.
          CALL mixtureproperties( regions(ireg),ibc,iec,.false. )
        ELSE
          CALL mixtureproperties( regions(ireg),ibc,iec,.true. )
        ENDIF
#endif

#ifdef RFLO
      ENDIF  ! region on this processor and active
#endif
    ENDDO    ! iReg

#ifdef RFLO
! - facilitate global non-dummy communication ---------------------------------

    CALL globalcommunicationmp( regions )

#ifdef GENX
! - send fluids density at interface

    CALL com_call_function( global%genxHandleBc,2,alpha,1 )

    DO ireg=1,global%nRegions
      IF (regions(ireg)%procid==global%myProcid .AND. &   ! region active and
          regions(ireg)%active==active) THEN              ! on my processor
        IF (regions(ireg)%mixtInput%externalBc) THEN
          CALL rflo_sendboundaryvaluesalpha( regions(ireg) )
        ENDIF
      ENDIF
    ENDDO

! - get BC values at the interface; set BC values in dummy cells

    CALL com_call_function( global%genxHandleBc,2,alpha,2 )
#endif

    DO ireg=1,global%nRegions
      IF (regions(ireg)%procid==global%myProcid .AND. &   ! region active and
          regions(ireg)%active==active) THEN              ! on my processor
        IF (regions(ireg)%mixtInput%externalBc) THEN
          CALL rflo_getboundaryvalues( regions(ireg) )
        ENDIF
        IF (global%solverType == solv_implicit) THEN
          CALL updatetbc( regions(ireg),time,global%dtMin,.true. )
        ENDIF
        CALL rflo_boundaryconditionsset( regions,ireg )
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_rflo_ransbndconditionsset( regions,ireg )
        ENDIF
#endif
      ENDIF
    ENDDO

! - receive variables from other processors; send BC values to external code

    DO ireg=1,global%nRegions
      IF (regions(ireg)%procid==global%myProcid .AND. &   ! region active and
          regions(ireg)%active==active) THEN              ! on my processor
        CALL rflo_boundaryconditionsrecv( regions,ireg )
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_rflo_ransbndconditionsrecv( regions,ireg )
        ENDIF
#endif
      ENDIF
    ENDDO

! - clear send requests

    DO ireg=1,global%nRegions
      IF (regions(ireg)%procid==global%myProcid .AND. &   ! region active and
          regions(ireg)%active==active) THEN              ! on my processor
        CALL rflo_clearsendrequests( regions,ireg,.false. )
#ifdef TURB
        IF (flowmodel == flow_navst .AND. &
            regions(ireg)%mixtInput%turbModel /= turb_model_none) THEN
          CALL turb_rflo_ransclearsendrequests( regions,ireg )
        ENDIF
#endif
      ENDIF
    ENDDO

#ifdef PERI
!    DO iReg=1,global%nRegions
!      IF (regions(iReg)%procid==global%myProcid .AND. &   ! region active and
!          regions(iReg)%active==ACTIVE) THEN              ! on my processor
!        IF (regions(iReg)%periInput%flowKind /= OFF) THEN
!          CALL PERI_CoMeanCorrection( regions(iReg) ) 
!        ENDIF
!      ENDIF
!    ENDDO
#endif

#endif



#ifdef RFLU
    CALL rflu_mpi_copywrapper(regions)

    DO ireg = 1,global%nRegionsLocal
      pregion => regions(ireg)

      CALL rflu_mpi_recvwrapper(pregion)    
      CALL rflu_setvars(pregion,pregion%grid%nCells+1,pregion%grid%nCellsTot) 
    END DO ! iReg

    DO ireg = 1,global%nRegionsLocal
      pregion => regions(ireg)

      CALL rflu_mpi_clearrequestwrapper(pregion)
    END DO ! iReg
#endif


  ENDDO   ! istage

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

  CALL deregisterfunction( global )

END SUBROUTINE explicitmultistage

!******************************************************************************
!
! RCS Revision history:
!
! $Log: ExplicitMultistage.F90,v $
! Revision 1.15  2008/12/06 08:44:09  mtcampbe
! Updated license.
!
! Revision 1.14  2008/11/19 22:17:22  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.13  2006/03/26 20:21:14  haselbac
! Changed to wrappers bcos of GL model
!
! Revision 1.12  2005/12/03 19:45:05  haselbac
! Apparent bug fix: Separated call to RFLU_MPI_ClearRequestWrapper into separate loop
!
! Revision 1.11  2005/10/31 21:09:34  haselbac
! Changed specModel and SPEC_MODEL_NONE
!
! Revision 1.10  2005/10/05 13:47:23  haselbac
! Adapted to new module for cell grads
!
! Revision 1.9  2005/07/13 17:26:20  haselbac
! Bug fix: Adapted to changes to RFLU_LimitGradCells
!
! Revision 1.8  2005/05/16 20:38:42  haselbac
! Renamed RFLU_ZeroDummyCells, moved computation of dt into RFLU_TimeStepping
!
! Revision 1.7  2005/04/29 00:06:09  haselbac
! Added routines to clear send requests
!
! Revision 1.6  2005/04/15 15:06:01  haselbac
! Converted to MPI
!
! Revision 1.5  2005/04/06 01:59:42  wasistho
! mv call to PERI_CoMeanCorrection to after BC treatment, commented for now
!
! Revision 1.4  2005/03/10 02:08:25  wasistho
! commented PERI_coMeanCorrection temporarily for testing
!
! Revision 1.3  2005/03/07 05:05:18  wasistho
! install hybrid DESSA turbulence model
!
! Revision 1.2  2004/12/28 22:49:15  wasistho
! moved RFLO_Bcond* and RFLO_BoundaryCond* routines into RFLO_ModBoundaryConditions
!
! Revision 1.1  2004/12/01 16:48:30  haselbac
! Initial revision after changing case
!
! Revision 1.68  2004/11/17 23:44:33  wasistho
! used generic RK-update for rocturb
!
! Revision 1.67  2004/11/14 19:34:26  haselbac
! Replaced call to mixtureProperties by RFLU_SetVars
!
! Revision 1.66  2004/10/19 19:25:34  haselbac
! Adapted calls to time step routines
!
! Revision 1.65  2004/07/26 19:08:59  wasistho
! add RFLO_CheckValidity
!
! Revision 1.64  2004/03/27 03:02:15  wasistho
! added ifdef RFLO within ifdef TURB
!
! Revision 1.63  2004/03/20 00:26:35  wasistho
! set turb_rflo_ransNumericalDiss to turb_ransNumerical..
!
! Revision 1.62  2004/03/19 02:39:55  wasistho
! renamed TURB_RFLO_RansZeroDummyCells to TURB_RansZeroDummyCells
!
! Revision 1.61  2004/03/11 03:32:13  wasistho
! changed rocturb nomenclature
!
! Revision 1.60  2004/03/05 22:09:00  jferry
! created global variables for peul, plag, and inrt use
!
! Revision 1.59  2004/03/02 21:49:21  jferry
! Added inrtUsed flag to mixture data structure
!
! Revision 1.58  2004/02/26 21:13:55  wasistho
! changed TURB_ransEmsInit to TURB_emsInit
!
! Revision 1.57  2004/01/29 22:52:41  haselbac
! Added info argument to cell-gradient routine call
!
! Revision 1.56  2003/12/04 03:22:59  haselbac
! Added second-order scheme, viscous fluxes, validity check
!
! Revision 1.55  2003/11/25 21:01:39  haselbac
! Added routine to update dummy cells
!
! Revision 1.54  2003/11/20 16:40:35  mdbrandy
! Backing out RocfluidMP changes from 11-17-03
!
! Revision 1.51  2003/10/27 04:49:29  wasistho
! replace ransCentralDiss by ransNumericalDiss.
!
! Revision 1.50  2003/10/21 03:57:42  wasistho
! put turb and peul spectralradii before its smoothingcoef
!
! Revision 1.49  2003/10/16 20:12:27  wasistho
! completed incorporation of RaNS/DES
!
! Revision 1.48  2003/10/03 20:12:11  wasistho
! initial installation of turbModel SA and DES
!
! Revision 1.47  2003/09/26 21:45:34  fnajjar
! Modified ModInterfaces calls to proper physical modules
!
! Revision 1.46  2003/08/28 20:32:48  wasistho
! excluced ModInterfacesTurbulence,Radiation,Periodic from ModInterfaces
!
! Revision 1.45  2003/08/13 02:19:42  wasistho
! added call to radiation source term routine
!
! Revision 1.44  2003/07/22 01:52:40  haselbac
! Bug fix: dual-time stepping only for RFLO, lead to core dump
!
! Revision 1.43  2003/07/03 21:48:44  jblazek
! Implemented dual-time stepping.
!
! Revision 1.42  2003/05/15 02:57:02  jblazek
! Inlined index function.
!
! Revision 1.41  2003/05/09 17:01:03  jiao
! Renamed the COM_call_function_handlers to COM_call_function.
!
! Revision 1.40  2003/04/05 02:01:08  wasistho
! regions to region in PERI_solutionUpdate
!
! Revision 1.39  2003/03/29 03:26:21  wasistho
! install ROCPERI
!
! Revision 1.38  2003/03/15 16:17:58  haselbac
! Added zeroing of dummies (?) and changed FEM to IDXL
!
! Revision 1.37  2003/03/05 20:41:21  jiao
! ACH: Split update inbuff calls into get correct dependency of rhof on mdot
!
! Revision 1.36  2003/03/04 22:12:34  jferry
! Initial import of Rocinteract
!
! Revision 1.35  2003/02/11 22:52:50  jferry
! Initial import of Rocsmoke
!
! Revision 1.34  2002/12/06 22:29:25  jblazek
! Corrected bug for geometry exchange between minimal patches.
!
! Revision 1.33  2002/10/27 18:47:32  haselbac
! Added call to RFLU_CheckPositivity
!
! Revision 1.32  2002/10/17 06:50:42  jiao
! Changed 0. to 0._RFREAL.
!
! Revision 1.31  2002/10/12 20:21:25  jblazek
! Rearranged UpdateTbc (after getting values from GenX);
! corrected bug in externalBc.
!
! Revision 1.30  2002/10/05 18:34:08  haselbac
! Moved ifdef CHARM inside ifdef RFLU
!
! Revision 1.29  2002/09/20 22:22:35  jblazek
! Finalized integration into GenX.
!
! Revision 1.27  2002/09/09 13:57:49  haselbac
! added viscous routines, bug fix for iec, mixtInput under regions
!
! Revision 1.26  2002/09/05 17:40:20  jblazek
! Variable global moved into regions().
!
! Revision 1.25  2002/08/29 21:53:18  jblazek
! Added support for moving grids.
!
! Revision 1.24  2002/08/16 21:33:47  jblazek
! Changed interface to MixtureProperties.
!
! Revision 1.23  2002/08/15 19:48:05  jblazek
! Implemented grid deformation capability.
!
! Revision 1.22  2002/07/25 14:50:52  haselbac
! Added RFLU_ModFEM and call for update of dummy cells
!
! Revision 1.21  2002/07/23 20:28:19  wasistho
! #ifdef RFLO around ViscousFluxes removed
!
! Revision 1.20  2002/07/22 17:01:23  jblazek
! Removed MPI_Barrier at the end of stage loop.
!
! Revision 1.19  2002/07/05 23:20:46  jblazek
! Corrected bug in perfgasDependentVars.F90; did some cosmetics.
!
! Revision 1.18  2002/06/27 15:46:27  haselbac
! Added FEM calls for communication - not activated yet
!
! Revision 1.17  2002/06/14 20:11:48  haselbac
! Deleted ModLocal, renamed local%nRegions to global%nRegionsLocal
!
! Revision 1.16  2002/06/07 16:40:36  jblazek
! Grid & solution for all regions in one file.
!
! Revision 1.15  2002/05/28 14:14:39  haselbac
! Enclosed viscous fluxes within RFLO conditional statement
!
! Revision 1.14  2002/05/28 13:44:05  haselbac
! Cosmetic changes only?
!
! Revision 1.13  2002/05/21 01:51:56  wasistho
! add viscous terms
!
! Revision 1.12  2002/05/04 16:31:49  haselbac
! Added RFLU statements
!
! Revision 1.11  2002/04/12 17:36:23  jblazek
! Added timer.
!
! Revision 1.10  2002/04/01 19:36:07  jblazek
! Added routine to clear send requests.
!
! Revision 1.9  2002/03/30 00:50:49  jblazek
! Cleaned up with flint.
!
! Revision 1.8  2002/03/29 23:15:22  jblazek
! Corrected bug in MPI send.
!
! Revision 1.7  2002/03/18 22:25:45  jblazek
! Finished multiblock and MPI.
!
! Revision 1.6  2002/02/21 23:25:05  jblazek
! Blocks renamed as regions.
!
! Revision 1.5  2002/02/09 01:47:01  jblazek
! Added multi-probe option, residual smoothing, physical time step.
!
! Revision 1.4  2002/01/28 23:55:22  jblazek
! Added flux computation (central scheme).
!
! Revision 1.3  2002/01/23 23:37:32  jblazek
! All blocks passed to time integration routines.
!
! Revision 1.2  2002/01/23 03:51:24  jblazek
! Added low-level time-stepping routines.
!
! Revision 1.1  2002/01/16 22:03:34  jblazek
! Added time-stepping routines.
!
!******************************************************************************