RFLU_ModPartitionRegionUtils.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: Suite of utility routines to partition region.
!
! Description: None.
!
! Notes: None.
!
! ******************************************************************************
!
! $Id: RFLU_ModPartitionRegionUtils.F90,v 1.8 2008/12/06 08:45:03 mtcampbe Exp $
!
! Copyright: (c) 2005-2007 by the University of Illinois
!
! ******************************************************************************

MODULE rflu_modpartitionregionutils

  USE modglobal, ONLY: t_global
  USE moddatatypes
  USE modparameters
  USE moderror
  USE modbndpatch, ONLY: t_patch
  USE moddatastruct, ONLY: t_level,t_region
  USE modgrid, ONLY: t_grid
  USE modmpi

  USE modsortsearch

  IMPLICIT NONE
  
  PRIVATE
  PUBLIC :: rflu_part_addvirtualcellsinv1, & 
            rflu_part_addvirtualcellsinv2

! ******************************************************************************
! Declarations and definitions
! ******************************************************************************

  CHARACTER(CHRLEN) :: & 
    RCSIdentString = '$RCSfile: RFLU_ModPartitionRegionUtils.F90,v $ $Revision: 1.8 $'


! ******************************************************************************
! Routines
! ******************************************************************************

  CONTAINS






! ******************************************************************************
!
! Purpose: Add virtual cells for inviscid fluxes based on first-order scheme.
!
! Description: None.
!
! Input:
!   pRegion             Pointer to region
!   pRegionSerial       Pointer to serial region        
!   vc                  List of virtual cells (empty)
!   nCellsVirtMax       Maximum allowable number of virtual cells  
!
! Output: 
!   vc                  List of virtual cells 
!   nCellsVirt          Number of virtual cells
!
! Notes: None.
!
! ******************************************************************************

  SUBROUTINE rflu_part_addvirtualcellsinv1(pRegion,pRegionSerial,vc, &
                                           ncellsvirtmax,ncellsvirt)                                   

    IMPLICIT NONE
      
! ******************************************************************************
!   Declarations and definitions
! ******************************************************************************

! ==============================================================================
!   Arguments
! ==============================================================================

    INTEGER, INTENT(IN) :: ncellsvirtmax
    INTEGER, INTENT(OUT) :: ncellsvirt
    INTEGER, INTENT(INOUT) :: vc(ncellsvirtmax)    
    TYPE(t_region), POINTER :: pregion,pregionserial

! ==============================================================================
!   Locals
! ==============================================================================

    INTEGER :: c1,c2,errorflag,icg,icl,ifg,ifl,iflbeg,iflend,iloc
    TYPE(t_grid), POINTER :: pgrid,pgridserial
    TYPE(t_global), POINTER :: global     
      
! ******************************************************************************
!   Start
! ******************************************************************************

    global => pregionserial%global

    CALL registerfunction(global,'RFLU_PART_AddVirtualCellsInv1',&
  'RFLU_ModPartitionRegionUtils.F90')

    IF ( global%verbLevel > verbose_none ) THEN
      WRITE(stdout,'(A,1X,A)') solver_name, & 
            'Adding virtual cells for inviscid first-order stencil...' 
    END IF ! global%verbLevel

! ******************************************************************************
!   Set pointers and initialize variables
! ******************************************************************************

    pgrid       => pregion%grid
    pgridserial => pregionserial%grid    

    ncellsvirt = 0

    iflbeg = pgridserial%avfCSRInfo(pregion%iRegionGlobal)

    IF ( pregion%iRegionGlobal /= global%nRegionsLocal ) THEN
      iflend = pgridserial%avfCSRInfo(pregion%iRegionGlobal+1)-1
    ELSE 
      iflend = 2*pgridserial%nFacesCut
    END IF ! pRegion%iRegionGlobal
   
