region.F90       coverage:  84.21 %func     48.14 %block


     1) module Region_module
     2)  
     3)   use Geometry_module
     4)   
     5)   use PFLOTRAN_Constants_module
     6) 
     7)   implicit none
     8) 
     9)   private
    10) 
    11) #include "petsc/finclude/petscsys.h"
    12) 
    13)   PetscInt, parameter, public :: STRUCTURED_GRID_REGION = 1
    14)   PetscInt, parameter, public :: UNSTRUCTURED_GRID_REGION = 2
    15) 
    16)   PetscInt, parameter, public :: DEFINED_BY_BLOCK = 1
    17)   PetscInt, parameter, public :: DEFINED_BY_COORD = 2
    18)   PetscInt, parameter, public :: DEFINED_BY_CELL_IDS = 3
    19)   PetscInt, parameter, public :: DEFINED_BY_CELL_AND_FACE_IDS = 4
    20)   PetscInt, parameter, public :: DEFINED_BY_VERTEX_IDS = 5
    21)   PetscInt, parameter, public :: DEFINED_BY_SIDESET_UGRID = 6
    22)   PetscInt, parameter, public :: DEFINED_BY_FACE_UGRID_EXP = 7
    23)   PetscInt, parameter, public :: DEFINED_BY_POLY_BOUNDARY_FACE = 8
    24)   PetscInt, parameter, public :: DEFINED_BY_POLY_CELL_CENTER = 9
    25)   PetscInt, parameter, public :: DEFINED_BY_CARTESIAN_BOUNDARY = 10
    26) 
    27)   type, public :: block_type        
    28)     PetscInt :: i1,i2,j1,j2,k1,k2    
    29)     type(block_type), pointer :: next
    30)   end type block_type
    31)  
    32)   type, public :: region_type
    33)     PetscInt :: id
    34)     PetscInt :: def_type
    35)     PetscBool :: hdf5_ugrid_kludge  !TODO(geh) tear this out!!!!!
    36)     character(len=MAXWORDLENGTH) :: name
    37)     character(len=MAXSTRINGLENGTH) :: filename
    38)     PetscInt :: i1,i2,j1,j2,k1,k2
    39)     type(point3d_type), pointer :: coordinates(:)
    40)     PetscInt :: iface
    41)     PetscInt :: num_cells
    42)     PetscInt, pointer :: cell_ids(:)
    43)     PetscInt, pointer :: faces(:)
    44)     !TODO(geh): Tear anything to do with structured/unstructured grids other
    45)     !           than cell id ane face id out of region.
    46)     PetscInt, pointer :: vertex_ids(:,:) ! For Unstructured mesh
    47)     PetscInt :: num_verts              ! For Unstructured mesh
    48)     PetscInt :: grid_type  ! To identify whether region is applicable to a Structured or Unstructred mesh
    49)     type(region_sideset_type), pointer :: sideset
    50)     type(region_explicit_face_type), pointer :: explicit_faceset
    51)     type(polygonal_volume_type), pointer :: polygonal_volume
    52)     type(region_type), pointer :: next
    53)   end type region_type
    54)   
    55)   type, public :: region_ptr_type
    56)     type(region_type), pointer :: ptr
    57)   end type region_ptr_type
    58)   
    59)   type, public :: region_list_type
    60)     PetscInt :: num_regions
    61)     type(region_type), pointer :: first
    62)     type(region_type), pointer :: last
    63)     type(region_type), pointer :: array(:)
    64)   end type region_list_type
    65) 
    66)   type, public :: region_sideset_type
    67)     PetscInt :: nfaces
    68)     PetscInt, pointer :: face_vertices(:,:)
    69)   end type region_sideset_type
    70)   
    71)   type, public :: region_explicit_face_type
    72)     type(point3d_type), pointer :: face_centroids(:)
    73)     PetscReal, pointer :: face_areas(:)
    74)   end type region_explicit_face_type
    75)   
    76)   interface RegionCreate
    77)     module procedure RegionCreateWithBlock
    78)     module procedure RegionCreateWithList
    79)     module procedure RegionCreateWithNothing
    80)     module procedure RegionCreateWithRegion    
    81)   end interface RegionCreate
    82)   
    83)   interface RegionReadFromFile
    84)     module procedure RegionReadFromFileId
    85)     module procedure RegionReadFromFilename
    86)     module procedure RegionReadSideSet
    87)     module procedure RegionReadExplicitFaceSet
    88)   end interface RegionReadFromFile
    89)   
    90)   public :: RegionCreate, RegionDestroy, RegionAddToList, RegionReadFromFile, &
    91)             RegionInitList, RegionDestroyList, RegionGetPtrFromList, & 
    92)             RegionRead, RegionReadSideSet, RegionCreateSideset, &
    93)             RegionInputRecord
    94)   
    95) contains
    96) 
    97) ! ************************************************************************** !
    98) 
    99) function RegionCreateWithNothing()
   100)   ! 
   101)   ! Creates a region with no arguments
   102)   ! 
   103)   ! Author: Glenn Hammond
   104)   ! Date: 10/23/07
   105)   ! 
   106) 
   107)   implicit none
   108)   
   109)   type(region_type), pointer :: RegionCreateWithNothing
   110)   
   111)   type(region_type), pointer :: region
   112)   
   113)   allocate(region)
   114)   region%id = 0
   115)   region%def_type = 0
   116)   region%hdf5_ugrid_kludge = PETSC_FALSE
   117)   region%name = ""
   118)   region%filename = ""
   119)   region%i1 = 0
   120)   region%i2 = 0
   121)   region%j1 = 0
   122)   region%j2 = 0
   123)   region%k1 = 0
   124)   region%k2 = 0
   125)   region%iface = 0
   126)   region%num_cells = 0
   127)   ! By default it is assumed that the region is applicable to strucutred grid,
   128)   ! unless explicitly stated in pflotran input file
   129)   region%grid_type = STRUCTURED_GRID_REGION
   130)   region%num_verts = 0
   131)   nullify(region%coordinates)
   132)   nullify(region%cell_ids)
   133)   nullify(region%faces)
   134)   nullify(region%vertex_ids)
   135)   nullify(region%sideset)
   136)   nullify(region%explicit_faceset)
   137)   nullify(region%polygonal_volume)
   138)   nullify(region%next)
   139)   
   140)   RegionCreateWithNothing => region
   141) 
   142) end function RegionCreateWithNothing
   143) 
   144) ! ************************************************************************** !
   145) 
   146) function RegionCreateSideset()
   147)   ! 
   148)   ! Creates a sideset
   149)   ! 
   150)   ! Author: Glenn Hammond
   151)   ! Date: 12/19/11
   152)   ! 
   153) 
   154)   implicit none
   155)   
   156)   type(region_sideset_type), pointer :: RegionCreateSideset
   157)   
   158)   type(region_sideset_type), pointer :: sideset
   159)   
   160)   allocate(sideset)
   161)   sideset%nfaces = 0
   162)   nullify(sideset%face_vertices)
   163)   
   164)   RegionCreateSideset => sideset
   165) 
   166) end function RegionCreateSideset
   167) 
   168) ! ************************************************************************** !
   169) 
   170) function RegionCreateExplicitFaceSet()
   171)   ! 
   172)   ! Creates a sideset
   173)   ! 
   174)   ! Author: Glenn Hammond
   175)   ! Date: 12/19/11
   176)   ! 
   177) 
   178)   implicit none
   179)   
   180)   type(region_explicit_face_type), pointer :: RegionCreateExplicitFaceSet
   181)   
   182)   type(region_explicit_face_type), pointer :: explicit_faceset
   183)   
   184)   allocate(explicit_faceset)
   185)   nullify(explicit_faceset%face_centroids)
   186)   nullify(explicit_faceset%face_areas)
   187)   
   188)   RegionCreateExplicitFaceSet => explicit_faceset
   189) 
   190) end function RegionCreateExplicitFaceSet
   191) 
   192) ! ************************************************************************** !
   193) 
   194) function RegionCreateWithBlock(i1,i2,j1,j2,k1,k2)
   195)   ! 
   196)   ! Creates a region with i,j,k indices for arguments
   197)   ! 
   198)   ! Author: Glenn Hammond
   199)   ! Date: 10/23/07
   200)   ! 
   201) 
   202)   implicit none
   203)   
   204)   PetscInt :: i1, i2, j1, j2, k1, k2
   205)   
   206)   type(region_type), pointer :: RegionCreateWithBlock
   207) 
   208)   type(region_type), pointer :: region
   209)   
   210)   region => RegionCreateWithNothing()
   211)   region%i1 = i1
   212)   region%i2 = i2
   213)   region%j1 = j2
   214)   region%j2 = j2
   215)   region%k1 = k1
   216)   region%k2 = k2
   217)   region%num_cells = (abs(i2-i1)+1)*(abs(j2-j1)+1)* &
   218)                                     (abs(k2-k1)+1)
   219)                                     
   220)   RegionCreateWithBlock => region                                    
   221) 
   222) end function RegionCreateWithBlock
   223) 
   224) ! ************************************************************************** !
   225) 
   226) function RegionCreateWithList(list)
   227)   ! 
   228)   ! RegionCreate: Creates a region from a list of cells
   229)   ! 
   230)   ! Author: Glenn Hammond
   231)   ! Date: 10/23/07
   232)   ! 
   233) 
   234)   implicit none
   235)   
   236)   PetscInt :: list(:)
   237)   
   238)   type(region_type), pointer :: RegionCreateWithList
   239)   
   240)   type(region_type), pointer :: region
   241) 
   242)   region => RegionCreateWithNothing()
   243)   region%num_cells = size(list)
   244)   allocate(region%cell_ids(region%num_cells))
   245)   region%cell_ids = list
   246)   
   247)   RegionCreateWithList => region
   248) 
   249) end function RegionCreateWithList
   250) 
   251) ! ************************************************************************** !
   252) 
   253) function RegionCreateWithRegion(region)
   254)   ! 
   255)   ! Creates a copy of a region
   256)   ! 
   257)   ! Author: Glenn Hammond
   258)   ! Date: 02/22/08
   259)   ! 
   260) 
   261)   use Grid_Unstructured_Cell_module
   262) 
   263)   implicit none
   264)   
   265)   type(region_type), pointer :: RegionCreateWithRegion
   266)   type(region_type), pointer :: region
   267)   
   268)   type(region_type), pointer :: new_region
   269)   PetscInt :: icount, temp_int
   270)   
   271)   new_region => RegionCreateWithNothing()
   272)   
   273)   new_region%id = region%id
   274)   new_region%def_type = region%def_type
   275)   new_region%hdf5_ugrid_kludge = region%hdf5_ugrid_kludge
   276)   new_region%name = region%name
   277)   new_region%filename = region%filename
   278)   new_region%i1 = region%i1
   279)   new_region%i2 = region%i2
   280)   new_region%j1 = region%j1
   281)   new_region%j2 = region%j2
   282)   new_region%k1 = region%k1
   283)   new_region%k2 = region%k2
   284)   new_region%iface = region%iface
   285)   new_region%num_cells = region%num_cells
   286)   new_region%num_verts = region%num_verts
   287)   new_region%grid_type = region%grid_type
   288)   if (associated(region%coordinates)) then
   289)     call GeometryCopyCoordinates(region%coordinates, &
   290)                                  new_region%coordinates)
   291)   endif
   292)   if (associated(region%cell_ids)) then
   293)     allocate(new_region%cell_ids(new_region%num_cells))
   294)     new_region%cell_ids(1:new_region%num_cells) = &
   295)       region%cell_ids(1:region%num_cells)
   296)   endif
   297)   if (associated(region%faces)) then
   298)     allocate(new_region%faces(new_region%num_cells))
   299)     new_region%faces(1:new_region%num_cells) = &
   300)       region%faces(1:region%num_cells)
   301)   endif
   302)   if (associated(region%vertex_ids)) then
   303)     allocate(new_region%vertex_ids(0:MAX_VERT_PER_FACE,1:new_region%num_verts))
   304)     new_region%vertex_ids(0:MAX_VERT_PER_FACE,1:new_region%num_verts) = &
   305)     region%vertex_ids(0:MAX_VERT_PER_FACE,1:new_region%num_verts)
   306)   endif
   307)   if (associated(region%sideset)) then
   308)     new_region%sideset => RegionCreateSideSet()
   309)     new_region%sideset%nfaces = region%sideset%nfaces
   310)     allocate(new_region%sideset%face_vertices( &
   311)                size(region%sideset%face_vertices,1), &
   312)                size(region%sideset%face_vertices,2)))
   313)     new_region%sideset%face_vertices = region%sideset%face_vertices
   314)   endif
   315)   if (associated(region%explicit_faceset)) then
   316)     new_region%explicit_faceset => RegionCreateExplicitFaceSet()
   317)     allocate(new_region%explicit_faceset%face_centroids( &
   318)                size(region%explicit_faceset%face_centroids)))
   319)     new_region%explicit_faceset%face_centroids = &
   320)       region%explicit_faceset%face_centroids
   321)     do icount = 1, size(region%explicit_faceset%face_centroids)
   322)       new_region%explicit_faceset%face_centroids(icount)%x = &
   323)         region%explicit_faceset%face_centroids(icount)%x
   324)       new_region%explicit_faceset%face_centroids(icount)%y = &
   325)         region%explicit_faceset%face_centroids(icount)%y
   326)       new_region%explicit_faceset%face_centroids(icount)%z = &
   327)         region%explicit_faceset%face_centroids(icount)%z
   328)     enddo
   329)     allocate(new_region%explicit_faceset%face_areas( &
   330)                size(region%explicit_faceset%face_areas)))
   331)     new_region%explicit_faceset%face_areas = &
   332)       region%explicit_faceset%face_areas
   333)   endif
   334)   if (associated(region%polygonal_volume)) then
   335)     new_region%polygonal_volume => GeometryCreatePolygonalVolume()
   336)     if (associated(region%polygonal_volume%xy_coordinates)) then
   337)       call GeometryCopyCoordinates(region%polygonal_volume%xy_coordinates, &
   338)                                    new_region%polygonal_volume%xy_coordinates)
   339)     endif
   340)     if (associated(region%polygonal_volume%xz_coordinates)) then
   341)       call GeometryCopyCoordinates(region%polygonal_volume%xz_coordinates, &
   342)                                    new_region%polygonal_volume%xz_coordinates)
   343)     endif
   344)     if (associated(region%polygonal_volume%yz_coordinates)) then
   345)       call GeometryCopyCoordinates(region%polygonal_volume%yz_coordinates, &
   346)                                    new_region%polygonal_volume%yz_coordinates)
   347)     endif
   348)   endif
   349)   
   350)   RegionCreateWithRegion => new_region
   351)   
   352) end function RegionCreateWithRegion
   353) 
   354) ! ************************************************************************** !
   355) 
   356) subroutine RegionInitList(list)
   357)   ! 
   358)   ! Initializes a region list
   359)   ! 
   360)   ! Author: Glenn Hammond
   361)   ! Date: 10/29/07
   362)   ! 
   363) 
   364)   implicit none
   365) 
   366)   type(region_list_type) :: list
   367)   
   368)   nullify(list%first)
   369)   nullify(list%last)
   370)   nullify(list%array)
   371)   list%num_regions = 0
   372) 
   373) end subroutine RegionInitList
   374) 
   375) ! ************************************************************************** !
   376) 
   377) subroutine RegionAddToList(new_region,list)
   378)   ! 
   379)   ! Adds a new region to a region list
   380)   ! 
   381)   ! Author: Glenn Hammond
   382)   ! Date: 10/29/07
   383)   ! 
   384) 
   385)   implicit none
   386)   
   387)   type(region_type), pointer :: new_region
   388)   type(region_list_type) :: list
   389)   
   390)   list%num_regions = list%num_regions + 1
   391)   new_region%id = list%num_regions
   392)   if (.not.associated(list%first)) list%first => new_region
   393)   if (associated(list%last)) list%last%next => new_region
   394)   list%last => new_region
   395)   
   396) end subroutine RegionAddToList
   397) 
   398) ! ************************************************************************** !
   399) 
   400) subroutine RegionRead(region,input,option)
   401)   ! 
   402)   ! Reads a region from the input file
   403)   ! 
   404)   ! Author: Glenn Hammond
   405)   ! Date: 02/20/08
   406)   ! 
   407) 
   408)   use Input_Aux_module
   409)   use String_module
   410)   use Option_module
   411)   use Grid_Structured_module
   412)   
   413)   implicit none
   414)   
   415)   type(option_type) :: option
   416)   type(region_type) :: region
   417)   type(input_type), pointer :: input
   418)   
   419)   character(len=MAXWORDLENGTH) :: keyword, word
   420) 
   421)   input%ierr = 0
   422)   do
   423)   
   424)     call InputReadPflotranString(input,option)
   425)     if (InputError(input)) exit
   426)     if (InputCheckExit(input,option)) exit
   427)     
   428)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   429)     call InputErrorMsg(input,option,'keyword','REGION')
   430)     call StringToUpper(keyword)   
   431) 
   432)     select case(trim(keyword))
   433)     
   434)       case('BLOCK')
   435)         region%def_type = DEFINED_BY_BLOCK
   436)         call InputReadInt(input,option,region%i1)
   437)         if (InputError(input)) then
   438)           input%ierr = 0
   439)           call InputReadPflotranString(input,option)
   440)           call InputReadStringErrorMsg(input,option,'REGION')
   441)           call InputReadInt(input,option,region%i1) 
   442)         endif
   443)         call InputErrorMsg(input,option,'i1','REGION')
   444)         call InputReadInt(input,option,region%i2)
   445)         call InputErrorMsg(input,option,'i2','REGION')
   446)         call InputReadInt(input,option,region%j1)
   447)         call InputErrorMsg(input,option,'j1','REGION')
   448)         call InputReadInt(input,option,region%j2)
   449)         call InputErrorMsg(input,option,'j2','REGION')
   450)         call InputReadInt(input,option,region%k1)
   451)         call InputErrorMsg(input,option,'k1','REGION')
   452)         call InputReadInt(input,option,region%k2)
   453)         call InputErrorMsg(input,option,'k2','REGION')
   454)       case('CARTESIAN_BOUNDARY')
   455)         region%def_type = DEFINED_BY_CARTESIAN_BOUNDARY
   456)         call InputReadWord(input,option,word,PETSC_TRUE)
   457)         call InputErrorMsg(input,option,'cartesian boundary face','REGION')
   458)         call StringToUpper(word)
   459)         select case(word)
   460)           case('WEST')
   461)             region%iface = WEST_FACE
   462)           case('EAST')
   463)             region%iface = EAST_FACE
   464)           case('NORTH')
   465)             region%iface = NORTH_FACE
   466)           case('SOUTH')
   467)             region%iface = SOUTH_FACE
   468)           case('BOTTOM')
   469)             region%iface = BOTTOM_FACE
   470)           case('TOP')
   471)             region%iface = TOP_FACE
   472)           case default
   473)             option%io_buffer = 'Cartesian boundary face "' // trim(word) // &
   474)               '" not recognized.'
   475)             call printErrMsg(option)
   476)         end select
   477)       case('COORDINATE')
   478)         region%def_type = DEFINED_BY_COORD
   479)         allocate(region%coordinates(1))
   480)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%x) 
   481)         if (InputError(input)) then
   482)           input%ierr = 0
   483)           call InputReadPflotranString(input,option)
   484)           call InputReadStringErrorMsg(input,option,'REGION')
   485)           call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%x)
   486)         endif
   487)         call InputErrorMsg(input,option,'x-coordinate','REGION')
   488)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%y)
   489)         call InputErrorMsg(input,option,'y-coordinate','REGION')
   490)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%z)
   491)         call InputErrorMsg(input,option,'z-coordinate','REGION')
   492)       case('COORDINATES')
   493)         region%def_type = DEFINED_BY_COORD
   494)         call GeometryReadCoordinates(input,option,region%name, &
   495)                                      region%coordinates)
   496)       case('POLYGON')
   497)         if (.not.associated(region%polygonal_volume)) then
   498)           region%polygonal_volume => GeometryCreatePolygonalVolume()
   499)         endif
   500)         do
   501)           call InputReadPflotranString(input,option)
   502)           if (InputError(input)) exit
   503)           if (InputCheckExit(input,option)) exit
   504)           call InputReadWord(input,option,word,PETSC_TRUE)
   505)           call InputErrorMsg(input,option,'keyword','REGION')
   506)           call StringToUpper(word)   
   507)           select case(trim(word))
   508)             case('TYPE')
   509)               call InputReadWord(input,option,word,PETSC_TRUE)
   510)               call InputErrorMsg(input,option,'polygon type','REGION')
   511)               call StringToUpper(word)
   512)               select case(word)
   513)                 case('BOUNDARY_FACES_IN_VOLUME')
   514)                   region%def_type = DEFINED_BY_POLY_BOUNDARY_FACE
   515)                 case('CELL_CENTERS_IN_VOLUME')
   516)                   region%def_type = DEFINED_BY_POLY_CELL_CENTER
   517)                 case default
   518)                   option%io_buffer = 'REGION->POLYGON->"' // trim(word) // &
   519)                     '" not recognized.'
   520)                   call printErrMsg(option)
   521)               end select
   522)             case('XY')
   523)               call GeometryReadCoordinates(input,option,region%name, &
   524)                                          region%polygonal_volume%xy_coordinates)
   525)             case('XZ')
   526)               call GeometryReadCoordinates(input,option,region%name, &
   527)                                          region%polygonal_volume%xz_coordinates)
   528)             case('YZ')
   529)               call GeometryReadCoordinates(input,option,region%name, &
   530)                                          region%polygonal_volume%yz_coordinates)
   531)             case default
   532)               option%io_buffer = 'Keyword not recognized for REGION POLYGON.'
   533)               call printErrMsg(option)
   534)           end select
   535)         enddo
   536)       case('FILE')
   537)         call InputReadNChars(input,option,region%filename,MAXSTRINGLENGTH,PETSC_TRUE)
   538)         call InputErrorMsg(input,option,'filename','REGION')
   539)         call InputReadWord(input,option,word,PETSC_TRUE)
   540)         if (input%ierr == 0) then
   541)           if (StringCompareIgnoreCase(word,'OLD_FORMAT')) then
   542)             region%hdf5_ugrid_kludge = PETSC_TRUE
   543)           endif
   544)         endif
   545)       case('LIST')
   546)         option%io_buffer = 'REGION LIST currently not implemented'
   547)         call printErrMsg(option)
   548)       case('FACE')
   549)         call InputReadWord(input,option,word,PETSC_TRUE)
   550)         call InputErrorMsg(input,option,'face','REGION')
   551)         call StringToUpper(word)
   552)         select case(word)
   553)           case('WEST')
   554)             region%iface = WEST_FACE
   555)           case('EAST')
   556)             region%iface = EAST_FACE
   557)           case('NORTH')
   558)             region%iface = NORTH_FACE
   559)           case('SOUTH')
   560)             region%iface = SOUTH_FACE
   561)           case('BOTTOM')
   562)             region%iface = BOTTOM_FACE
   563)           case('TOP')
   564)             region%iface = TOP_FACE
   565)           case default
   566)             option%io_buffer = 'FACE "' // trim(word) // &
   567)               '" not recognized.'
   568)             call printErrMsg(option)
   569)         end select
   570)       case('GRID','SURF_GRID')
   571)         call InputReadWord(input,option,word,PETSC_TRUE)
   572)         call InputErrorMsg(input,option,'GRID','REGION')
   573)         call StringToUpper(word)
   574)         select case(trim(word))
   575)           case('STRUCTURED')
   576)             region%grid_type = STRUCTURED_GRID_REGION
   577)           case('UNSTRUCTURED')
   578)             region%grid_type = UNSTRUCTURED_GRID_REGION
   579)           case default
   580)             option%io_buffer = 'REGION keyword: GRID = '//trim(word)//'not supported yet'
   581)           call printErrMsg(option)
   582)         end select
   583)       case default
   584)         call InputKeywordUnrecognized(keyword,'REGION',option)
   585)     end select
   586)   enddo
   587)  
   588) end subroutine RegionRead
   589) 
   590) ! ************************************************************************** !
   591) 
   592) subroutine RegionReadFromFilename(region,option,filename)
   593)   ! 
   594)   ! Reads a list of cells from a file named filename
   595)   ! 
   596)   ! Author: Glenn Hammond
   597)   ! Date: 10/29/07
   598)   ! 
   599) 
   600)   use Input_Aux_module
   601)   use Option_module
   602)   use Utility_module
   603)   
   604)   implicit none
   605)   
   606)   type(region_type) :: region
   607)   type(option_type) :: option
   608)   type(input_type), pointer :: input
   609)   character(len=MAXSTRINGLENGTH) :: filename
   610)   
   611)   input => InputCreate(IUNIT_TEMP,filename,option)
   612)   call RegionReadFromFileId(region,input,option)          
   613)   call InputDestroy(input)         
   614) 
   615) end subroutine RegionReadFromFilename
   616) 
   617) ! ************************************************************************** !
   618) 
   619) subroutine RegionReadFromFileId(region,input,option)
   620)   ! 
   621)   ! Reads a list of cells from an open file
   622)   ! 
   623)   ! Author: Glenn Hammond
   624)   ! Date: 10/29/07
   625)   ! 
   626) 
   627)   use Input_Aux_module
   628)   use Option_module
   629)   use Utility_module
   630)   use Logging_module
   631)   use Grid_Unstructured_Cell_module
   632)   
   633)   implicit none
   634)   
   635)   type(region_type) :: region
   636)   type(option_type) :: option
   637)   type(input_type), pointer :: input
   638)   
   639)   character(len=MAXWORDLENGTH) :: word
   640)   character(len=1) :: backslash
   641) 
   642)   PetscInt, pointer :: temp_int_array(:)
   643)   PetscInt, pointer :: cell_ids_p(:)
   644)   PetscInt, pointer :: face_ids_p(:)
   645)   PetscInt, pointer :: vert_id_0_p(:)
   646)   PetscInt, pointer :: vert_id_1_p(:)
   647)   PetscInt, pointer :: vert_id_2_p(:)
   648)   PetscInt, pointer :: vert_id_3_p(:)
   649)   PetscInt, pointer :: vert_id_4_p(:)
   650)   PetscInt :: max_size
   651)   PetscInt :: count
   652)   PetscInt :: temp_int
   653)   PetscInt :: input_data_type
   654)   PetscInt :: ii
   655)   PetscInt :: istart
   656)   PetscInt :: iend
   657)   PetscInt :: remainder
   658)   PetscErrorCode :: ierr
   659) 
   660)   PetscInt, parameter :: CELL_IDS_ONLY = 1
   661)   PetscInt, parameter :: CELL_IDS_WITH_FACE_IDS = 2
   662)   PetscInt, parameter :: VERTEX_IDS = 3
   663) 
   664)   call PetscLogEventBegin(logging%event_region_read_ascii,ierr);CHKERRQ(ierr)
   665)   
   666)   !TODO(geh): clean and optimize this subroutine
   667)   
   668)   max_size = 1000
   669)   backslash = achar(92)  ! 92 = "\" Some compilers choke on \" thinking it
   670)                           ! is a double quote as in c/c++
   671)   
   672)   allocate(temp_int_array(max_size))
   673)   allocate(cell_ids_p(max_size))
   674)   allocate(face_ids_p(max_size))
   675)   allocate(vert_id_0_p(max_size))
   676)   allocate(vert_id_0_p(max_size))
   677)   allocate(vert_id_1_p(max_size))
   678)   allocate(vert_id_2_p(max_size))
   679)   allocate(vert_id_3_p(max_size))
   680)   allocate(vert_id_4_p(max_size))
   681)   
   682)   temp_int_array = 0
   683)   cell_ids_p = 0
   684)   face_ids_p = 0
   685)   vert_id_0_p = 0
   686)   vert_id_1_p = -1
   687)   vert_id_2_p = -1
   688)   vert_id_3_p = -1
   689)   vert_id_4_p = -1
   690)   
   691)   
   692)   count = 0
   693) 
   694)   ! Determine if region definition in the input data is one of the following:
   695)   !  1) Contains cell ids only : Only ONE entry per line
   696)   !  2) Contains cell ids and face ids: TWO entries per line
   697)   !  3) Contains vertex ids that make up the face: MORE than two entries per
   698)   !     line
   699)   count = 0
   700)   call InputReadPflotranString(input, option)
   701)   do 
   702)     call InputReadInt(input, option, temp_int)
   703)     if (InputError(input)) exit
   704)     count = count + 1
   705)     temp_int_array(count) = temp_int
   706)   enddo
   707) 
   708)   if (count == 1) then
   709)     !
   710)     ! Input data contains only cell ids
   711)     !
   712)     input_data_type = CELL_IDS_ONLY
   713)     cell_ids_p(1) = temp_int_array(1)
   714)     count = 1
   715)     region%def_type = DEFINED_BY_CELL_IDS
   716) 
   717)     ! Read the data
   718)     do
   719)       call InputReadPflotranString(input, option)
   720)       if (InputError(input)) exit
   721)       call InputReadInt(input, option, temp_int)
   722)       if (.not.InputError(input)) then
   723)         count = count + 1
   724)         cell_ids_p(count) = temp_int
   725)       endif
   726)       if (count+1 > max_size) then ! resize temporary array
   727)         call reallocateIntArray(cell_ids_p, max_size)
   728)       endif
   729)     enddo
   730) 
   731)     ! Depending on processor rank, save only a portion of data
   732)     region%num_cells = count/option%mycommsize
   733)       remainder = count - region%num_cells*option%mycommsize
   734)     if (option%myrank < remainder) region%num_cells = region%num_cells + 1
   735)     istart = 0
   736)     iend   = 0
   737)     call MPI_Exscan(region%num_cells, istart, ONE_INTEGER_MPI, MPIU_INTEGER, &
   738)                     MPI_SUM, option%mycomm, ierr)
   739)     call MPI_Scan(region%num_cells, iend, ONE_INTEGER_MPI, MPIU_INTEGER, &
   740)                    MPI_SUM, option%mycomm, ierr)
   741) 
   742)     ! Allocate memory and save the data
   743)     region%num_cells = iend - istart
   744)     allocate(region%cell_ids(region%num_cells))
   745)     region%cell_ids(1:region%num_cells) = cell_ids_p(istart+1:iend)
   746)     deallocate(cell_ids_p)
   747) 
   748)   else if (count == 2) then
   749)     !
   750)     ! Input data contains cell ids + face ids
   751)     !
   752)     input_data_type = CELL_IDS_WITH_FACE_IDS
   753)     cell_ids_p(1) = temp_int_array(1)
   754)     face_ids_p(1) = temp_int_array(2)
   755)     count = 1 ! reset the counter to represent the num of rows read
   756)     region%def_type = DEFINED_BY_CELL_AND_FACE_IDS
   757) 
   758)     ! Read the data
   759)     do
   760)       call InputReadPflotranString(input, option)
   761)       if (InputError(input)) exit
   762)       call InputReadInt(input, option, temp_int)
   763)       if (InputError(input)) exit
   764)       count = count + 1
   765)       cell_ids_p(count) = temp_int
   766) 
   767)       call InputReadInt(input,option,temp_int)
   768)       if (InputError(input)) then
   769)         option%io_buffer = 'ERROR while reading the region "' // &
   770)           trim(region%name) // '" from file'
   771)         call printErrMsg(option)
   772)       endif
   773)       face_ids_p(count) = temp_int
   774)       if (count+1 > max_size) then ! resize temporary array
   775)         call reallocateIntArray(cell_ids_p, max_size)
   776)         call reallocateIntArray(face_ids_p, max_size)
   777)       endif
   778)     enddo
   779) 
   780)     ! Depending on processor rank, save only a portion of data
   781)     region%num_cells = count/option%mycommsize
   782)       remainder = count - region%num_cells*option%mycommsize
   783)     if (option%myrank < remainder) region%num_cells = region%num_cells + 1
   784)     istart = 0
   785)     iend   = 0
   786)     call MPI_Exscan(region%num_cells,istart,ONE_INTEGER_MPI,MPIU_INTEGER, &
   787)                     MPI_SUM,option%mycomm,ierr)
   788)     call MPI_Scan(region%num_cells,iend,ONE_INTEGER_MPI,MPIU_INTEGER, &
   789)                    MPI_SUM,option%mycomm,ierr)
   790) 
   791)     ! Allocate memory and save the data
   792)     allocate(region%cell_ids(region%num_cells))
   793)     allocate(region%faces(region%num_cells))
   794)     region%cell_ids(1:region%num_cells) = cell_ids_p(istart + 1:iend)
   795)     region%faces(1:region%num_cells) = face_ids_p(istart + 1:iend)
   796)     deallocate(cell_ids_p)
   797)     deallocate(face_ids_p)
   798) 
   799)   else
   800)     !
   801)     ! Input data contains vertices
   802)     !
   803)     input_data_type = VERTEX_IDS
   804)     vert_id_0_p(1) = temp_int_array(1)
   805)     vert_id_1_p(1) = temp_int_array(2)
   806)     vert_id_2_p(1) = temp_int_array(3)
   807)     vert_id_3_p(1) = temp_int_array(4)
   808)     if (vert_id_0_p(1) == 4 ) vert_id_4_p(1) = temp_int_array(5)
   809)     count = 1 ! reset the counter to represent the num of rows read
   810)     region%def_type = DEFINED_BY_VERTEX_IDS
   811) 
   812)     ! Read the data
   813)     do
   814)       call InputReadPflotranString(input,option)
   815)       if (InputError(input)) exit
   816)       call InputReadInt(input,option,temp_int)
   817)       if (InputError(input)) exit
   818)       count = count + 1
   819)       vert_id_0_p(count) = temp_int
   820) 
   821)       vert_id_4_p(count) = UNINITIALIZED_INTEGER
   822)       do ii = 1, vert_id_0_p(count)
   823)         call InputReadInt(input,option,temp_int)
   824)         if (InputError(input)) then
   825)           option%io_buffer = 'ERROR while reading the region "' // &
   826)             trim(region%name) // '" from file'
   827)           call printErrMsg(option)
   828)         endif
   829) 
   830)         select case(ii)
   831)           case(1)
   832)             vert_id_1_p(count) = temp_int
   833)           case(2)
   834)             vert_id_2_p(count) = temp_int
   835)           case(3)
   836)             vert_id_3_p(count) = temp_int
   837)           case(4)
   838)             vert_id_4_p(count) = temp_int
   839)         end select
   840) 
   841)         if (count+1 > max_size) then ! resize temporary array
   842)           call reallocateIntArray(vert_id_0_p,max_size)
   843)           call reallocateIntArray(vert_id_1_p,max_size)
   844)           call reallocateIntArray(vert_id_2_p,max_size)
   845)           call reallocateIntArray(vert_id_3_p,max_size)
   846)           call reallocateIntArray(vert_id_4_p,max_size)
   847)         endif
   848)       enddo
   849)     enddo
   850) 
   851)     ! Depending on processor rank, save only a portion of data
   852)     region%num_verts = count/option%mycommsize
   853)       remainder = count - region%num_verts*option%mycommsize
   854)     if (option%myrank < remainder) region%num_verts = region%num_verts + 1
   855)     istart = 0
   856)     iend   = 0
   857)     call MPI_Exscan(region%num_verts,istart,ONE_INTEGER_MPI,MPIU_INTEGER, &
   858)                     MPI_SUM,option%mycomm,ierr)
   859)     call MPI_Scan(region%num_verts,iend,ONE_INTEGER_MPI,MPIU_INTEGER, &
   860)                    MPI_SUM,option%mycomm,ierr)
   861) 
   862)     ! Allocate memory and save the data
   863)     region%num_verts = iend - istart
   864)     allocate(region%vertex_ids(0:MAX_VERT_PER_FACE,1:region%num_verts))
   865)     region%vertex_ids(0,1:region%num_verts) = vert_id_0_p(istart + 1: iend)
   866)     region%vertex_ids(1,1:region%num_verts) = vert_id_1_p(istart + 1: iend)
   867)     region%vertex_ids(2,1:region%num_verts) = vert_id_2_p(istart + 1: iend)
   868)     region%vertex_ids(3,1:region%num_verts) = vert_id_3_p(istart + 1: iend)
   869)     region%vertex_ids(4,1:region%num_verts) = vert_id_4_p(istart + 1: iend)
   870)     deallocate(vert_id_0_p)
   871)     deallocate(vert_id_1_p)
   872)     deallocate(vert_id_2_p)
   873)     deallocate(vert_id_3_p)
   874)     deallocate(vert_id_4_p)
   875) 
   876)   endif
   877)   
   878) #if 0  
   879)   count = 1
   880)   do
   881)     call InputReadPflotranString(input,option)
   882)     if (InputError(input)) exit
   883)     call InputReadInt(input,option,temp_int)
   884)     if (.not.InputError(input)) then
   885)       count = count + 1
   886)       temp_int_array(count) = temp_int
   887)       write(*,*) count,temp_int
   888)     endif
   889)     if (count+1 > max_size) then ! resize temporary array
   890)       call reallocateIntArray(temp_int_array,max_size)
   891)     endif
   892)   enddo
   893) 
   894)   if (count > 0) then
   895)     region%num_cells = count
   896)     allocate(region%cell_ids(count))
   897)     region%cell_ids(1:count) = temp_int_array(1:count)
   898)   else
   899)     region%num_cells = 0
   900)     nullify(region%cell_ids)
   901)   endif
   902) #endif
   903)   deallocate(temp_int_array) 
   904) 
   905)   call PetscLogEventEnd(logging%event_region_read_ascii,ierr);CHKERRQ(ierr)
   906) 
   907) end subroutine RegionReadFromFileId
   908) 
   909) ! ************************************************************************** !
   910) 
   911) subroutine RegionReadSideSet(sideset,filename,option)
   912)   ! 
   913)   ! Reads an unstructured grid sideset
   914)   ! 
   915)   ! Author: Glenn Hammond
   916)   ! Date: 12/19/11
   917)   ! 
   918) 
   919)   use Input_Aux_module
   920)   use Option_module
   921)   use String_module
   922)   
   923)   implicit none
   924)   
   925)   type(region_sideset_type) :: sideset
   926)   character(len=MAXSTRINGLENGTH) :: filename
   927)   type(option_type) :: option
   928)   
   929)   type(input_type), pointer :: input
   930)   character(len=MAXSTRINGLENGTH) :: string, hint
   931)   character(len=MAXWORDLENGTH) :: word
   932)   PetscInt :: num_faces_local_save
   933)   PetscInt :: num_faces_local
   934)   PetscInt :: num_to_read
   935)   PetscInt, parameter :: max_nvert_per_face = 4
   936)   PetscInt, allocatable :: temp_int_array(:,:)
   937) 
   938)   PetscInt :: iface, ivertex, irank, num_vertices
   939)   PetscInt :: remainder
   940)   PetscErrorCode :: ierr
   941)   PetscMPIInt :: status_mpi(MPI_STATUS_SIZE)
   942)   PetscMPIInt :: int_mpi
   943)   PetscInt :: fileid
   944)   
   945)   fileid = 86
   946)   input => InputCreate(fileid,filename,option)
   947) 
   948) ! Format of sideset file
   949) ! type: T=triangle, Q=quadrilateral
   950) ! vertn(Q) = 4
   951) ! vertn(T) = 3
   952) ! -----------------------------------------------------------------
   953) ! num_faces  (integer)
   954) ! type vert1 vert2 ... vertn  ! for face 1 (integers)
   955) ! type vert1 vert2 ... vertn  ! for face 2
   956) ! ...
   957) ! ...
   958) ! type vert1 vert2 ... vertn  ! for face num_faces
   959) ! -----------------------------------------------------------------
   960) 
   961)   hint = 'Unstructured Sideset'
   962) 
   963)   call InputReadPflotranString(input,option)
   964)   string = 'unstructured sideset'
   965)   call InputReadStringErrorMsg(input,option,hint)  
   966) 
   967)   ! read num_faces
   968)   call InputReadInt(input,option,sideset%nfaces)
   969)   call InputErrorMsg(input,option,'number of faces',hint)
   970) 
   971)   ! divide faces across ranks
   972)   num_faces_local = sideset%nfaces/option%mycommsize 
   973)   num_faces_local_save = num_faces_local
   974)   remainder = sideset%nfaces - num_faces_local*option%mycommsize
   975)   if (option%myrank < remainder) num_faces_local = &
   976)                                  num_faces_local + 1
   977) 
   978)   ! allocate array to store vertices for each faces
   979)   allocate(sideset%face_vertices(max_nvert_per_face, &
   980)                                  num_faces_local))
   981)   sideset%face_vertices = UNINITIALIZED_INTEGER
   982) 
   983)   ! for now, read all faces from ASCII file through io_rank and communicate
   984)   ! to other ranks
   985)   if (option%myrank == option%io_rank) then
   986)     allocate(temp_int_array(max_nvert_per_face, &
   987)                             num_faces_local_save+1))
   988)     ! read for other processors
   989)     do irank = 0, option%mycommsize-1
   990)       temp_int_array = UNINITIALIZED_INTEGER
   991)       num_to_read = num_faces_local_save
   992)       if (irank < remainder) num_to_read = num_to_read + 1
   993) 
   994)       do iface = 1, num_to_read
   995)         ! read in the vertices defining the cell face
   996)         call InputReadPflotranString(input,option)
   997)         call InputReadStringErrorMsg(input,option,hint)  
   998)         call InputReadWord(input,option,word,PETSC_TRUE)
   999)         call InputErrorMsg(input,option,'face type',hint)
  1000)         call StringToUpper(word)
  1001)         select case(word)
  1002)           case('Q')
  1003)             num_vertices = 4
  1004)           case('T')
  1005)             num_vertices = 3
  1006)           case('L')
  1007)             num_vertices = 2
  1008)           case default
  1009)             option%io_buffer = 'Unknown face type: ' // trim(word)
  1010)         end select
  1011)         do ivertex = 1, num_vertices
  1012)           call InputReadInt(input,option,temp_int_array(ivertex,iface))
  1013)           call InputErrorMsg(input,option,'vertex id',hint)
  1014)         enddo
  1015)       enddo
  1016)       ! if the faces reside on io_rank
  1017)       if (irank == option%io_rank) then
  1018) #if UGRID_DEBUG
  1019)         write(string,*) num_faces_local
  1020)         string = trim(adjustl(string)) // ' faces stored on p0'
  1021)         print *, trim(string)
  1022) #endif
  1023)         sideset%nfaces = num_faces_local
  1024)         sideset%face_vertices(:,1:num_faces_local) = &
  1025)           temp_int_array(:,1:num_faces_local)
  1026)       else
  1027)         ! otherwise communicate to other ranks
  1028) #if UGRID_DEBUG
  1029)         write(string,*) num_to_read
  1030)         write(word,*) irank
  1031)         string = trim(adjustl(string)) // ' faces sent from p0 to p' // &
  1032)                  trim(adjustl(word))
  1033)         print *, trim(string)
  1034) #endif
  1035)         int_mpi = num_to_read*max_nvert_per_face
  1036)         call MPI_Send(temp_int_array,int_mpi,MPIU_INTEGER,irank, &
  1037)                       num_to_read,option%mycomm,ierr)
  1038)       endif
  1039)     enddo
  1040)     deallocate(temp_int_array)
  1041)   else
  1042)     ! other ranks post the recv
  1043) #if UGRID_DEBUG
  1044)         write(string,*) num_faces_local
  1045)         write(word,*) option%myrank
  1046)         string = trim(adjustl(string)) // ' faces received from p0 at p' // &
  1047)                  trim(adjustl(word))
  1048)         print *, trim(string)
  1049) #endif
  1050)     sideset%nfaces = num_faces_local
  1051)     int_mpi = num_faces_local*max_nvert_per_face
  1052)     call MPI_Recv(sideset%face_vertices,int_mpi, &
  1053)                   MPIU_INTEGER,option%io_rank, &
  1054)                   MPI_ANY_TAG,option%mycomm,status_mpi,ierr)
  1055)   endif
  1056) 
  1057) !  unstructured_grid%nlmax = num_faces_local
  1058) !  unstructured_grid%num_vertices_local = num_vertices_local
  1059) 
  1060)   call InputDestroy(input)
  1061) 
  1062) end subroutine RegionReadSideSet
  1063) 
  1064) ! ************************************************************************** !
  1065) 
  1066) subroutine RegionReadExplicitFaceSet(explicit_faceset,cell_ids,filename,option)
  1067)   ! 
  1068)   ! Reads an unstructured grid explicit region
  1069)   ! 
  1070)   ! Author: Glenn Hammond
  1071)   ! Date: 05/18/12
  1072)   ! 
  1073) 
  1074)   use Input_Aux_module
  1075)   use Option_module
  1076)   use String_module
  1077)   
  1078)   implicit none
  1079)   
  1080)   type(region_explicit_face_type), pointer :: explicit_faceset
  1081)   PetscInt, pointer :: cell_ids(:)
  1082)   character(len=MAXSTRINGLENGTH) :: filename
  1083)   type(option_type) :: option
  1084)   
  1085)   type(input_type), pointer :: input
  1086)   character(len=MAXSTRINGLENGTH) :: string, hint
  1087)   character(len=MAXWORDLENGTH) :: word
  1088)   PetscInt :: fileid
  1089)   
  1090)   PetscInt :: num_connections
  1091)   PetscInt :: iconn
  1092)   
  1093)   explicit_faceset => RegionCreateExplicitFaceSet()
  1094)   
  1095)   fileid = 86
  1096)   input => InputCreate(fileid,filename,option)
  1097)   
  1098) ! Format of explicit unstructured grid file
  1099) ! id_ = integer
  1100) ! x_, y_, z_, area_ = real
  1101) ! definitions
  1102) ! id_ = id of grid cell
  1103) ! x_ = x coordinate of cell face
  1104) ! y_ = y coordinate of cell face
  1105) ! z_ = z coordinate of cell face
  1106) ! area_ = area of grid cell face
  1107) ! -----------------------------------------------------------------
  1108) ! CONNECTIONS <integer>   integer = # connections (M)
  1109) ! id_1 x_1 y_1 z_1 area_1
  1110) ! id_2 x_2 y_2 z_2 area_2
  1111) ! ...
  1112) ! ...
  1113) ! id_M x_M y_M z_M area_M
  1114) ! -----------------------------------------------------------------
  1115) 
  1116)   do
  1117)     call InputReadPflotranString(input,option)
  1118)     if (InputError(input)) exit
  1119) 
  1120)     call InputReadWord(input,option,word,PETSC_FALSE)
  1121)     call StringToUpper(word)
  1122)     hint = trim(word)
  1123)   
  1124)     select case(word)
  1125)       case('CONNECTIONS')
  1126)         hint = 'Explicit Unstructured Grid CONNECTIONS'
  1127)         call InputReadInt(input,option,num_connections)
  1128)         call InputErrorMsg(input,option,'number of connections',hint)
  1129)         
  1130)         allocate(cell_ids(num_connections))
  1131)         cell_ids = 0
  1132)         allocate(explicit_faceset%face_areas(num_connections))
  1133)         explicit_faceset%face_areas = 0    
  1134)         allocate(explicit_faceset%face_centroids(num_connections))
  1135)         do iconn = 1, num_connections
  1136)           explicit_faceset%face_centroids(iconn)%x = 0.d0
  1137)           explicit_faceset%face_centroids(iconn)%y = 0.d0
  1138)           explicit_faceset%face_centroids(iconn)%z = 0.d0
  1139)         enddo
  1140)         do iconn = 1, num_connections
  1141)           call InputReadPflotranString(input,option)
  1142)           call InputReadStringErrorMsg(input,option,hint)  
  1143)           call InputReadInt(input,option,cell_ids(iconn))
  1144)           call InputErrorMsg(input,option,'cell id',hint)
  1145)           call InputReadDouble(input,option, &
  1146)                                explicit_faceset%face_centroids(iconn)%x)
  1147)           call InputErrorMsg(input,option,'face x coordinate',hint)
  1148)           call InputReadDouble(input,option, &
  1149)                                explicit_faceset%face_centroids(iconn)%y)
  1150)           call InputErrorMsg(input,option,'face y coordinate',hint)
  1151)           call InputReadDouble(input,option, &
  1152)                                explicit_faceset%face_centroids(iconn)%z)
  1153)           call InputErrorMsg(input,option,'face z coordinate',hint)
  1154)           call InputReadDouble(input,option, &
  1155)                                explicit_faceset%face_areas(iconn))
  1156)           call InputErrorMsg(input,option,'face area',hint)
  1157)         enddo
  1158)       case default
  1159)         call InputKeywordUnrecognized(word, &
  1160)                'REGION (explicit unstructured grid)',option)
  1161)     end select
  1162)   enddo
  1163) 
  1164)   call InputDestroy(input)
  1165) 
  1166) end subroutine RegionReadExplicitFaceSet
  1167) 
  1168) ! ************************************************************************** !
  1169) 
  1170) function RegionGetPtrFromList(region_name,region_list)
  1171)   ! 
  1172)   ! Returns a pointer to the region matching region_name
  1173)   ! 
  1174)   ! Author: Glenn Hammond
  1175)   ! Date: 11/01/07
  1176)   ! 
  1177) 
  1178)   use String_module
  1179) 
  1180)   implicit none
  1181)   
  1182)   type(region_type), pointer :: RegionGetPtrFromList
  1183)   character(len=MAXWORDLENGTH) :: region_name
  1184)   PetscInt :: length
  1185)   type(region_list_type) :: region_list
  1186) 
  1187)   type(region_type), pointer :: region
  1188)     
  1189)   nullify(RegionGetPtrFromList)
  1190)   region => region_list%first
  1191)   
  1192)   do 
  1193)     if (.not.associated(region)) exit
  1194)     length = len_trim(region_name)
  1195)     if (length == len_trim(region%name) .and. &
  1196)         StringCompare(region%name,region_name,length)) then
  1197)       RegionGetPtrFromList => region
  1198)       return
  1199)     endif
  1200)     region => region%next
  1201)   enddo
  1202)   
  1203) end function RegionGetPtrFromList
  1204) 
  1205) ! **************************************************************************** !
  1206) 
  1207) subroutine RegionInputRecord(region_list)
  1208)   ! 
  1209)   ! Prints ingested region information to the input record file
  1210)   ! 
  1211)   ! Author: Jenn Frederick
  1212)   ! Date: 03/30/2016
  1213)   ! 
  1214)   use Grid_Structured_module
  1215) 
  1216)   implicit none
  1217) 
  1218)   type(region_list_type), pointer :: region_list
  1219)   
  1220)   type(region_type), pointer :: cur_region
  1221)   character(len=MAXWORDLENGTH) :: word1, word2
  1222)   character(len=MAXSTRINGLENGTH) :: string
  1223)   PetscInt :: k
  1224)   PetscInt :: id = INPUT_RECORD_UNIT
  1225)   character(len=10) :: sFormat, iFormat
  1226)   
  1227)   sFormat = '(ES14.7)'
  1228)   iFormat = '(I10)'
  1229) 
  1230)   write(id,'(a)') ' '
  1231)   write(id,'(a)') '---------------------------------------------------------&
  1232)                   &-----------------------'
  1233)   write(id,'(a29)',advance='no') '---------------------------: '
  1234)   write(id,'(a)') 'REGIONS'
  1235)   
  1236)   cur_region => region_list%first
  1237)   do
  1238)     if (.not.associated(cur_region)) exit
  1239)     write(id,'(a29)',advance='no') 'region: '
  1240)     write(id,'(a)') adjustl(trim(cur_region%name))
  1241)     if (len_trim(cur_region%filename) > 0) then
  1242)       write(id,'(a29)',advance='no') 'from file: '
  1243)       write(id,'(a)') adjustl(trim(cur_region%filename)) 
  1244)     endif
  1245)     
  1246)     select case (cur_region%def_type)
  1247)     !--------------------------------
  1248)       case (DEFINED_BY_BLOCK)
  1249)         write(id,'(a29)',advance='no') 'defined by: '
  1250)         write(id,'(a)') 'BLOCK'
  1251)         write(id,'(a29)',advance='no') 'I indices: '
  1252)         write(word1,iFormat) cur_region%i1
  1253)         write(word2,iFormat) cur_region%i2
  1254)         write(id,'(a)') adjustl(trim(word1)) // ' ' // adjustl(trim(word2))
  1255)         write(id,'(a29)',advance='no') 'J indices: '
  1256)         write(word1,iFormat) cur_region%j1
  1257)         write(word2,iFormat) cur_region%j2
  1258)         write(id,'(a)') adjustl(trim(word1)) // ' ' // adjustl(trim(word2))
  1259)         write(id,'(a29)',advance='no') 'K indices: '
  1260)         write(word1,iFormat) cur_region%k1
  1261)         write(word2,iFormat) cur_region%k2
  1262)         write(id,'(a)') adjustl(trim(word1)) // ' ' // adjustl(trim(word2))
  1263)     !--------------------------------
  1264)       case (DEFINED_BY_CARTESIAN_BOUNDARY)
  1265)         write(id,'(a29)',advance='no') 'defined by: '
  1266)         write(id,'(a)') 'CARTESIAN BOUNDARY'
  1267)     !--------------------------------
  1268)       case (DEFINED_BY_COORD)
  1269)         write(id,'(a29)',advance='no') 'defined by: '
  1270)         write(id,'(a)') 'COORDINATE(S)'
  1271)         write(id,'(a29)',advance='no') 'X coordinate(s): '
  1272)         string = ''
  1273)         do k = 1,size(cur_region%coordinates)
  1274)          write(word1,sFormat) cur_region%coordinates(k)%x
  1275)          string = adjustl(trim(string)) // ' ' // adjustl(trim(word1))
  1276)         enddo
  1277)         write(id,'(a)') adjustl(trim(string)) // ' m'
  1278)         write(id,'(a29)',advance='no') 'Y coordinate(s): '
  1279)         string = ''
  1280)         do k = 1,size(cur_region%coordinates)
  1281)           write(word1,sFormat) cur_region%coordinates(k)%y
  1282)           string = adjustl(trim(string)) // ' ' // adjustl(trim(word1))
  1283)         enddo
  1284)         write(id,'(a)') adjustl(trim(string)) // ' m'
  1285)         write(id,'(a29)',advance='no') 'Z coordinate(s): '
  1286)         string = ''
  1287)         do k = 1,size(cur_region%coordinates)
  1288)           write(word1,sFormat) cur_region%coordinates(k)%z
  1289)           string = adjustl(trim(string)) // ' ' // adjustl(trim(word1))
  1290)         enddo
  1291)         write(id,'(a)') adjustl(trim(string)) // ' m'
  1292)     !--------------------------------
  1293)       case (DEFINED_BY_CELL_AND_FACE_IDS)
  1294)         write(id,'(a29)',advance='no') 'defined by: '
  1295)         write(id,'(a)') 'CELL AND FACE IDS'
  1296)     !--------------------------------
  1297)       case (DEFINED_BY_CELL_IDS)
  1298)         write(id,'(a29)',advance='no') 'defined by: '
  1299)         write(id,'(a)') 'CELL IDS'
  1300)     !--------------------------------
  1301)       case (DEFINED_BY_VERTEX_IDS)
  1302)         write(id,'(a29)',advance='no') 'defined by: '
  1303)         write(id,'(a)') 'VERTEX IDS'
  1304)     !--------------------------------
  1305)       case (DEFINED_BY_FACE_UGRID_EXP)
  1306)         write(id,'(a29)',advance='no') 'defined by: '
  1307)         write(id,'(a)') 'FACE UNSTRUCTURED GRID EXPLICIT'
  1308)     !--------------------------------
  1309)       case (DEFINED_BY_POLY_BOUNDARY_FACE)
  1310)         write(id,'(a29)',advance='no') 'defined by: '
  1311)         write(id,'(a)') 'POLYGON BOUNDARY FACES IN VOLUME'
  1312)     !--------------------------------
  1313)       case (DEFINED_BY_POLY_CELL_CENTER)
  1314)         write(id,'(a29)',advance='no') 'defined by: '
  1315)         write(id,'(a)') 'POLYGON CELL CENTERS IN VOLUME'
  1316)     !--------------------------------
  1317)     end select
  1318)     
  1319)     if (cur_region%iface /= 0) then
  1320)       write(id,'(a29)',advance='no') 'face: '
  1321)       select case (cur_region%iface)
  1322)         case (WEST_FACE)
  1323)           write(id,'(a)') 'west'
  1324)         case (EAST_FACE)
  1325)           write(id,'(a)') 'east'
  1326)         case (NORTH_FACE)
  1327)           write(id,'(a)') 'north'
  1328)         case (SOUTH_FACE)
  1329)           write(id,'(a)') 'south'
  1330)         case (BOTTOM_FACE)
  1331)           write(id,'(a)') 'bottom'
  1332)         case (TOP_FACE)
  1333)           write(id,'(a)') 'top'
  1334)       end select
  1335)     endif
  1336)     
  1337)     write(id,'(a29)') '---------------------------: '
  1338)     cur_region => cur_region%next
  1339)   enddo
  1340)   
  1341) end subroutine RegionInputRecord
  1342) 
  1343) ! **************************************************************************** !
  1344) 
  1345) subroutine RegionDestroySideset(sideset)
  1346)   ! 
  1347)   ! Deallocates a unstructured grid side set
  1348)   ! 
  1349)   ! Author: Glenn Hammond
  1350)   ! Date: 11/01/09
  1351)   ! 
  1352) 
  1353)   implicit none
  1354)   
  1355)   type(region_sideset_type), pointer :: sideset
  1356)   
  1357)   if (.not.associated(sideset)) return
  1358)   
  1359)   if (associated(sideset%face_vertices)) deallocate(sideset%face_vertices)
  1360)   nullify(sideset%face_vertices)
  1361)   
  1362)   deallocate(sideset)
  1363)   nullify(sideset)
  1364)   
  1365) end subroutine RegionDestroySideset
  1366) 
  1367) ! ************************************************************************** !
  1368) 
  1369) subroutine RegionDestroyExplicitFaceSet(explicit_faceset)
  1370)   ! 
  1371)   ! Deallocates a unstructured grid explicit grid
  1372)   ! 
  1373)   ! Author: Glenn Hammond
  1374)   ! Date: 05/18/12
  1375)   ! 
  1376) 
  1377)   use Utility_module, only : DeallocateArray
  1378) 
  1379)   implicit none
  1380)   
  1381)   type(region_explicit_face_type), pointer :: explicit_faceset
  1382)   
  1383)   if (.not.associated(explicit_faceset)) return
  1384)   
  1385)   if (associated(explicit_faceset%face_centroids)) &
  1386)     deallocate(explicit_faceset%face_centroids)
  1387)   nullify(explicit_faceset%face_centroids)
  1388)   call DeallocateArray(explicit_faceset%face_areas)
  1389)   
  1390)   deallocate(explicit_faceset)
  1391)   nullify(explicit_faceset)
  1392)   
  1393) end subroutine RegionDestroyExplicitFaceSet
  1394) 
  1395) ! ************************************************************************** !
  1396) 
  1397) subroutine RegionDestroyList(region_list)
  1398)   ! 
  1399)   ! Deallocates a list of regions
  1400)   ! 
  1401)   ! Author: Glenn Hammond
  1402)   ! Date: 11/01/07
  1403)   ! 
  1404) 
  1405)   implicit none
  1406)   
  1407)   type(region_list_type), pointer :: region_list
  1408)   
  1409)   type(region_type), pointer :: region, prev_region
  1410)   
  1411)   if (.not.associated(region_list)) return
  1412)   
  1413)   region => region_list%first
  1414)   do 
  1415)     if (.not.associated(region)) exit
  1416)     prev_region => region
  1417)     region => region%next
  1418)     call RegionDestroy(prev_region)
  1419)   enddo
  1420)   
  1421)   region_list%num_regions = 0
  1422)   nullify(region_list%first)
  1423)   nullify(region_list%last)
  1424)   if (associated(region_list%array)) deallocate(region_list%array)
  1425)   nullify(region_list%array)
  1426)   
  1427)   deallocate(region_list)
  1428)   nullify(region_list)
  1429) 
  1430) end subroutine RegionDestroyList
  1431) 
  1432) ! ************************************************************************** !
  1433) 
  1434) subroutine RegionDestroy(region)
  1435)   ! 
  1436)   ! Deallocates a region
  1437)   ! 
  1438)   ! Author: Glenn Hammond
  1439)   ! Date: 10/23/07
  1440)   ! 
  1441)   use Utility_module, only : DeallocateArray
  1442)   
  1443)   implicit none
  1444)   
  1445)   type(region_type), pointer :: region
  1446)   
  1447)   if (.not.associated(region)) return
  1448)   
  1449)   call DeallocateArray(region%cell_ids)
  1450)   call DeallocateArray(region%faces)
  1451)   if (associated(region%coordinates)) deallocate(region%coordinates)
  1452)   nullify(region%coordinates)
  1453)   call RegionDestroySideset(region%sideset)
  1454)   call RegionDestroyExplicitFaceSet(region%explicit_faceset)
  1455)   call GeometryDestroyPolygonalVolume(region%polygonal_volume)
  1456)   
  1457)   if (associated(region%vertex_ids)) deallocate(region%vertex_ids)
  1458)   nullify(region%vertex_ids)
  1459)   
  1460)   nullify(region%next)
  1461) 
  1462)   deallocate(region)
  1463)   nullify(region)
  1464) 
  1465) end subroutine RegionDestroy
  1466) 
  1467) end module Region_module

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody