srcsink_sandbox.F90       coverage:  91.67 %func     78.21 %block


     1) module SrcSink_Sandbox_module
     2) 
     3)   use SrcSink_Sandbox_Base_class
     4)   use SrcSink_Sandbox_WIPP_Gas_class
     5)   use SrcSink_Sandbox_Mass_Rate_class
     6)   use SrcSink_Sandbox_Downreg_class
     7)   use SrcSink_Sandbox_WIPP_Well_class
     8)   use PFLOTRAN_Constants_module
     9) 
    10)   implicit none
    11)   
    12)   private
    13)   
    14) #include "petsc/finclude/petscsys.h"
    15) 
    16)   class(srcsink_sandbox_base_type), pointer, public :: ss_sandbox_list
    17)   PetscBool :: print_mass_balance
    18) 
    19)   interface SSSandboxRead
    20)     module procedure SSSandboxRead1
    21)     module procedure SSSandboxRead2
    22)   end interface
    23)   
    24)   interface SSSandboxDestroyList
    25)     module procedure SSSandboxDestroyList1
    26)     module procedure SSSandboxDestroyList2
    27)   end interface
    28)   
    29)   public :: SSSandboxInit, &
    30)             SSSandboxRead, &
    31)             SSSandboxSetup, &
    32)             SSSandboxUpdate, &
    33)             SSSandbox, &
    34)             SSSandboxDestroyList
    35) 
    36) contains
    37) 
    38) ! ************************************************************************** !
    39) 
    40) subroutine SSSandboxInit(option)
    41)   ! 
    42)   ! Initializes the sandbox list
    43)   ! 
    44)   ! Author: Glenn Hammond
    45)   ! Date: 04/11/14
    46)   ! 
    47)   use Option_module
    48)   implicit none
    49)   type(option_type) :: option
    50) 
    51)   if (associated(ss_sandbox_list)) then
    52)     call SSSandboxDestroyList()
    53)   endif
    54)   nullify(ss_sandbox_list)
    55)   print_mass_balance = PETSC_FALSE
    56) 
    57) end subroutine SSSandboxInit
    58) 
    59) ! ************************************************************************** !
    60) 
    61) subroutine SSSandboxRead1(input,option)
    62)   ! 
    63)   ! Reads input deck for source/sink sandbox parameters
    64)   ! 
    65)   ! Author: Glenn Hammond
    66)   ! Date: 04/11/14
    67)   ! 
    68) 
    69)   use Option_module
    70)   use String_module
    71)   use Input_Aux_module
    72)   use Utility_module
    73)   
    74)   implicit none
    75)   
    76)   type(input_type), pointer :: input
    77)   type(option_type) :: option
    78) 
    79)   call SSSandboxRead(ss_sandbox_list,input,option)
    80) 
    81) end subroutine SSSandboxRead1
    82) 
    83) ! ************************************************************************** !
    84) 
    85) subroutine SSSandboxRead2(local_sandbox_list,input,option)
    86)   ! 
    87)   ! Reads input deck for src/sink sandbox parameters
    88)   ! 
    89)   ! Author: Glenn Hammond
    90)   ! Date: 04/11/14
    91)   ! 
    92) 
    93)   use Option_module
    94)   use String_module
    95)   use Input_Aux_module
    96)   use Utility_module
    97)   
    98)   implicit none
    99)   
   100)   class(srcsink_sandbox_base_type), pointer :: local_sandbox_list  
   101)   type(input_type), pointer :: input
   102)   type(option_type) :: option
   103) 
   104)   character(len=MAXSTRINGLENGTH) :: string
   105)   character(len=MAXWORDLENGTH) :: word
   106)   class(srcsink_sandbox_base_type), pointer :: new_sandbox, cur_sandbox
   107)   
   108)   nullify(new_sandbox)
   109)   do 
   110)     call InputReadPflotranString(input,option)
   111)     if (InputError(input)) exit
   112)     if (InputCheckExit(input,option)) exit
   113) 
   114)     call InputReadWord(input,option,word,PETSC_TRUE)
   115)     call InputErrorMsg(input,option,'keyword','SOURCE_SINK_SANDBOX')
   116)     call StringToUpper(word)   
   117) 
   118)     select case(trim(word))
   119)       case('WIPP-WELL')
   120)         new_sandbox => WIPPWellCreate()
   121)       case('WIPP-GAS_GENERATION')
   122)         new_sandbox => WIPPGasGenerationCreate()
   123)       case('MASS_RATE')
   124)         new_sandbox => MassRateCreate()
   125)       case('MASS_RATE_DOWNREGULATED')
   126)         new_sandbox => DownregCreate()
   127)       case('MASS_BALANCE')
   128)         print_mass_balance = PETSC_TRUE
   129)       case default
   130)         call InputKeywordUnrecognized(word,'SRCSINK_SANDBOX',option)
   131)     end select
   132)     
   133)     if (associated(new_sandbox)) then
   134)       call new_sandbox%ReadInput(input,option)
   135)       if (.not.associated(local_sandbox_list)) then
   136)         local_sandbox_list => new_sandbox
   137)       else
   138)         cur_sandbox => local_sandbox_list
   139)         do
   140)           if (.not.associated(cur_sandbox%next)) exit
   141)           cur_sandbox => cur_sandbox%next
   142)         enddo
   143)         cur_sandbox%next => new_sandbox
   144)       endif
   145)     endif
   146)     nullify(new_sandbox)
   147)   enddo
   148)   
   149) end subroutine SSSandboxRead2
   150) 
   151) 
   152) ! ************************************************************************** !
   153) 
   154) subroutine SSSandboxSetup(grid,option,output_option)
   155)   ! 
   156)   ! Calls all the initialization routines for all source/sinks in
   157)   ! the sandbox list
   158)   ! 
   159)   ! Author: Glenn Hammond
   160)   ! Date: 04/11/14
   161)   ! 
   162) 
   163)   use Option_module
   164)   use Output_Aux_module
   165)   use Grid_module
   166)   
   167)   implicit none
   168)   
   169)   type(grid_type) :: grid
   170)   type(option_type) :: option
   171)   type(output_option_type) :: output_option
   172)   
   173)   class(srcsink_sandbox_base_type), pointer :: cur_sandbox  
   174)   class(srcsink_sandbox_base_type), pointer :: prev_sandbox  
   175)   class(srcsink_sandbox_base_type), pointer :: next_sandbox  
   176)   PetscBool :: exists
   177) 
   178)   ! sandbox source/sinks
   179)   cur_sandbox => ss_sandbox_list
   180)   nullify(prev_sandbox)
   181)   do
   182)     if (.not.associated(cur_sandbox)) exit
   183)     next_sandbox => cur_sandbox%next
   184)     call cur_sandbox%Setup(grid,option)
   185)     ! destory if not on process
   186)     if (.not.Initialized(cur_sandbox%local_cell_id)) then
   187)       if (associated(prev_sandbox)) then
   188)         prev_sandbox%next => next_sandbox
   189)       else
   190)         ss_sandbox_list => next_sandbox
   191)       endif
   192)       nullify(cur_sandbox%next)
   193)       call SSSandboxDestroy(cur_sandbox)
   194)     else
   195)       if (print_mass_balance) then
   196)         allocate(cur_sandbox%instantaneous_mass_rate(option%nflowdof))
   197)         cur_sandbox%instantaneous_mass_rate = 0.d0
   198)         allocate(cur_sandbox%cumulative_mass(option%nflowdof))
   199)         cur_sandbox%cumulative_mass = 0.d0
   200)       endif    
   201)     endif
   202)     if (associated(cur_sandbox)) prev_sandbox => cur_sandbox
   203)     cur_sandbox => next_sandbox
   204)   enddo
   205)   
   206)   if (print_mass_balance) then
   207)     call SSSandboxOutputHeader(ss_sandbox_list,grid,option,output_option)
   208)   endif
   209) 
   210) end subroutine SSSandboxSetup
   211) 
   212) ! ************************************************************************** !
   213) 
   214) subroutine SSSandbox(residual,Jacobian,compute_derivative, &
   215)                      grid,material_auxvars,option)
   216)   ! 
   217)   ! Evaluates source/sink term storing residual and/or Jacobian
   218)   ! 
   219)   ! Author: Glenn Hammond
   220)   ! Date: 04/11/14
   221)   ! 
   222) 
   223)   use Option_module
   224)   use Grid_module
   225)   use Material_Aux_class, only: material_auxvar_type
   226)   
   227)   implicit none
   228)   
   229) #include "petsc/finclude/petscvec.h"
   230) #include "petsc/finclude/petscvec.h90"
   231) #include "petsc/finclude/petscmat.h"
   232) #include "petsc/finclude/petscmat.h90"
   233) 
   234)   PetscBool :: compute_derivative
   235)   Vec :: residual
   236)   Mat :: Jacobian
   237)   class(material_auxvar_type), pointer :: material_auxvars(:)
   238)   type(grid_type) :: grid
   239)   type(option_type) :: option
   240)   
   241)   PetscReal, pointer :: r_p(:)
   242)   PetscReal :: res(option%nflowdof)
   243)   PetscReal :: Jac(option%nflowdof,option%nflowdof)
   244)   class(srcsink_sandbox_base_type), pointer :: cur_srcsink
   245)   PetscInt :: i, local_id, ghosted_id, istart, iend
   246)   PetscReal :: aux_real(0)
   247)   PetscErrorCode :: ierr
   248)   
   249)   if (.not.compute_derivative) then
   250)     call VecGetArrayF90(residual,r_p,ierr);CHKERRQ(ierr)
   251)   endif
   252)   
   253)   cur_srcsink => ss_sandbox_list
   254)   do
   255)     if (.not.associated(cur_srcsink)) exit
   256)     local_id = cur_srcsink%local_cell_id
   257)     ghosted_id = grid%nL2G(local_id)
   258)     res = 0.d0
   259)     Jac = 0.d0
   260)     call cur_srcsink%Evaluate(res,Jac,compute_derivative, &
   261)                               material_auxvars(ghosted_id), &
   262)                               aux_real,option)
   263)     if (compute_derivative) then
   264)       call MatSetValuesBlockedLocal(Jacobian,1,ghosted_id-1,1, &
   265)                                     ghosted_id-1,Jac,ADD_VALUES, &
   266)                                     ierr);CHKERRQ(ierr)
   267)     else
   268)       iend = local_id*option%nflowdof
   269)       istart = iend - option%nflowdof + 1
   270)       r_p(istart:iend) = r_p(istart:iend) + res
   271)     endif
   272)     cur_srcsink => cur_srcsink%next
   273)   enddo
   274)   
   275)   if (.not.compute_derivative) then
   276)     call VecRestoreArrayF90(residual,r_p,ierr);CHKERRQ(ierr)
   277)   endif
   278) 
   279) end subroutine SSSandbox
   280) 
   281) ! ************************************************************************** !
   282) 
   283) subroutine SSSandboxUpdate(sandbox_list,time,option,output_option)
   284)   ! 
   285)   ! Updates datasets associated with a sandbox, if they exist
   286)   ! 
   287)   ! Author: Glenn Hammond
   288)   ! Date: 06/22/15
   289)   ! 
   290)   use Option_module
   291)   use Output_Aux_module
   292) 
   293)   implicit none
   294) 
   295)   class(srcsink_sandbox_base_type), pointer :: sandbox_list
   296)   PetscReal :: time
   297)   type(option_type) :: option
   298)   type(output_option_type) :: output_option
   299) 
   300)   class(srcsink_sandbox_base_type), pointer :: cur_sandbox
   301)   
   302)   cur_sandbox => sandbox_list
   303)   do
   304)     if (.not.associated(cur_sandbox)) exit
   305)     call cur_sandbox%Update(time,option)
   306)     cur_sandbox => cur_sandbox%next
   307)   enddo 
   308)   
   309)   if (print_mass_balance) then
   310)     call SSSandboxOutput(sandbox_list,option,output_option)
   311)   endif
   312) 
   313) end subroutine SSSandboxUpdate
   314) 
   315) ! ************************************************************************** !
   316) 
   317) function SSSandboxOutputFilename(option)
   318)   ! 
   319)   ! Generates filename for source/sink sandbox output
   320)   ! 
   321)   ! Author: Glenn Hammond
   322)   ! Date: 02/23/16
   323) 
   324)   use Option_module
   325) 
   326)   implicit none
   327)   
   328)   type(option_type) :: option
   329)   character(len=MAXSTRINGLENGTH) :: SSSandboxOutputFilename
   330)   character(len=MAXWORDLENGTH) :: word
   331) 
   332)   write(word,'(i6)') option%myrank
   333)   SSSandboxOutputFilename = trim(option%global_prefix) // &
   334)                             trim(option%group_prefix) // &
   335)                             '-ss_mass-' // trim(adjustl(word)) // '.dat'
   336)   
   337) end function SSSandboxOutputFilename  
   338) 
   339) ! ************************************************************************** !
   340) 
   341) subroutine SSSandboxOutputHeader(sandbox_list,grid,option,output_option)
   342)   ! 
   343)   ! Writes header for source/sink sandbox output
   344)   ! 
   345)   ! Author: Glenn Hammond
   346)   ! Date: 02/23/16
   347) 
   348)   use Option_module
   349)   use Output_Aux_module
   350)   use Grid_module
   351)   use Utility_module, only : BestFloat
   352)   
   353)   implicit none
   354)   
   355)   class(srcsink_sandbox_base_type), pointer :: sandbox_list
   356)   type(grid_type) :: grid
   357)   type(option_type) :: option
   358)   type(output_option_type) :: output_option
   359) 
   360)   class(srcsink_sandbox_base_type), pointer :: cur_srcsink
   361)   character(len=MAXSTRINGLENGTH) :: cell_string
   362)   character(len=MAXWORDLENGTH) :: x_string, y_string, z_string
   363)   character(len=MAXWORDLENGTH) :: units_string, variable_string
   364)   character(len=MAXSTRINGLENGTH) :: filename
   365)   PetscInt :: fid
   366)   PetscInt :: icolumn, i
   367)   PetscInt :: local_id, ghosted_id
   368)   
   369)   filename = SSSandboxOutputFilename(option)
   370)   open(unit=IUNIT_TEMP,file=filename,action="write",status="replace")  
   371)   
   372)   if (output_option%print_column_ids) then
   373)     icolumn = 1
   374)   else
   375)     icolumn = -1
   376)   endif 
   377)   
   378)   write(IUNIT_TEMP,'(a)',advance="no") ' "Time [' // &
   379)     trim(output_option%tunit) // ']"'
   380) 
   381)   cur_srcsink => sandbox_list
   382)   do
   383)     if (.not.associated(cur_srcsink)) exit
   384)     local_id = cur_srcsink%local_cell_id
   385)     ghosted_id = grid%nL2G(local_id)
   386) 
   387)     ! cell natural id
   388)     write(cell_string,*) grid%nG2A(ghosted_id)
   389)     cell_string = ' (' // trim(adjustl(cell_string)) // ')'
   390)     ! coordinate of cell
   391)     x_string = BestFloat(grid%x(ghosted_id),1.d4,1.d-2)
   392)     y_string = BestFloat(grid%y(ghosted_id),1.d4,1.d-2)
   393)     z_string = BestFloat(grid%z(ghosted_id),1.d4,1.d-2)
   394)     cell_string = trim(cell_string) // &
   395)              ' (' // trim(adjustl(x_string)) // &
   396)              ' ' // trim(adjustl(y_string)) // &
   397)              ' ' // trim(adjustl(z_string)) // ')'
   398)     select case(option%iflowmode)
   399)       case(RICHARDS_MODE,G_MODE)
   400)         variable_string = ' Water'
   401)         ! cumulative
   402)         units_string = 'kg'
   403)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   404)                                  cell_string,icolumn)
   405)         ! instantaneous
   406)         units_string = 'kg/' // trim(adjustl(output_option%tunit))
   407)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   408)                                  cell_string,icolumn)
   409)     end select
   410)     select case(option%iflowmode)
   411)       case(G_MODE)
   412)         variable_string = ' Gas Component'
   413)         ! cumulative
   414)         units_string = 'kg'
   415)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   416)                                  cell_string,icolumn)
   417)         ! instantaneous
   418)         units_string = 'kg/' // trim(adjustl(output_option%tunit))
   419)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   420)                                  cell_string,icolumn)
   421)         variable_string = ' Energy'
   422)         ! cumulative
   423)         units_string = 'MJ'
   424)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   425)                                  cell_string,icolumn)
   426)         ! instantaneous
   427)         units_string = 'MJ/' // trim(adjustl(output_option%tunit))
   428)         call OutputWriteToHeader(IUNIT_TEMP,variable_string,units_string, &
   429)                                  cell_string,icolumn)
   430)     end select
   431)     cur_srcsink => cur_srcsink%next
   432)   enddo
   433)   
   434)   close(IUNIT_TEMP)
   435)   
   436) end subroutine SSSandboxOutputHeader
   437) 
   438) ! ************************************************************************** !
   439) 
   440) subroutine SSSandboxOutput(sandbox_list,option,output_option)
   441)   ! 
   442)   ! Writes output for for source/sink sandbox
   443)   ! 
   444)   ! Author: Glenn Hammond
   445)   ! Date: 02/23/16
   446) 
   447)   use Option_module
   448)   use Output_Aux_module
   449)   use General_Aux_module, only : fmw_comp
   450) 
   451)   implicit none
   452)   
   453)   class(srcsink_sandbox_base_type), pointer :: sandbox_list
   454)   type(option_type) :: option
   455)   type(output_option_type) :: output_option
   456)   
   457)   class(srcsink_sandbox_base_type), pointer :: cur_srcsink
   458)   character(len=MAXSTRINGLENGTH) :: filename
   459)   PetscInt :: i
   460)   PetscReal :: flow_dof_scale(3)  
   461)   
   462)   if (.not.associated(sandbox_list)) return
   463) 
   464)   flow_dof_scale = 1.d0
   465)   select case(option%iflowmode)
   466)     case(RICHARDS_MODE)
   467)       flow_dof_scale(1) = FMWH2O
   468)     case(TH_MODE)
   469)       flow_dof_scale(1) = FMWH2O
   470)     case(MIS_MODE)
   471)       flow_dof_scale(1) = FMWH2O
   472)       flow_dof_scale(2) = FMWGLYC
   473)     case(G_MODE)
   474)       flow_dof_scale(1) = fmw_comp(1)
   475)       flow_dof_scale(2) = fmw_comp(2)
   476)     case(MPH_MODE,FLASH2_MODE,IMS_MODE)
   477)       flow_dof_scale(1) = FMWH2O
   478)       flow_dof_scale(2) = FMWCO2
   479)   end select  
   480)   
   481) 100 format(100es16.8)
   482) 
   483)   filename = SSSandboxOutputFilename(option)
   484)   open(unit=IUNIT_TEMP,file=filename,action="write",status="old", &
   485)        position="append")
   486) 
   487)   ! this time is set at the end of the reactive transport step
   488)   write(IUNIT_TEMP,100,advance="no") option%time / output_option%tconv
   489)   
   490)   cur_srcsink => sandbox_list
   491)   do
   492)     if (.not.associated(cur_srcsink)) exit
   493)     do i = 1, option%nflowdof
   494)       write(IUNIT_TEMP,100,advance="no") &
   495)         cur_srcsink%cumulative_mass(i)*flow_dof_scale(i), &
   496)         cur_srcsink%instantaneous_mass_rate(i)*flow_dof_scale(i)* &
   497)           output_option%tconv
   498)     enddo
   499)     cur_srcsink => cur_srcsink%next
   500)   enddo
   501)   close(IUNIT_TEMP)
   502)   
   503) end subroutine SSSandboxOutput
   504) 
   505) ! ************************************************************************** !
   506) 
   507) subroutine SSSandboxDestroy(sandbox)
   508)   ! 
   509)   ! Destroys arbitrary sandbox list
   510)   ! 
   511)   ! Author: Glenn Hammond
   512)   ! Date: 04/11/14
   513)   ! 
   514) 
   515)   implicit none
   516) 
   517)   class(srcsink_sandbox_base_type), pointer :: sandbox
   518) 
   519)   if (.not.associated(sandbox)) return
   520) 
   521)   call sandbox%Destroy()
   522)   deallocate(sandbox)
   523)   nullify(sandbox)
   524) 
   525) end subroutine SSSandboxDestroy
   526) 
   527) 
   528) ! ************************************************************************** !
   529) 
   530) subroutine SSSandboxDestroyList1()
   531)   ! 
   532)   ! Destroys master sandbox list
   533)   ! 
   534)   ! Author: Glenn Hammond
   535)   ! Date: 04/11/14
   536)   ! 
   537) 
   538)   implicit none
   539) 
   540)   call SSSandboxDestroyList(ss_sandbox_list)
   541)   
   542) end subroutine SSSandboxDestroyList1
   543) 
   544) ! ************************************************************************** !
   545) 
   546) subroutine SSSandboxDestroyList2(local_sandbox_list)
   547)   ! 
   548)   ! Destroys arbitrary sandbox list
   549)   ! 
   550)   ! Author: Glenn Hammond
   551)   ! Date: 04/11/14
   552)   ! 
   553) 
   554)   implicit none
   555) 
   556)   class(srcsink_sandbox_base_type), pointer :: local_sandbox_list
   557) 
   558)   class(srcsink_sandbox_base_type), pointer :: cur_sandbox, prev_sandbox
   559)   
   560)   ! sandbox source/sinks
   561)   cur_sandbox => local_sandbox_list
   562)   do
   563)     if (.not.associated(cur_sandbox)) exit
   564)     prev_sandbox => cur_sandbox%next
   565)     call SSSandboxDestroy(cur_sandbox)
   566)     cur_sandbox => prev_sandbox
   567)   enddo  
   568) 
   569) end subroutine SSSandboxDestroyList2
   570) 
   571) end module SrcSink_Sandbox_module

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