! ******************************************************************************
!  Loop over actual-virtual faces
! ******************************************************************************
    
    DO ifl = iflbeg,iflend        
      ifg = pgridserial%avfCSR(ifl)

      DO icl = 1,2
        icg = pgridserial%f2c(icl,ifg)

        IF ( pgridserial%sc2r(icg) /= pregion%iRegionGlobal ) THEN
          IF ( ncellsvirt > 0 ) THEN 
            CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg,iloc)
          ELSE 
            iloc = element_not_found
          END IF ! nCellsVirt

          IF ( iloc == element_not_found ) THEN 
            ncellsvirt = ncellsvirt + 1

            IF ( ncellsvirt <= ncellsvirtmax ) THEN 
              vc(ncellsvirt) = icg
            ELSE 
              CALL errorstop(global,err_exceed_dimens,__line__,'vc')
            END IF ! nCellsVirt
              
            IF ( ncellsvirt > 1 ) THEN    
              CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt)              
            END IF ! nCellsVirt 
          END IF ! iLoc
        END IF ! pGridSerial%sc2r 
      END DO ! icl     
    END DO ! ifl 
                          
! ******************************************************************************
!   End
! ******************************************************************************

    IF ( global%myProcid == masterproc .AND. &
         global%verbLevel > verbose_none ) THEN   
      WRITE(stdout,'(A,1X,A)') solver_name, & 
            'Adding virtual cells for inviscid first-order stencil done.' 
    END IF ! global%myProcid

    CALL deregisterfunction(global)

  END SUBROUTINE rflu_part_addvirtualcellsinv1








! ******************************************************************************
!
! Purpose: Add virtual cells for inviscid fluxes based on higher-order scheme.
!
! Description: Building list of virtual cells proceeds in several steps:
!   1. Build list of vertices from list of actual-virtual faces
!   2. Build list of source cells adjacent to actual-virtual faces which are 
!      in same region as that for which virtual cells are constructed. One 
!      layer of source cells is sufficient. If more layers of source cells 
!      are specified here, this should not change the number of virtual cells.
!   3. Build stencils for source cells and add stencil members to list of 
!      virtual cells if not in same region
!   4. Loop over existing virtual cells and add layers.
!
! Input:
!   pRegion             Pointer to region
!   pRegionSerial       Pointer to serial region        
!   vc                  List of virtual cells (empty)
!   nCellsVirtMax       Maximum allowable number of virtual cells  
!
! Output: 
!   vc                  List of virtual cells 
!   nCellsVirt          Number of virtual cells
!
! Notes: None.
!
! ******************************************************************************

  SUBROUTINE rflu_part_addvirtualcellsinv2(pRegion,pRegionSerial,vc, &
                                           ncellsvirtmax,ncellsvirt)                                   

    USE rflu_modstencilscells
    USE rflu_modtopologyutils, ONLY: rflu_buildfacevertlist, & 
                                     rflu_buildvertcellnghblist

    IMPLICIT NONE
      
! ******************************************************************************
!   Declarations and definitions
! ******************************************************************************

! ==============================================================================
!   Arguments
! ==============================================================================

    INTEGER, INTENT(IN) :: ncellsvirtmax
    INTEGER, INTENT(OUT) :: ncellsvirt
    INTEGER, INTENT(INOUT) :: vc(ncellsvirtmax)    
    TYPE(t_region), POINTER :: pregion,pregionserial

! ==============================================================================
!   Locals
! ==============================================================================

    INTEGER, PARAMETER :: loop_counter_limit = 5
    INTEGER :: errorflag,i,icg,icg2,idimens,idir,iflbeg,iflend,ilayer,iloc,j, &
              loopcounter,nlayers,nvert,nvertest,scdim,vcnewdim,vcnewdimmax, &
              vcolddim
    INTEGER, DIMENSION(:), ALLOCATABLE :: avv,sc,vcnew,vcold
    TYPE(t_grid), POINTER :: pgrid,pgridserial
    TYPE(t_global), POINTER :: global     
      
! ******************************************************************************
!   Start
! ******************************************************************************

    global => pregionserial%global

    CALL registerfunction(global,'RFLU_PART_AddVirtualCellsInv2',&
  'RFLU_ModPartitionRegionUtils.F90')

    IF ( global%myProcid == masterproc .AND. & 
         global%verbLevel > verbose_none ) THEN
      WRITE(stdout,'(A,1X,A)') solver_name, & 
            'Adding virtual cells for inviscid higher-order stencil...' 
    END IF ! global%myProcid

! ******************************************************************************
!   Set pointers and initialize variables
! ******************************************************************************

    pgrid       => pregion%grid
    pgridserial => pregionserial%grid    

    ncellsvirt = 0

    loopcounter = 0

! ******************************************************************************
!   Build list of interface vertices
! ******************************************************************************

