geomechanics_region.F90       coverage:  90.91 %func     75.79 %block


     1) module Geomechanics_Region_module
     2)  
     3)   use Geometry_module
     4)   use PFLOTRAN_Constants_module
     5)   
     6)   implicit none
     7) 
     8)   private
     9) 
    10) #include "petsc/finclude/petscsys.h"
    11)  
    12)   type, public :: gm_region_type
    13)     PetscInt :: id
    14)     character(len=MAXWORDLENGTH) :: name
    15)     character(len=MAXSTRINGLENGTH) :: filename
    16)     type(point3d_type), pointer :: coordinates(:)
    17)     PetscInt :: num_verts 
    18)     PetscInt, pointer :: vertex_ids(:) 
    19)     type(gm_region_type), pointer :: next
    20)   end type gm_region_type
    21)   
    22)   type, public :: gm_region_ptr_type
    23)     type(gm_region_type), pointer :: ptr
    24)   end type gm_region_ptr_type
    25)   
    26)   type, public :: gm_region_list_type
    27)     PetscInt :: num_regions
    28)     type(gm_region_type), pointer :: first
    29)     type(gm_region_type), pointer :: last
    30)     type(gm_region_type), pointer :: array(:)
    31)   end type gm_region_list_type
    32) 
    33)   interface GeomechRegionCreate
    34)     module procedure GeomechRegionCreateWithList
    35)     module procedure GeomechRegionCreateWithNothing
    36)     module procedure GeomechRegionCreateWithGeomechRegion
    37)   end interface GeomechRegionCreate
    38)   
    39)   interface GeomechRegionReadFromFile
    40)     module procedure GeomechRegionReadFromFileId
    41)     module procedure GeomechRegionReadFromFilename
    42)   end interface GeomechRegionReadFromFile
    43)   
    44)    public :: GeomechRegionCreate, &
    45)              GeomechRegionDestroy, &
    46)              GeomechRegionAddToList, &
    47)              GeomechRegionReadFromFile, &
    48)              GeomechRegionDestroyList, &
    49)              GeomechRegionRead, &
    50)              GeomechRegionInitList, &
    51)              GeomechRegionGetPtrFromList
    52)              
    53)  contains
    54) 
    55) ! ************************************************************************** !
    56) 
    57) function GeomechRegionCreateWithNothing()
    58)   ! 
    59)   ! Creates a region with no arguments for
    60)   ! geomechanics
    61)   ! 
    62)   ! Author: Satish Karra, LANL
    63)   ! Date: 06/06/13
    64)   ! 
    65) 
    66)   implicit none
    67)   
    68)   type(gm_region_type), pointer :: GeomechRegionCreateWithNothing
    69)   
    70)   type(gm_region_type), pointer :: region
    71)   
    72)   allocate(region)
    73)   region%id = 0
    74)   region%name = ""
    75)   region%filename = ""
    76)   region%num_verts = 0
    77)   nullify(region%coordinates)
    78)   nullify(region%vertex_ids)
    79)   nullify(region%next)
    80)   
    81)   GeomechRegionCreateWithNothing => region
    82) 
    83) end function GeomechRegionCreateWithNothing
    84) 
    85) ! ************************************************************************** !
    86) 
    87) function GeomechRegionCreateWithList(list)
    88)   ! 
    89)   ! Creates a region from a list of vertices
    90)   ! 
    91)   ! Author: Satish Karra, LANL
    92)   ! Date: 06/06/13
    93)   ! 
    94) 
    95)   implicit none
    96)   
    97)   PetscInt :: list(:)
    98)   
    99)   type(gm_region_type), pointer :: GeomechRegionCreateWithList
   100)   
   101)   type(gm_region_type), pointer :: region
   102) 
   103)   region => GeomechRegionCreateWithNothing()
   104)   region%num_verts = size(list)
   105)   allocate(region%vertex_ids(region%num_verts))
   106)   region%vertex_ids = list
   107)   
   108)   GeomechRegionCreateWithList => region
   109) 
   110) end function GeomechRegionCreateWithList
   111) 
   112) ! ************************************************************************** !
   113) 
   114) function GeomechRegionCreateWithGeomechRegion(region)
   115)   ! 
   116)   ! Creates a copy of a region
   117)   ! 
   118)   ! Author: Satish Karra, LANL
   119)   ! Date: 06/06/13
   120)   ! 
   121) 
   122)   use Grid_Unstructured_Cell_module
   123) 
   124)   implicit none
   125)   
   126)   type(gm_region_type), pointer :: GeomechRegionCreateWithGeomechRegion
   127)   type(gm_region_type), pointer :: region
   128)   
   129)   type(gm_region_type), pointer :: new_region
   130)   PetscInt :: icount, temp_int
   131)   
   132)   new_region => GeomechRegionCreateWithNothing()
   133)   
   134)   new_region%id = region%id
   135)   new_region%name = region%name
   136)   new_region%filename = region%filename
   137)   new_region%num_verts = region%num_verts
   138)   if (associated(region%coordinates)) then
   139)     call GeometryCopyCoordinates(region%coordinates, &
   140)                                  new_region%coordinates)
   141)   endif
   142)   if (associated(region%vertex_ids)) then
   143)     allocate(new_region%vertex_ids(new_region%num_verts))
   144)     new_region%vertex_ids(1:new_region%num_verts) = &
   145)     region%vertex_ids(1:new_region%num_verts)
   146)   endif
   147)  
   148)   GeomechRegionCreateWithGeomechRegion => new_region
   149)   
   150) end function GeomechRegionCreateWithGeomechRegion
   151) 
   152) ! ************************************************************************** !
   153) 
   154) subroutine GeomechRegionInitList(list)
   155)   ! 
   156)   ! Initializes a region list
   157)   ! 
   158)   ! Author: Satish Karra, LANL
   159)   ! Date: 06/06/13
   160)   ! 
   161) 
   162)   implicit none
   163) 
   164)   type(gm_region_list_type) :: list
   165)   
   166)   nullify(list%first)
   167)   nullify(list%last)
   168)   nullify(list%array)
   169)   list%num_regions = 0
   170) 
   171) end subroutine GeomechRegionInitList
   172) 
   173) ! ************************************************************************** !
   174) 
   175) subroutine GeomechRegionAddToList(new_region,list)
   176)   ! 
   177)   ! Adds a new region to a region list
   178)   ! 
   179)   ! Author: Satish Karra, LANL
   180)   ! Date: 06/06/13
   181)   ! 
   182) 
   183)   implicit none
   184)   
   185)   type(gm_region_type), pointer :: new_region
   186)   type(gm_region_list_type) :: list
   187)   
   188)   list%num_regions = list%num_regions + 1
   189)   new_region%id = list%num_regions
   190)   if (.not.associated(list%first)) list%first => new_region
   191)   if (associated(list%last)) list%last%next => new_region
   192)   list%last => new_region
   193)   
   194) end subroutine GeomechRegionAddToList
   195) 
   196) ! ************************************************************************** !
   197) 
   198) subroutine GeomechRegionRead(region,input,option)
   199)   ! 
   200)   ! Reads a region from the input file
   201)   ! 
   202)   ! Author: Satish Karra, LANL
   203)   ! Date: 06/06/13
   204)   ! 
   205) 
   206)   use Input_Aux_module
   207)   use String_module
   208)   use Option_module
   209)   
   210)   implicit none
   211)   
   212)   type(option_type) :: option
   213)   type(gm_region_type) :: region
   214)   type(input_type), pointer :: input
   215)   
   216)   character(len=MAXWORDLENGTH) :: keyword, word
   217)  
   218)   input%ierr = 0
   219)   do
   220)   
   221)     call InputReadPflotranString(input,option)
   222)     if (InputError(input)) exit
   223)     if (InputCheckExit(input,option)) exit
   224)     
   225)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   226)     call InputErrorMsg(input,option,'keyword','GEOMECHANICS_REGION')
   227)     call StringToUpper(keyword)   
   228) 
   229)     select case(trim(keyword))
   230)       case('COORDINATE')
   231)         allocate(region%coordinates(1))
   232)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%x) 
   233)         if (InputError(input)) then
   234)           input%ierr = 0
   235)           call InputReadPflotranString(input,option)
   236)           call InputReadStringErrorMsg(input,option,'GEOMECHANICS_REGION')
   237)           call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%x)
   238)         endif
   239)         call InputErrorMsg(input,option,'x-coordinate','GEOMECHANICS_REGION')
   240)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%y)
   241)         call InputErrorMsg(input,option,'y-coordinate','GEOMECHANICS_REGION')
   242)         call InputReadDouble(input,option,region%coordinates(ONE_INTEGER)%z)
   243)         call InputErrorMsg(input,option,'z-coordinate','GEOMECHANICS_REGION')
   244)       case('COORDINATES')
   245)         call GeometryReadCoordinates(input,option,region%name, &
   246)                                      region%coordinates)
   247)       case('FILE')
   248)         call InputReadNChars(input,option,region%filename, &
   249)             MAXSTRINGLENGTH,PETSC_TRUE)
   250)         call InputErrorMsg(input,option,'filename','GEOMECHANICS_REGION')
   251)         call GeomechRegionReadFromFilename(region,option,region%filename)
   252)       case('LIST')
   253)         option%io_buffer = 'GEOMECHANICS_REGION LIST currently not implemented'
   254)         call printErrMsg(option)
   255)       case default
   256)         call InputKeywordUnrecognized(keyword,'GEOMECHANICS_REGION',option)
   257)     end select
   258)   enddo
   259)  
   260) end subroutine GeomechRegionRead
   261) 
   262) ! ************************************************************************** !
   263) 
   264) subroutine GeomechRegionReadFromFilename(region,option,filename)
   265)   ! 
   266)   ! Reads a list of vertex ids from a file named
   267)   ! filename
   268)   ! 
   269)   ! Author: Satish Karra, LANL
   270)   ! Date: 06/06/13
   271)   ! 
   272) 
   273)   use Input_Aux_module
   274)   use Option_module
   275)   use Utility_module
   276)   
   277)   implicit none
   278)   
   279)   type(gm_region_type) :: region
   280)   type(option_type) :: option
   281)   type(input_type), pointer :: input
   282)   character(len=MAXSTRINGLENGTH) :: filename
   283)   
   284)   input => InputCreate(IUNIT_TEMP,filename,option)
   285)   call GeomechRegionReadFromFileId(region,input,option)          
   286)   call InputDestroy(input)         
   287) 
   288) end subroutine GeomechRegionReadFromFilename
   289) 
   290) ! ************************************************************************** !
   291) 
   292) subroutine GeomechRegionReadFromFileId(region,input,option)
   293)   ! 
   294)   ! Reads a list of vertex ids from an open file
   295)   ! 
   296)   ! Author: Satish Karra, LANL
   297)   ! Date: 06/06/13
   298)   ! 
   299) 
   300)   use Input_Aux_module
   301)   use Option_module
   302)   use Utility_module
   303)   use Logging_module
   304)   use Grid_Unstructured_Cell_module
   305)   
   306)   implicit none
   307)   
   308)   type(gm_region_type) :: region
   309)   type(option_type) :: option
   310)   type(input_type), pointer :: input
   311)   
   312)   character(len=MAXWORDLENGTH) :: word
   313)   character(len=1) :: backslash
   314)   character(len=MAXSTRINGLENGTH) :: string, string1
   315) 
   316)   PetscInt, pointer :: temp_int_array(:)
   317)   PetscInt, pointer :: vertex_ids(:)
   318)   PetscInt :: max_size
   319)   PetscInt :: count
   320)   PetscInt :: temp_int
   321)   PetscInt :: input_data_type
   322)   PetscInt :: ii
   323)   PetscInt :: istart
   324)   PetscInt :: iend
   325)   PetscInt :: remainder
   326)   PetscErrorCode :: ierr
   327) 
   328)   max_size = 1000
   329)   backslash = achar(92)  ! 92 = "\" Some compilers choke on \" thinking it
   330)                           ! is a double quote as in c/c++
   331)   
   332)   allocate(temp_int_array(max_size))
   333)   allocate(vertex_ids(max_size))
   334)   
   335)   temp_int_array = 0
   336)   vertex_ids = 0
   337)   
   338)   count = 0
   339)   call InputReadPflotranString(input, option)
   340)   do 
   341)     call InputReadInt(input, option, temp_int)
   342)     if (InputError(input)) exit
   343)     count = count + 1
   344)     temp_int_array(count) = temp_int
   345)   enddo
   346) 
   347)   if (count == 1) then
   348)     !
   349)     ! Input data contains only cell ids
   350)     !
   351)     vertex_ids(1) = temp_int_array(1)
   352)     count = 1
   353) 
   354)     ! Read the data
   355)     do
   356)       call InputReadPflotranString(input, option)
   357)       if (InputError(input)) exit
   358)       call InputReadInt(input, option, temp_int)
   359)       if (.not.InputError(input)) then
   360)         count = count + 1
   361)         vertex_ids(count) = temp_int
   362)       endif
   363)       if (count+1 > max_size) then ! resize temporary array
   364)         call reallocateIntArray(vertex_ids, max_size)
   365)       endif
   366)     enddo
   367) 
   368)     ! Depending on processor rank, save only a portion of data
   369)     region%num_verts = count/option%mycommsize
   370)       remainder = count - region%num_verts*option%mycommsize
   371)     if (option%myrank < remainder) region%num_verts = region%num_verts + 1
   372)     istart = 0
   373)     iend   = 0
   374)     call MPI_Exscan(region%num_verts, istart, ONE_INTEGER_MPI, MPIU_INTEGER, &
   375)                     MPI_SUM, option%mycomm, ierr)
   376)     call MPI_Scan(region%num_verts, iend, ONE_INTEGER_MPI, MPIU_INTEGER, &
   377)                    MPI_SUM, option%mycomm, ierr)
   378) 
   379)     ! Allocate memory and save the data
   380)     region%num_verts = iend - istart
   381)     allocate(region%vertex_ids(region%num_verts))
   382)     region%vertex_ids(1:region%num_verts) = vertex_ids(istart+1:iend)
   383)     deallocate(vertex_ids)
   384)   else
   385)    option%io_buffer = 'Provide one vertex_id per line, GEOMECHANICS_REGION.'
   386)    call printErrMsg(option) 
   387)   endif
   388)   
   389)   deallocate(temp_int_array)
   390)   
   391) #ifdef GEOMECH_DEBUG
   392)   write(string,*) option%myrank
   393)   write(string1,*) region%name
   394)   string = 'geomech_region_' // trim(adjustl(string1)) // '_vertex_ids' &
   395)     // trim(adjustl(string)) // '.out'
   396)   open(unit=86,file=trim(string))
   397)   do ii = 1, region%num_verts
   398)     write(86,'(i5)') region%vertex_ids(ii)
   399)   enddo  
   400)   close(86)
   401) #endif    
   402)     
   403) end subroutine GeomechRegionReadFromFileId
   404) 
   405) ! ************************************************************************** !
   406) 
   407) function GeomechRegionGetPtrFromList(region_name,region_list)
   408)   ! 
   409)   ! Returns a pointer to the region matching region_name
   410)   ! 
   411)   ! Author: Satish Karra, LANL
   412)   ! Date: 06/06/13
   413)   ! 
   414) 
   415)   use String_module
   416) 
   417)   implicit none
   418)   
   419)   type(gm_region_type), pointer :: GeomechRegionGetPtrFromList
   420)   character(len=MAXWORDLENGTH) :: region_name
   421)   PetscInt :: length
   422)   type(gm_region_list_type) :: region_list
   423) 
   424)   type(gm_region_type), pointer :: region
   425)     
   426)   nullify(GeomechRegionGetPtrFromList)
   427)   region => region_list%first
   428)   
   429)   do 
   430)     if (.not.associated(region)) exit
   431)     length = len_trim(region_name)
   432)     if (length == len_trim(region%name) .and. &
   433)         StringCompare(region%name,region_name,length)) then
   434)       GeomechRegionGetPtrFromList => region
   435)       return
   436)     endif
   437)     region => region%next
   438)   enddo
   439)   
   440) end function GeomechRegionGetPtrFromList
   441) 
   442) ! ************************************************************************** !
   443) 
   444) subroutine GeomechRegionDestroyList(region_list)
   445)   ! 
   446)   ! Deallocates a list of regions
   447)   ! 
   448)   ! Author: Satish Karra, LANL
   449)   ! Date: 06/06/13
   450)   ! 
   451) 
   452)   implicit none
   453)   
   454)   type(gm_region_list_type), pointer :: region_list
   455)   
   456)   type(gm_region_type), pointer :: region, prev_region
   457)   
   458)   if (.not.associated(region_list)) return
   459)   
   460)   region => region_list%first
   461)   do 
   462)     if (.not.associated(region)) exit
   463)     prev_region => region
   464)     region => region%next
   465)     call GeomechRegionDestroy(prev_region)
   466)   enddo
   467)   
   468)   region_list%num_regions = 0
   469)   nullify(region_list%first)
   470)   nullify(region_list%last)
   471)   if (associated(region_list%array)) deallocate(region_list%array)
   472)   nullify(region_list%array)
   473)   
   474)   deallocate(region_list)
   475)   nullify(region_list)
   476) 
   477) end subroutine GeomechRegionDestroyList
   478) 
   479) ! ************************************************************************** !
   480) 
   481) subroutine GeomechRegionDestroy(region)
   482)   ! 
   483)   ! Deallocates a region
   484)   ! 
   485)   ! Author: Satish Karra, LANL
   486)   ! Date: 06/06/13
   487)   ! 
   488) 
   489)   implicit none
   490)   
   491)   type(gm_region_type), pointer :: region
   492)   
   493)   if (.not.associated(region)) return
   494)   
   495)   if (associated(region%vertex_ids)) deallocate(region%vertex_ids)
   496)   nullify(region%vertex_ids)
   497)   if (associated(region%coordinates)) deallocate(region%coordinates)
   498)   nullify(region%coordinates)
   499)   
   500)   nullify(region%next)
   501) 
   502)   deallocate(region)
   503)   nullify(region)
   504)   
   505) end subroutine GeomechRegionDestroy
   506) 
   507) end module Geomechanics_Region_module

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