observation.F90       coverage:  77.78 %func     67.86 %block


     1) module Observation_module
     2) 
     3)   use Region_module
     4)   use Connection_module
     5)   
     6)   use PFLOTRAN_Constants_module
     7) 
     8)   implicit none
     9)   
    10)   private
    11)   
    12) #include "petsc/finclude/petscsys.h"
    13) 
    14)   PetscInt, parameter, public :: OBSERVATION_SCALAR = 1
    15)   PetscInt, parameter, public :: OBSERVATION_FLUX = 2
    16)   PetscInt, parameter, public :: OBSERVATION_AT_CELL_CENTER = 1
    17)   PetscInt, parameter, public :: OBSERVATION_AT_COORDINATE = 2
    18) 
    19)   type, public :: observation_type
    20)     ! all added variables must be included in ObservationCreateFromObservation
    21)     PetscInt :: id
    22)     PetscInt :: itype
    23)     PetscBool :: print_velocities
    24)     PetscBool :: print_secondary_data(5)          ! first entry is for temp., second is for conc. and third is for mineral vol frac., fourth is rate, fifth is SI
    25)     PetscBool :: at_cell_center
    26)     character(len=MAXWORDLENGTH) :: name
    27)     character(len=MAXWORDLENGTH) :: linkage_name
    28)     type(connection_set_type), pointer :: connection_set
    29)     type(region_type), pointer :: region
    30)     type(observation_type), pointer :: next
    31)   end type observation_type
    32)   
    33)   type, public :: observation_list_type
    34)     PetscInt :: num_observations
    35)     type(observation_type), pointer :: first
    36)     type(observation_type), pointer :: last
    37)     type(observation_type), pointer :: array(:)
    38)   end type observation_list_type
    39) 
    40)   public :: ObservationCreate, ObservationDestroy, ObservationRead, &
    41)             ObservationAddToList, ObservationInitList, ObservationDestroyList, &
    42)             ObservationGetPtrFromList, ObservationRemoveFromList
    43) 
    44)   interface ObservationCreate
    45)     module procedure ObservationCreate1
    46)     module procedure ObservationCreateFromObservation
    47)   end interface
    48)     
    49) contains
    50) 
    51) ! ************************************************************************** !
    52) 
    53) function ObservationCreate1()
    54)   ! 
    55)   ! Create object that stores observation regions
    56)   ! 
    57)   ! Author: Glenn Hammond
    58)   ! Date: 02/11/08
    59)   ! 
    60) 
    61)   implicit none
    62)   
    63)   type(observation_type), pointer :: ObservationCreate1
    64)   
    65)   type(observation_type), pointer :: observation
    66)   
    67)   allocate(observation)
    68)   
    69)   observation%name = ""
    70)   observation%linkage_name = ""
    71)   observation%id = 0
    72)   observation%itype = OBSERVATION_SCALAR
    73)   observation%print_velocities = PETSC_FALSE
    74)   observation%at_cell_center = PETSC_TRUE
    75)   observation%print_secondary_data = PETSC_FALSE
    76)   nullify(observation%region)
    77)   nullify(observation%next)
    78)   
    79)   ObservationCreate1 => observation
    80) 
    81) end function ObservationCreate1
    82) 
    83) ! ************************************************************************** !
    84) 
    85) function ObservationCreateFromObservation(observation)
    86)   ! 
    87)   ! Create object that stores observation regions
    88)   ! 
    89)   ! Author: Glenn Hammond
    90)   ! Date: 02/11/08
    91)   ! 
    92) 
    93)   implicit none
    94)   
    95)   type(observation_type), pointer :: ObservationCreateFromObservation
    96)   type(observation_type), pointer :: observation
    97) 
    98)   type(observation_type), pointer :: new_observation
    99)   
   100)   new_observation => ObservationCreate1()
   101)   
   102)   new_observation%name = observation%name
   103)   new_observation%linkage_name = observation%linkage_name
   104)   new_observation%id = observation%id
   105)   new_observation%itype = observation%itype
   106)   new_observation%print_velocities = observation%print_velocities
   107)   new_observation%at_cell_center = observation%at_cell_center
   108)   new_observation%print_secondary_data = &
   109)        observation%print_secondary_data
   110)   ! keep these null for now to catch bugs
   111)   nullify(new_observation%region)
   112)   nullify(new_observation%next)
   113)   
   114)   ObservationCreateFromObservation => new_observation
   115) 
   116) end function ObservationCreateFromObservation
   117) 
   118) ! ************************************************************************** !
   119) 
   120) subroutine ObservationRead(observation,input,option)
   121)   ! 
   122)   ! Reads observation data from the input file
   123)   ! 
   124)   ! Author: Glenn Hammond
   125)   ! Date: 02/11/08
   126)   ! 
   127) 
   128)   use Input_Aux_module
   129)   use String_module
   130)   use Option_module
   131)   
   132)   implicit none
   133)   
   134)   type(observation_type) :: observation
   135)   type(input_type), pointer :: input
   136)   type(option_type) :: option
   137)   
   138)   character(len=MAXWORDLENGTH) :: keyword
   139)   
   140)   input%ierr = 0
   141)   do
   142)   
   143)     call InputReadPflotranString(input,option)
   144)     
   145)     if (InputCheckExit(input,option)) exit  
   146) 
   147)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   148)     call InputErrorMsg(input,option,'keyword','OBSERVATION')   
   149)       
   150)     select case(trim(keyword))
   151)     
   152)       case('BOUNDARY_CONDITION')
   153)         call InputReadWord(input,option,observation%linkage_name,PETSC_TRUE)
   154)         call InputErrorMsg(input,option,'boundary condition name','OBSERVATION')
   155)         option%flow%store_fluxes = PETSC_TRUE
   156)         option%transport%store_fluxes = PETSC_TRUE
   157)         observation%itype = OBSERVATION_FLUX
   158)       case('REGION')
   159)         call InputReadWord(input,option,observation%linkage_name,PETSC_TRUE)
   160)         call InputErrorMsg(input,option,'region name','OBSERVATION')
   161)         observation%itype = OBSERVATION_SCALAR
   162)       case('VELOCITY')
   163)         observation%print_velocities = PETSC_TRUE
   164)       case('SECONDARY_TEMPERATURE')
   165)         if (option%use_mc) then
   166)           observation%print_secondary_data(1) = PETSC_TRUE
   167)         else
   168)           option%io_buffer = 'Keyword SECONDARY_TEMPERATURE can only be used' // &
   169)                              ' MULTIPLE_CONTINUUM keyword'
   170)           call printErrMsg(option)
   171)         endif
   172)       case('SECONDARY_CONCENTRATION')
   173)         if (option%use_mc) then
   174)           observation%print_secondary_data(2) = PETSC_TRUE
   175)         else
   176)           option%io_buffer = 'Keyword SECONDARY_CONCENTRATION can only be used' // &
   177)                              ' MULTIPLE_CONTINUUM keyword'
   178)           call printErrMsg(option)
   179)         endif
   180)       case('SECONDARY_MINERAL_VOLFRAC')
   181)         if (option%use_mc) then
   182)           observation%print_secondary_data(3) = PETSC_TRUE
   183)         else
   184)           option%io_buffer = 'Keyword SECONDARY_MINERAL_VOLFRAC can ' // &
   185)                              'only be used MULTIPLE_CONTINUUM keyword'
   186)           call printErrMsg(option)
   187)         endif
   188)       case('SECONDARY_MINERAL_RATE')
   189)         if (option%use_mc) then
   190)           observation%print_secondary_data(4) = PETSC_TRUE
   191)         else
   192)           option%io_buffer = 'Keyword SECONDARY_MINERAL_RATE can ' // &
   193)                              'only be used MULTIPLE_CONTINUUM keyword'
   194)           call printErrMsg(option)
   195)         endif
   196)       case('SECONDARY_MINERAL_SI')
   197)         if (option%use_mc) then
   198)           observation%print_secondary_data(5) = PETSC_TRUE
   199)         else
   200)           option%io_buffer = 'Keyword SECONDARY_MINERAL_SI can ' // &
   201)                              'only be used MULTIPLE_CONTINUUM keyword'
   202)           call printErrMsg(option)
   203)         endif
   204)       case('AT_CELL_CENTER')
   205)         observation%at_cell_center = PETSC_TRUE
   206)       case('AT_COORDINATE')
   207)         observation%at_cell_center = PETSC_FALSE
   208)       case default
   209)         call InputKeywordUnrecognized(keyword,'OBSERVATION',option)
   210)     end select 
   211)   
   212)   enddo  
   213) 
   214) end subroutine ObservationRead
   215) 
   216) ! ************************************************************************** !
   217) 
   218) subroutine ObservationInitList(list)
   219)   ! 
   220)   ! Initializes a observation list
   221)   ! 
   222)   ! Author: Glenn Hammond
   223)   ! Date: 02/11/08
   224)   ! 
   225) 
   226)   implicit none
   227) 
   228)   type(observation_list_type) :: list
   229)   
   230)   nullify(list%first)
   231)   nullify(list%last)
   232)   nullify(list%array)
   233)   list%num_observations = 0
   234) 
   235) end subroutine ObservationInitList
   236) 
   237) ! ************************************************************************** !
   238) 
   239) subroutine ObservationAddToList(new_observation,list)
   240)   ! 
   241)   ! Adds a new observation to a observation list
   242)   ! 
   243)   ! Author: Glenn Hammond
   244)   ! Date: 02/11/08
   245)   ! 
   246) 
   247)   implicit none
   248)   
   249)   type(observation_type), pointer :: new_observation
   250)   type(observation_list_type) :: list
   251)   
   252)   list%num_observations = list%num_observations + 1
   253)   new_observation%id = list%num_observations
   254)   if (.not.associated(list%first)) list%first => new_observation
   255)   if (associated(list%last)) list%last%next => new_observation
   256)   list%last => new_observation
   257)   
   258) end subroutine ObservationAddToList
   259) 
   260) ! ************************************************************************** !
   261) 
   262) subroutine ObservationRemoveFromList(observation,list)
   263)   ! 
   264)   ! Removes a observation from a observation list
   265)   ! 
   266)   ! Author: Glenn Hammond
   267)   ! Date: 02/11/08
   268)   ! 
   269) 
   270)   implicit none
   271)   
   272)   type(observation_type), pointer :: observation
   273)   type(observation_list_type) :: list
   274)   
   275)   type(observation_type), pointer :: cur_observation, prev_observation
   276)   
   277)   cur_observation => list%first
   278)   nullify(prev_observation)
   279)   
   280)   do
   281)     if (.not.associated(cur_observation)) exit
   282)     if (associated(cur_observation,observation)) then
   283)       if (associated(prev_observation)) then
   284)         prev_observation%next => cur_observation%next
   285)       else
   286)         list%first => cur_observation%next
   287)       endif
   288)       if (.not.associated(cur_observation%next)) then
   289)         list%last => prev_observation
   290)       endif
   291)       list%num_observations = list%num_observations-1
   292)       call ObservationDestroy(cur_observation)
   293)       return
   294)     endif
   295)     prev_observation => cur_observation
   296)     cur_observation => cur_observation%next
   297)   enddo
   298)   
   299) end subroutine ObservationRemoveFromList
   300) 
   301) ! ************************************************************************** !
   302) 
   303) function ObservationGetPtrFromList(observation_name,observation_list)
   304)   ! 
   305)   ! Returns a pointer to the observation matching &
   306)   ! observation_name
   307)   ! 
   308)   ! Author: Glenn Hammond
   309)   ! Date: 02/11/08
   310)   ! 
   311) 
   312)   use String_module
   313) 
   314)   implicit none
   315)   
   316)   type(observation_type), pointer :: ObservationGetPtrFromList
   317)   character(len=MAXWORDLENGTH) :: observation_name
   318)   type(observation_list_type) :: observation_list
   319)  
   320)   PetscInt :: length
   321)   type(observation_type), pointer :: observation
   322)     
   323)   nullify(ObservationGetPtrFromList)
   324)   observation => observation_list%first
   325)   
   326)   do 
   327)     if (.not.associated(observation)) exit
   328)     length = len_trim(observation_name)
   329)     if (length == len_trim(observation%name) .and. &
   330)         StringCompare(observation%name,observation_name, &
   331)                         length)) then
   332)       ObservationGetPtrFromList => observation
   333)       return
   334)     endif
   335)     observation => observation%next
   336)   enddo
   337)   
   338) end function ObservationGetPtrFromList
   339) 
   340) ! ************************************************************************** !
   341) 
   342) subroutine ObservationDestroyList(observation_list)
   343)   ! 
   344)   ! Deallocates a list of observations
   345)   ! 
   346)   ! Author: Glenn Hammond
   347)   ! Date: 02/11/08
   348)   ! 
   349) 
   350)   implicit none
   351)   
   352)   type(observation_list_type), pointer :: observation_list
   353)   
   354)   type(observation_type), pointer :: observation, prev_observation
   355)   
   356)   if (.not.associated(observation_list)) return
   357)   
   358)   observation => observation_list%first
   359)   do 
   360)     if (.not.associated(observation)) exit
   361)     prev_observation => observation
   362)     observation => observation%next
   363)     call ObservationDestroy(prev_observation)
   364)   enddo
   365)   
   366)   observation_list%num_observations = 0
   367)   nullify(observation_list%first)
   368)   nullify(observation_list%last)
   369)   if (associated(observation_list%array)) deallocate(observation_list%array)
   370)   nullify(observation_list%array)
   371)   
   372)   deallocate(observation_list)
   373)   nullify(observation_list)
   374) 
   375) end subroutine ObservationDestroyList
   376) 
   377) ! ************************************************************************** !
   378) 
   379) subroutine ObservationDestroy(observation)
   380)   ! 
   381)   ! Deallocates a observation
   382)   ! 
   383)   ! Author: Glenn Hammond
   384)   ! Date: 10/23/07
   385)   ! 
   386) 
   387)   implicit none
   388)   
   389)   type(observation_type), pointer :: observation
   390)   
   391)   PetscInt :: i
   392)   
   393)   if (.not.associated(observation)) return
   394)   
   395)   nullify(observation%region)
   396)   nullify(observation%connection_set)
   397)   deallocate(observation)
   398)   nullify(observation)
   399) 
   400) end subroutine ObservationDestroy
   401) 
   402) end module Observation_module

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