! ==============================================================================
!   Estimate number of interface vertices 
! ==============================================================================

    iflbeg = pgridserial%avfCSRInfo(pregion%iRegionGlobal)

    IF ( pregion%iRegionGlobal /= global%nRegionsLocal ) THEN
      iflend = pgridserial%avfCSRInfo(pregion%iRegionGlobal+1)-1
    ELSE 
      iflend = 2*pgridserial%nFacesCut
    END IF ! pRegion%iRegionGlobal
      
    nvertest = 2*(iflend - iflbeg + 1)  

    IF ( nvertest < 100 ) THEN ! Kludge
      nvertest = 100 + 4*nvertest
    ELSE IF ( nvertest < 1000 ) THEN 
      nvertest = 2*nvertest
    END IF ! nVertEst      

! ==============================================================================
!   Build list of interface vertices. NOTE this sometimes failed for some cases
!   because the estimated size was too small, so added capability of adaptively 
!   increasing size within some limits to prevent code from aborting.
! ==============================================================================

    emptyloop: DO 
      loopcounter = loopcounter + 1

      ALLOCATE(avv(nvertest),stat=errorflag)
      global%error = errorflag   
      IF ( global%error /= err_none ) THEN 
        CALL errorstop(global,err_allocate,__line__,'avv')
      END IF ! global%error      
      
      CALL rflu_buildfacevertlist(global,pgridserial, & 
                                  pgridserial%avfCSR(iflbeg:iflend), & 
                                  iflend-iflbeg+1,avv,nvertest,nvert,errorflag)
     
      IF ( errorflag /= err_none ) THEN 
        IF ( loopcounter <= loop_counter_limit ) THEN 
          DEALLOCATE(avv,stat=errorflag)
          global%error = errorflag   
          IF ( global%error /= err_none ) THEN 
            CALL errorstop(global,err_allocate,__line__,'avv')
          END IF ! global%error      
  
          nvertest = 2*nvertest

          global%warnCounter = global%warnCounter + 1

          IF ( global%myProcid == masterproc .AND. &
               global%verbLevel > verbose_low ) THEN
            WRITE(stdout,'(A,3X,A,1X,I2,1X,A)') solver_name, & 
                         '*** WARNING *** Attempt ',loopcounter, & 
                         'to build vertex list failed because array too small.' 
            WRITE(stdout,'(A,19X,A,1X,I7)') solver_name, & 
                         'Attempting again with array size:',nvertest
          END IF ! global%myProcid
        ELSE
          CALL errorstop(global,err_allocate_adaptive,__line__) 
        END IF ! loopCounter
      ELSE 
        EXIT emptyloop
      END IF ! errorFlag
    END DO emptyloop

! ******************************************************************************
!   Based on list of interface vertices, build list of source cells adjacent to 
!   actual-virtual faces for construction of virtual cells. One layer of source 
!   cells is sufficient. If more layers of source cells are specified here, this
!   should not change the number of virtual cells.
! ******************************************************************************
                               
    ALLOCATE(sc(ncellsvirtmax),stat=errorflag)
    global%error = errorflag   
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_allocate,__line__,'sc')
    END IF ! global%error             

    nlayers = 1    

    CALL rflu_buildvertcellnghblist(global,pgridserial,avv,nvert,nlayers, &
                                    pregion%iRegionGlobal,sc,ncellsvirtmax, &
                                    scdim)

    DEALLOCATE(avv,stat=errorflag)
    global%error = errorflag   
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_deallocate,__line__,'avv')
    END IF ! global%error 
                                
