fluid.F90       coverage:  100.00 %func     89.71 %block


     1) module Fluid_module
     2)  
     3)   use PFLOTRAN_Constants_module
     4) 
     5)   implicit none
     6) 
     7)   private
     8) 
     9) #include "petsc/finclude/petscsys.h"
    10)  
    11)   type, public :: fluid_property_type
    12)     PetscReal :: tort_bin_diff
    13)     PetscReal :: vap_air_diff_coef
    14)     PetscReal :: exp_binary_diff
    15)     PetscReal :: enh_binary_diff_coef
    16)     PetscReal :: diff_base
    17)     PetscReal :: diff_exp
    18)     character(len=MAXWORDLENGTH) :: phase_name
    19)     PetscInt :: phase_id
    20)     PetscReal :: diffusion_coefficient
    21)     PetscReal :: gas_diffusion_coefficient
    22)     PetscReal :: diffusion_activation_energy
    23)     PetscReal :: nacl_concentration
    24)     type(fluid_property_type), pointer :: next
    25)   end type fluid_property_type
    26)   
    27)   public :: FluidPropertyCreate, FluidPropertyDestroy, &
    28)             FluidPropertyRead, FluidPropertyAddToList
    29) 
    30) contains
    31) 
    32) ! ************************************************************************** !
    33) 
    34) function FluidPropertyCreate()
    35)   ! 
    36)   ! Creates a fluid property object
    37)   ! 
    38)   ! Author: Glenn Hammond
    39)   ! Date: 01/21/09
    40)   ! 
    41)   
    42)   implicit none
    43) 
    44)   type(fluid_property_type), pointer :: FluidPropertyCreate
    45)   
    46)   type(fluid_property_type), pointer :: fluid_property
    47)   
    48)   allocate(fluid_property)
    49)   fluid_property%tort_bin_diff = 0.d0
    50)   fluid_property%vap_air_diff_coef = 2.13d-5
    51)   fluid_property%exp_binary_diff = 0.d0
    52)   fluid_property%enh_binary_diff_coef = 0.d0
    53)   fluid_property%diff_base = 0.d0
    54)   fluid_property%diff_exp = 0.d0
    55)   fluid_property%phase_name = 'LIQUID'
    56)   fluid_property%phase_id = 0
    57)   fluid_property%diffusion_coefficient = 1.d-9
    58)   fluid_property%gas_diffusion_coefficient = 2.13D-5
    59)   fluid_property%diffusion_activation_energy = 0.d0
    60)   fluid_property%nacl_concentration = 0.d0
    61)   nullify(fluid_property%next)
    62)   FluidPropertyCreate => fluid_property
    63) 
    64) end function FluidPropertyCreate
    65) 
    66) ! ************************************************************************** !
    67) 
    68) subroutine FluidPropertyRead(fluid_property,input,option)
    69)   ! 
    70)   ! Reads in contents of a fluid property card
    71)   ! 
    72)   ! Author: Glenn Hammond
    73)   ! Date: 01/21/09
    74)   ! 
    75) 
    76)   use Option_module
    77)   use Input_Aux_module
    78)   use String_module
    79)   use Units_module
    80) 
    81)   implicit none
    82)   
    83)   type(fluid_property_type) :: fluid_property
    84)   type(input_type), pointer :: input
    85)   type(option_type) :: option
    86)   
    87)   character(len=MAXWORDLENGTH) :: keyword, word
    88)   character(len=MAXWORDLENGTH) :: internal_units
    89) 
    90)   input%ierr = 0
    91)   do
    92)   
    93)     call InputReadPflotranString(input,option)
    94) 
    95)     if (InputCheckExit(input,option)) exit  
    96) 
    97)     call InputReadWord(input,option,keyword,PETSC_TRUE)
    98)     call InputErrorMsg(input,option,'keyword','FLUID_PROPERTY')
    99)     call StringToUpper(keyword)   
   100)       
   101)     select case(trim(keyword))
   102)     
   103)       case('PHASE') 
   104)         call InputReadWord(input,option,fluid_property%phase_name,PETSC_TRUE)
   105)         call InputErrorMsg(input,option,'phase','FLUID_PROPERTY')
   106)       case('DIFFUSION_COEFFICIENT','LIQUID_DIFFUSION_COEFFICIENT') 
   107)         call InputReadDouble(input,option,fluid_property%diffusion_coefficient)
   108)         call InputErrorMsg(input,option,'diffusion coefficient', &
   109)                            'FLUID_PROPERTY')
   110)         call InputReadAndConvertUnits(input, &
   111)                                       fluid_property%diffusion_coefficient, &
   112)                          'm^2/sec','FLUID_PROPERTY,diffusion_coeffient',option)
   113)       case('DIFFUSION_ACTIVATION_ENERGY') 
   114)         call InputReadDouble(input,option, &
   115)                              fluid_property%diffusion_activation_energy)
   116)         call InputErrorMsg(input,option,'diffusion activation energy', &
   117)                            'FLUID_PROPERTY')
   118)       case('GAS_DIFFUSION_COEFFICIENT') 
   119)         call InputReadDouble(input,option, &
   120)                              fluid_property%gas_diffusion_coefficient)
   121)         call InputErrorMsg(input,option,'gas diffusion coefficient', &
   122)                            'FLUID_PROPERTY')
   123)       case default
   124)         call InputKeywordUnrecognized(keyword,'FLUID_PROPERTY',option)
   125)     end select
   126)     
   127)   enddo  
   128) 
   129)   if (.not.(StringCompareIgnoreCase(fluid_property%phase_name,'LIQUID') .or. &
   130)             StringCompareIgnoreCase(fluid_property%phase_name,'GAS'))) then
   131)     option%io_buffer = 'PHASE in FLUID_PROPERTY should be LIQUID or GAS.'
   132)     call printErrMsg(option)
   133)   endif
   134) 
   135) 
   136) end subroutine FluidPropertyRead
   137) 
   138) ! ************************************************************************** !
   139) 
   140) subroutine FluidPropertyAddToList(fluid_property,list)
   141)   ! 
   142)   ! Adds a thermal property to linked list
   143)   ! 
   144)   ! Author: Glenn Hammond
   145)   ! Date: 01/21/09
   146)   ! 
   147) 
   148)   implicit none
   149)   
   150)   type(fluid_property_type), pointer :: fluid_property
   151)   type(fluid_property_type), pointer :: list
   152) 
   153)   type(fluid_property_type), pointer :: cur_fluid_property
   154)   
   155)   if (associated(list)) then
   156)     cur_fluid_property => list
   157)     ! loop to end of list
   158)     do
   159)       if (.not.associated(cur_fluid_property%next)) exit
   160)       cur_fluid_property => cur_fluid_property%next
   161)     enddo
   162)     cur_fluid_property%next => fluid_property
   163)   else
   164)     list => fluid_property
   165)   endif
   166)   
   167) end subroutine FluidPropertyAddToList
   168) 
   169) ! ************************************************************************** !
   170) 
   171) recursive subroutine FluidPropertyDestroy(fluid_property)
   172)   ! 
   173)   ! Destroys a fluid property
   174)   ! 
   175)   ! Author: Glenn Hammond
   176)   ! Date: 01/21/09
   177)   ! 
   178) 
   179)   implicit none
   180)   
   181)   type(fluid_property_type), pointer :: fluid_property
   182)   
   183)   if (.not.associated(fluid_property)) return
   184)   
   185)   call FluidPropertyDestroy(fluid_property%next)
   186) 
   187)   deallocate(fluid_property)
   188)   nullify(fluid_property)
   189)   
   190) end subroutine FluidPropertyDestroy
   191) 
   192) end module Fluid_module

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