! ******************************************************************************
!   Loop over list of source cells, build stencil for each cell, and add cells
!   in stencil as virtual cells if not in same region and not already added
! ******************************************************************************

    CALL rflu_setinfoc2cstencilwrapper(pregionserial, & 
                                       pregionserial%mixtInput%spaceOrder-1)        
    CALL rflu_createc2cstencilwrapper(pregionserial)   

    DO i = 1,scdim
      icg = sc(i)
      
      CALL rflu_buildc2cstencilwrapper(pregionserial,icg,constr_none)
      
      SELECT CASE ( pregionserial%mixtInput%stencilDimensCells ) 
        CASE ( 1 ) 
          DO idimens = 1,pregionserial%mixtInput%dimens
            idir = xcoord - 1 + idimens 

            DO j = 1,pgridserial%c2cs1D(idir,icg)%nCellMembs
              icg2 = pgridserial%c2cs1D(idir,icg)%cellMembs(j)

              IF ( pgridserial%sc2r(icg2) /= pregion%iRegionGlobal ) THEN
                IF ( ncellsvirt > 0 ) THEN
                  CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)
                ELSE
                  iloc = element_not_found
                END IF ! nCellsVirt

                IF ( iloc == element_not_found ) THEN
                  ncellsvirt = ncellsvirt + 1
  
                  vc(ncellsvirt) = icg2
  
                  IF ( ncellsvirt > 1 ) THEN
                    CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt)
                  END IF ! nCellsVirt 
                END IF ! iLoc
              END IF ! pGridSerial%sc2r                      
            END DO ! j
          END DO ! iDimens
        CASE ( 2,3 ) 
          DO j = 1,pgridserial%c2cs(icg)%nCellMembs
            icg2 = pgridserial%c2cs(icg)%cellMembs(j)
                                
            IF ( pgridserial%sc2r(icg2) /= pregion%iRegionGlobal ) THEN
              IF ( ncellsvirt > 0 ) THEN 
                CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)
              ELSE 
                iloc = element_not_found
              END IF ! nCellsVirt

              IF ( iloc == element_not_found ) THEN 
                ncellsvirt = ncellsvirt + 1
              
                vc(ncellsvirt) = icg2
                  
                IF ( ncellsvirt > 1 ) THEN    
                  CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt)              
                END IF ! nCellsVirt 
              END IF ! iLoc
            END IF ! pGridSerial%sc2r                      
          END DO ! j      
        CASE default 
          CALL errorstop(global,err_reached_default,__line__)
      END SELECT ! pRegionSerial%mixtInput%stencilDimensCells
    END DO ! i      
          
    DEALLOCATE(sc,stat=errorflag)
    global%error = errorflag   
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_deallocate,__line__,'sc')
    END IF ! global%error             

! ******************************************************************************
!   Loop over current list of virtual cells and add layers of cells. For the 
!   current scheme, two additional layers are required.
! ******************************************************************************

! ==============================================================================
!   Allocate temporary memory 
! ==============================================================================

    vcolddim    = ncellsvirt
    vcnewdimmax = ncellsvirtmax

    ALLOCATE(vcold(vcolddim),stat=errorflag)
    global%error = errorflag
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_allocate,__line__,'vcOld')
    END IF ! global%error

    DO i = 1,vcolddim
      vcold(i) = vc(i)
    END DO ! i

    ALLOCATE(vcnew(vcnewdimmax),stat=errorflag)
    global%error = errorflag
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_allocate,__line__,'vcNew')
    END IF ! global%error

! ==============================================================================
!   Loop over layers
! ==============================================================================

    nlayers = 2

    DO ilayer = 1,nlayers
      vcnewdim = 0            
         
      DO i = 1,vcolddim
        icg = vcold(i)

        CALL rflu_buildc2cstencilwrapper(pregionserial,icg,constr_none)
        
        SELECT CASE ( pregionserial%mixtInput%stencilDimensCells ) 
          CASE ( 1 ) 
            DO idimens = 1,pregionserial%mixtInput%dimens
              idir = xcoord - 1 + idimens 
  
              DO j = 1,pgridserial%c2cs1D(idir,icg)%nCellMembs
                icg2 = pgridserial%c2cs1D(idir,icg)%cellMembs(j)
                
                IF ( pgridserial%sc2r(icg2) /= pregion%iRegionGlobal ) THEN
                  CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)

                  IF ( iloc == element_not_found ) THEN
                    IF ( ncellsvirt < ncellsvirtmax ) THEN  
                      ncellsvirt = ncellsvirt + 1              
                      vc(ncellsvirt) = icg2
                    ELSE 
                      CALL errorstop(global,err_exceed_dimens,__line__,'vc')          
                    END IF ! nCellsVirt

                    IF ( vcnewdim < vcnewdimmax ) THEN 
                      vcnewdim = vcnewdim + 1              
                      vcnew(vcnewdim) = icg2
                    ELSE
                      CALL errorstop(global,err_exceed_dimens,__line__,'vcNew')         
                    END IF ! vcNewDim 

                    CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt) 
                  END IF ! iLoc
                END IF ! pGridSerial%sc2r                
              END DO ! j
            END DO ! iDimens
          CASE ( 2,3 )       
            DO j = 1,pgridserial%c2cs(icg)%nCellMembs
              icg2 = pgridserial%c2cs(icg)%cellMembs(j)

              IF ( pgridserial%sc2r(icg2) /= pregion%iRegionGlobal ) THEN
                CALL binarysearchinteger(vc(1:ncellsvirt),ncellsvirt,icg2,iloc)

                IF ( iloc == element_not_found ) THEN
                  IF ( ncellsvirt < ncellsvirtmax ) THEN  
                    ncellsvirt = ncellsvirt + 1              
                    vc(ncellsvirt) = icg2
                  ELSE 
                    CALL errorstop(global,err_exceed_dimens,__line__,'vc')          
                  END IF ! nCellsVirt

                  IF ( vcnewdim < vcnewdimmax ) THEN 
                    vcnewdim = vcnewdim + 1              
                    vcnew(vcnewdim) = icg2
                  ELSE
                    CALL errorstop(global,err_exceed_dimens,__line__,'vcNew')         
                  END IF ! vcNewDim 

                  CALL quicksortinteger(vc(1:ncellsvirt),ncellsvirt) 
                END IF ! iLoc
              END IF ! pGridSerial%sc2r
            END DO ! j
          CASE default
            CALL errorstop(global,err_reached_default,__line__)
        END SELECT ! pRegionSerial%mixtInput%stencilDimensCells                      
      END DO ! i
             
      DEALLOCATE(vcold,stat=errorflag)
      global%error = errorflag
      IF ( global%error /= err_none ) THEN 
        CALL errorstop(global,err_deallocate,__line__,'vcOld')
      END IF ! global%error
      
      IF ( ilayer /= nlayers ) THEN 
        vcolddim = vcnewdim
      
        ALLOCATE(vcold(vcolddim),stat=errorflag)
        global%error = errorflag
        IF ( global%error /= err_none ) THEN 
          CALL errorstop(global,err_allocate,__line__,'vcOld')
        END IF ! global%error    
      
        DO i = 1,vcolddim
          vcold(i) = vcnew(i)
        END DO ! i
      END IF ! iLayer      
    END DO ! iLayer   

! ==============================================================================
!   Deallocate temporary memory 
! ==============================================================================

    DEALLOCATE(vcnew,stat=errorflag)
    global%error = errorflag
    IF ( global%error /= err_none ) THEN 
      CALL errorstop(global,err_deallocate,__line__,'vcNew')
    END IF ! global%error
           
    CALL rflu_destroyc2cstencilwrapper(pregionserial)          
                          
! ******************************************************************************
!   End
! ******************************************************************************

    IF ( global%myProcid == masterproc .AND. &
         global%verbLevel > verbose_none ) THEN   
      WRITE(stdout,'(A,1X,A)') solver_name, & 
            'Adding virtual cells for inviscid higher-order stencil done.' 
    END IF ! global%myProcid

    CALL deregisterfunction(global)

  END SUBROUTINE rflu_part_addvirtualcellsinv2







! ******************************************************************************
! End
! ******************************************************************************
  
END MODULE rflu_modpartitionregionutils


! ******************************************************************************
!
! RCS Revision history:
!
! $Log: RFLU_ModPartitionRegionUtils.F90,v $
! Revision 1.8  2008/12/06 08:45:03  mtcampbe
! Updated license.
!
! Revision 1.7  2008/11/19 22:18:14  mtcampbe
! Added Illinois Open Source License/Copyright
!
! Revision 1.6  2007/03/06 18:10:09  haselbac
! Adapted to building 1d stencils for higher-order scheme
!
! Revision 1.5  2006/01/06 22:17:20  haselbac
! Adapted to name changes
!
! Revision 1.4  2005/10/27 19:21:26  haselbac
! Adapted to changes in stencil routine names
!
! Revision 1.3  2005/10/05 14:25:11  haselbac
! Adapted to changes in stencil modules
!
! Revision 1.2  2005/06/14 17:48:14  haselbac
! Adapted RFLU_AddVirtualCellsInv2 to adaptive memory allocation
!
! Revision 1.1  2005/04/15 15:09:13  haselbac
! Initial revision
!
! Revision 1.1  2005/01/17 19:47:47  haselbac
! Initial revision
!
! ******************************************************************************