MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
flowmodelinterfacemodule Module Reference

Data Types

type  flowmodelinterfacetype
 

Functions/Subroutines

subroutine fmi_df (this, dis, idryinactive)
 Define the flow model interface. More...
 
subroutine fmi_ar (this, ibound)
 Allocate the package. More...
 
subroutine fmi_da (this)
 Deallocate variables. More...
 
subroutine allocate_scalars (this)
 Allocate scalars. More...
 
subroutine allocate_arrays (this, nodes)
 Allocate arrays. More...
 
subroutine source_options (this)
 @ brief Source input options for package More...
 
subroutine source_packagedata (this)
 @ brief Source input options for package More...
 
subroutine read_grid (this)
 Validate flow model grid. More...
 
subroutine initialize_bfr (this)
 Initialize the budget file reader. More...
 
subroutine advance_bfr (this)
 Advance the budget file reader. More...
 
subroutine finalize_bfr (this)
 Finalize the budget file reader. More...
 
subroutine initialize_hfr (this)
 Initialize the head file reader. More...
 
subroutine advance_hfr (this)
 Advance the head file reader. More...
 
subroutine finalize_hfr (this)
 Finalize the head file reader. More...
 
subroutine initialize_gwfterms_from_bfr (this)
 Initialize gwf terms from budget file. More...
 
subroutine initialize_gwfterms_from_gwfbndlist (this)
 Initialize gwf terms from a GWF exchange. More...
 
subroutine allocate_gwfpackages (this, ngwfterms)
 Allocate budget packages. More...
 
subroutine deallocate_gwfpackages (this)
 Deallocate memory in the gwfpackages array. More...
 
subroutine get_package_index (this, name, idx)
 Find the package index for the package with the given name. More...
 

Function/Subroutine Documentation

◆ advance_bfr()

subroutine flowmodelinterfacemodule::advance_bfr ( class(flowmodelinterfacetype this)
private

Advance the budget file reader by reading the next chunk of information for the current time step and stress period.

Definition at line 567 of file FlowModelInterface.f90.

568  ! -- modules
569  use tdismodule, only: kstp, kper
570  ! -- dummy
571  class(FlowModelInterfaceType) :: this
572  ! -- local
573  logical :: success
574  integer(I4B) :: n
575  integer(I4B) :: ipos
576  integer(I4B) :: nu, nr
577  integer(I4B) :: ip, i
578  logical :: readnext
579  ! -- format
580  character(len=*), parameter :: fmtkstpkper = &
581  "(1x,/1x,'FMI READING BUDGET TERMS &
582  &FOR KSTP ', i0, ' KPER ', i0)"
583  character(len=*), parameter :: fmtbudkstpkper = &
584  "(1x,/1x, 'FMI SETTING BUDGET TERMS &
585  &FOR KSTP ', i0, ' AND KPER ', &
586  &i0, ' TO BUDGET FILE TERMS FROM &
587  &KSTP ', i0, ' AND KPER ', i0)"
588  !
589  ! -- If the latest record read from the budget file is from a stress
590  ! -- period with only one time step, reuse that record (do not read a
591  ! -- new record) if the running model is still in that same stress period,
592  ! -- or if that record is the last one in the budget file.
593  readnext = .true.
594  if (kstp * kper > 1) then
595  if (this%bfr%header%kstp == 1) then
596  if (this%bfr%endoffile) then
597  readnext = .false.
598  else if (this%bfr%headernext%kper == kper + 1) then
599  readnext = .false.
600  end if
601  else if (this%bfr%endoffile) then
602  write (errmsg, '(4x,a)') 'REACHED END OF GWF BUDGET &
603  &FILE BEFORE READING SUFFICIENT BUDGET INFORMATION FOR THIS &
604  &GWT SIMULATION.'
605  call store_error(errmsg)
606  call store_error_unit(this%iubud)
607  end if
608  end if
609  !
610  ! -- Read the next record
611  if (readnext) then
612  !
613  ! -- Write the current time step and stress period
614  write (this%iout, fmtkstpkper) kstp, kper
615  !
616  ! -- loop through the budget terms for this stress period
617  ! i is the counter for gwf flow packages
618  ip = 1
619  do n = 1, this%bfr%nbudterms
620  call this%bfr%read_record(success, this%iout)
621  if (.not. success) then
622  write (errmsg, '(4x,a)') 'GWF BUDGET READ NOT SUCCESSFUL'
623  call store_error(errmsg)
624  call store_error_unit(this%iubud)
625  end if
626  !
627  ! -- Ensure kper is same between model and budget file
628  if (kper /= this%bfr%header%kper) then
629  write (errmsg, '(4x,a)') 'PERIOD NUMBER IN BUDGET FILE &
630  &DOES NOT MATCH PERIOD NUMBER IN TRANSPORT MODEL. IF THERE &
631  &IS MORE THAN ONE TIME STEP IN THE BUDGET FILE FOR A GIVEN &
632  &STRESS PERIOD, BUDGET FILE TIME STEPS MUST MATCH GWT MODEL &
633  &TIME STEPS ONE-FOR-ONE IN THAT STRESS PERIOD.'
634  call store_error(errmsg)
635  call store_error_unit(this%iubud)
636  end if
637  !
638  ! -- if budget file kstp > 1, then kstp must match
639  if (this%bfr%header%kstp > 1 .and. (kstp /= this%bfr%header%kstp)) then
640  write (errmsg, '(4x,a)') 'TIME STEP NUMBER IN BUDGET FILE &
641  &DOES NOT MATCH TIME STEP NUMBER IN TRANSPORT MODEL. IF THERE &
642  &IS MORE THAN ONE TIME STEP IN THE BUDGET FILE FOR A GIVEN STRESS &
643  &PERIOD, BUDGET FILE TIME STEPS MUST MATCH GWT MODEL TIME STEPS &
644  &ONE-FOR-ONE IN THAT STRESS PERIOD.'
645  call store_error(errmsg)
646  call store_error_unit(this%iubud)
647  end if
648  !
649  ! -- parse based on the type of data, and compress all user node
650  ! numbers into reduced node numbers
651  select type (h => this%bfr%header)
652  type is (budgetfileheadertype)
653  select case (trim(adjustl(h%budtxt)))
654  case ('FLOW-JA-FACE')
655  !
656  ! -- bfr%flowja contains only reduced connections so there is
657  ! a one-to-one match with this%gwfflowja
658  do ipos = 1, size(this%bfr%flowja)
659  this%gwfflowja(ipos) = this%bfr%flowja(ipos)
660  end do
661  case ('DATA-SPDIS')
662  do i = 1, h%nlist
663  nu = this%bfr%nodesrc(i)
664  nr = this%dis%get_nodenumber(nu, 0)
665  if (nr <= 0) cycle
666  this%gwfspdis(1, nr) = this%bfr%auxvar(1, i)
667  this%gwfspdis(2, nr) = this%bfr%auxvar(2, i)
668  this%gwfspdis(3, nr) = this%bfr%auxvar(3, i)
669  end do
670  case ('DATA-SAT')
671  do i = 1, h%nlist
672  nu = this%bfr%nodesrc(i)
673  nr = this%dis%get_nodenumber(nu, 0)
674  if (nr <= 0) cycle
675  this%gwfsat(nr) = this%bfr%auxvar(1, i)
676  end do
677  case ('STO-SS')
678  do nu = 1, this%dis%nodesuser
679  nr = this%dis%get_nodenumber(nu, 0)
680  if (nr <= 0) cycle
681  this%gwfstrgss(nr) = this%bfr%flow(nu)
682  end do
683  case ('STO-SY')
684  do nu = 1, this%dis%nodesuser
685  nr = this%dis%get_nodenumber(nu, 0)
686  if (nr <= 0) cycle
687  this%gwfstrgsy(nr) = this%bfr%flow(nu)
688  end do
689  case default
690  call this%gwfpackages(ip)%copy_values( &
691  h%nlist, &
692  this%bfr%nodesrc, &
693  this%bfr%flow, &
694  this%bfr%auxvar)
695  do i = 1, this%gwfpackages(ip)%nbound
696  nu = this%gwfpackages(ip)%nodelist(i)
697  nr = this%dis%get_nodenumber(nu, 0)
698  this%gwfpackages(ip)%nodelist(i) = nr
699  end do
700  ip = ip + 1
701  end select
702  end select
703  end do
704  else
705  !
706  ! -- write message to indicate that flows are being reused
707  write (this%iout, fmtbudkstpkper) kstp, kper, &
708  this%bfr%header%kstp, this%bfr%header%kper
709  !
710  ! -- set the flag to indicate that flows were not updated
711  this%iflowsupdated = 0
712  end if
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23
Here is the call graph for this function:

◆ advance_hfr()

subroutine flowmodelinterfacemodule::advance_hfr ( class(flowmodelinterfacetype this)
private

Definition at line 730 of file FlowModelInterface.f90.

731  ! modules
732  use tdismodule, only: kstp, kper
733  class(FlowModelInterfaceType) :: this
734  integer(I4B) :: nu, nr, i, ilay
735  integer(I4B) :: ncpl
736  real(DP) :: val
737  logical :: readnext
738  logical :: success
739  character(len=*), parameter :: fmtkstpkper = &
740  "(1x,/1x,'FMI READING HEAD FOR &
741  &KSTP ', i0, ' KPER ', i0)"
742  character(len=*), parameter :: fmthdskstpkper = &
743  "(1x,/1x, 'FMI SETTING HEAD FOR KSTP ', i0, ' AND KPER ', &
744  &i0, ' TO BINARY FILE HEADS FROM KSTP ', i0, ' AND KPER ', i0)"
745  !
746  ! -- If the latest record read from the head file is from a stress
747  ! -- period with only one time step, reuse that record (do not read a
748  ! -- new record) if the running model is still in that same stress period,
749  ! -- or if that record is the last one in the head file.
750  readnext = .true.
751  if (kstp * kper > 1) then
752  if (this%hfr%header%kstp == 1) then
753  if (this%hfr%endoffile) then
754  readnext = .false.
755  else if (this%hfr%headernext%kper == kper + 1) then
756  readnext = .false.
757  end if
758  else if (this%hfr%endoffile) then
759  write (errmsg, '(4x,a)') 'REACHED END OF GWF HEAD &
760  &FILE BEFORE READING SUFFICIENT HEAD INFORMATION FOR THIS &
761  &GWT SIMULATION.'
762  call store_error(errmsg)
763  call store_error_unit(this%iuhds)
764  end if
765  end if
766  !
767  ! -- Read the next record
768  if (readnext) then
769  !
770  ! -- write to list file that heads are being read
771  write (this%iout, fmtkstpkper) kstp, kper
772  !
773  ! -- loop through the layered heads for this time step
774  do ilay = 1, this%hfr%nlay
775  !
776  ! -- read next head chunk
777  call this%hfr%read_record(success, this%iout)
778  if (.not. success) then
779  write (errmsg, '(4x,a)') 'GWF HEAD READ NOT SUCCESSFUL'
780  call store_error(errmsg)
781  call store_error_unit(this%iuhds)
782  end if
783  !
784  ! -- Ensure kper is same between model and head file
785  if (kper /= this%hfr%header%kper) then
786  write (errmsg, '(4x,a)') 'PERIOD NUMBER IN HEAD FILE &
787  &DOES NOT MATCH PERIOD NUMBER IN TRANSPORT MODEL. IF THERE &
788  &IS MORE THAN ONE TIME STEP IN THE HEAD FILE FOR A GIVEN STRESS &
789  &PERIOD, HEAD FILE TIME STEPS MUST MATCH GWT MODEL TIME STEPS &
790  &ONE-FOR-ONE IN THAT STRESS PERIOD.'
791  call store_error(errmsg)
792  call store_error_unit(this%iuhds)
793  end if
794  !
795  ! -- if head file kstp > 1, then kstp must match
796  if (this%hfr%header%kstp > 1 .and. (kstp /= this%hfr%header%kstp)) then
797  write (errmsg, '(4x,a)') 'TIME STEP NUMBER IN HEAD FILE &
798  &DOES NOT MATCH TIME STEP NUMBER IN TRANSPORT MODEL. IF THERE &
799  &IS MORE THAN ONE TIME STEP IN THE HEAD FILE FOR A GIVEN STRESS &
800  &PERIOD, HEAD FILE TIME STEPS MUST MATCH GWT MODEL TIME STEPS &
801  &ONE-FOR-ONE IN THAT STRESS PERIOD.'
802  call store_error(errmsg)
803  call store_error_unit(this%iuhds)
804  end if
805  !
806  ! -- fill the head array for this layer and
807  ! compress into reduced form
808  ncpl = size(this%hfr%head)
809  do i = 1, ncpl
810  nu = (ilay - 1) * ncpl + i
811  nr = this%dis%get_nodenumber(nu, 0)
812  val = this%hfr%head(i)
813  if (nr > 0) this%gwfhead(nr) = val
814  end do
815  end do
816  else
817  write (this%iout, fmthdskstpkper) kstp, kper, &
818  this%hfr%header%kstp, this%hfr%header%kper
819  end if
Here is the call graph for this function:

◆ allocate_arrays()

subroutine flowmodelinterfacemodule::allocate_arrays ( class(flowmodelinterfacetype this,
integer(i4b), intent(in)  nodes 
)

Definition at line 240 of file FlowModelInterface.f90.

242  !modules
243  use constantsmodule, only: dzero
244  ! -- dummy
245  class(FlowModelInterfaceType) :: this
246  integer(I4B), intent(in) :: nodes
247  ! -- local
248  integer(I4B) :: n
249  !
250  ! -- Allocate ibdgwfsat0, which is an indicator array marking cells with
251  ! saturation greater than 0.0 with a value of 1
252  call mem_allocate(this%ibdgwfsat0, nodes, 'IBDGWFSAT0', this%memoryPath)
253  do n = 1, nodes
254  this%ibdgwfsat0(n) = 1
255  end do
256  !
257  ! -- Allocate differently depending on whether or not flows are
258  ! being read from a file.
259  if (this%flows_from_file) then
260  call mem_allocate(this%gwfflowja, this%dis%con%nja, &
261  'GWFFLOWJA', this%memoryPath)
262  call mem_allocate(this%gwfsat, nodes, 'GWFSAT', this%memoryPath)
263  call mem_allocate(this%gwfhead, nodes, 'GWFHEAD', this%memoryPath)
264  call mem_allocate(this%gwfspdis, 3, nodes, 'GWFSPDIS', this%memoryPath)
265  do n = 1, nodes
266  this%gwfsat(n) = done
267  this%gwfhead(n) = dzero
268  this%gwfspdis(:, n) = dzero
269  end do
270  do n = 1, size(this%gwfflowja)
271  this%gwfflowja(n) = dzero
272  end do
273  !
274  ! -- allocate and initialize storage arrays
275  if (this%igwfstrgss == 0) then
276  call mem_allocate(this%gwfstrgss, 1, 'GWFSTRGSS', this%memoryPath)
277  else
278  call mem_allocate(this%gwfstrgss, nodes, 'GWFSTRGSS', this%memoryPath)
279  end if
280  if (this%igwfstrgsy == 0) then
281  call mem_allocate(this%gwfstrgsy, 1, 'GWFSTRGSY', this%memoryPath)
282  else
283  call mem_allocate(this%gwfstrgsy, nodes, 'GWFSTRGSY', this%memoryPath)
284  end if
285  do n = 1, size(this%gwfstrgss)
286  this%gwfstrgss(n) = dzero
287  end do
288  do n = 1, size(this%gwfstrgsy)
289  this%gwfstrgsy(n) = dzero
290  end do
291  !
292  ! -- If there is no fmi package, then there are no flows at all or a
293  ! connected GWF model, so allocate gwfpackages to zero
294  if (this%inunit == 0) call this%allocate_gwfpackages(this%nflowpack)
295  end if
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65

◆ allocate_gwfpackages()

subroutine flowmodelinterfacemodule::allocate_gwfpackages ( class(flowmodelinterfacetype this,
integer(i4b), intent(in)  ngwfterms 
)

gwfpackages is an array of PackageBudget objects. This routine allocates gwfpackages to the proper size and initializes some member variables.

Definition at line 996 of file FlowModelInterface.f90.

997  ! -- modules
998  use constantsmodule, only: lenmempath
1000  ! -- dummy
1001  class(FlowModelInterfaceType) :: this
1002  integer(I4B), intent(in) :: ngwfterms
1003  ! -- local
1004  integer(I4B) :: n
1005  character(len=LENMEMPATH) :: memPath
1006  !
1007  ! -- direct allocate
1008  allocate (this%gwfpackages(ngwfterms))
1009  allocate (this%flowpacknamearray(ngwfterms))
1010  !
1011  ! -- mem_allocate
1012  call mem_allocate(this%igwfmvrterm, ngwfterms, 'IGWFMVRTERM', this%memoryPath)
1013  !
1014  ! -- initialize
1015  this%nflowpack = ngwfterms
1016  do n = 1, this%nflowpack
1017  this%igwfmvrterm(n) = 0
1018  this%flowpacknamearray(n) = ''
1019  !
1020  ! -- Create a mempath for each individual flow package data set
1021  ! of the form, MODELNAME/FMI-FTn
1022  write (mempath, '(a, i0)') trim(this%memoryPath)//'-FT', n
1023  call this%gwfpackages(n)%initialize(mempath)
1024  end do
integer(i4b), parameter lenmempath
maximum length of the memory path
Definition: Constants.f90:27

◆ allocate_scalars()

subroutine flowmodelinterfacemodule::allocate_scalars ( class(flowmodelinterfacetype this)

Definition at line 199 of file FlowModelInterface.f90.

200  ! -- modules
203  ! -- dummy
204  class(FlowModelInterfaceType) :: this
205  ! -- local
206  !
207  ! -- allocate scalars in NumericalPackageType
208  call this%NumericalPackageType%allocate_scalars()
209  !
210  ! -- Allocate
211  call mem_allocate(this%flows_from_file, 'FLOWS_FROM_FILE', this%memoryPath)
212  call mem_allocate(this%iflowsupdated, 'IFLOWSUPDATED', this%memoryPath)
213  call mem_allocate(this%igwfspdis, 'IGWFSPDIS', this%memoryPath)
214  call mem_allocate(this%igwfstrgss, 'IGWFSTRGSS', this%memoryPath)
215  call mem_allocate(this%igwfstrgsy, 'IGWFSTRGSY', this%memoryPath)
216  call mem_allocate(this%iubud, 'IUBUD', this%memoryPath)
217  call mem_allocate(this%iuhds, 'IUHDS', this%memoryPath)
218  call mem_allocate(this%iumvr, 'IUMVR', this%memoryPath)
219  call mem_allocate(this%iugrb, 'IUGRB', this%memoryPath)
220  call mem_allocate(this%nflowpack, 'NFLOWPACK', this%memoryPath)
221  call mem_allocate(this%idryinactive, "IDRYINACTIVE", this%memoryPath)
222  !
223  ! !
224  ! -- Initialize
225  this%flows_from_file = .true.
226  this%iflowsupdated = 1
227  this%igwfspdis = 0
228  this%igwfstrgss = 0
229  this%igwfstrgsy = 0
230  this%iubud = 0
231  this%iuhds = 0
232  this%iumvr = 0
233  this%iugrb = 0
234  this%nflowpack = 0
235  this%idryinactive = 1

◆ deallocate_gwfpackages()

subroutine flowmodelinterfacemodule::deallocate_gwfpackages ( class(flowmodelinterfacetype this)

Definition at line 1028 of file FlowModelInterface.f90.

1029  class(FlowModelInterfaceType) :: this
1030  integer(I4B) :: n
1031 
1032  do n = 1, this%nflowpack
1033  call this%gwfpackages(n)%da()
1034  end do

◆ finalize_bfr()

subroutine flowmodelinterfacemodule::finalize_bfr ( class(flowmodelinterfacetype this)

Definition at line 716 of file FlowModelInterface.f90.

717  class(FlowModelInterfaceType) :: this
718  call this%bfr%finalize()

◆ finalize_hfr()

subroutine flowmodelinterfacemodule::finalize_hfr ( class(flowmodelinterfacetype this)

Definition at line 823 of file FlowModelInterface.f90.

824  class(FlowModelInterfaceType) :: this
825  close (this%iuhds)

◆ fmi_ar()

subroutine flowmodelinterfacemodule::fmi_ar ( class(flowmodelinterfacetype this,
integer(i4b), dimension(:), pointer, contiguous  ibound 
)
private

Definition at line 138 of file FlowModelInterface.f90.

139  ! -- modules
140  ! -- dummy
141  class(FlowModelInterfaceType) :: this
142  integer(I4B), dimension(:), pointer, contiguous :: ibound
143  !
144  ! -- store pointers to arguments that were passed in
145  this%ibound => ibound
146  !
147  ! -- Allocate arrays
148  call this%allocate_arrays(this%dis%nodes)

◆ fmi_da()

subroutine flowmodelinterfacemodule::fmi_da ( class(flowmodelinterfacetype this)
private

Definition at line 153 of file FlowModelInterface.f90.

154  ! -- modules
156  ! -- dummy
157  class(FlowModelInterfaceType) :: this
158  ! -- todo: finalize hfr and bfr either here or in a finalize routine
159  !
160  ! -- deallocate any memory stored with gwfpackages
161  call this%deallocate_gwfpackages()
162  !
163  ! -- deallocate fmi arrays
164  deallocate (this%gwfpackages)
165  deallocate (this%flowpacknamearray)
166  call mem_deallocate(this%igwfmvrterm)
167  call mem_deallocate(this%ibdgwfsat0)
168  !
169  if (this%flows_from_file) then
170  call mem_deallocate(this%gwfstrgss)
171  call mem_deallocate(this%gwfstrgsy)
172  end if
173  !
174  ! -- special treatment, these could be from mem_checkin
175  call mem_deallocate(this%gwfhead, 'GWFHEAD', this%memoryPath)
176  call mem_deallocate(this%gwfsat, 'GWFSAT', this%memoryPath)
177  call mem_deallocate(this%gwfspdis, 'GWFSPDIS', this%memoryPath)
178  call mem_deallocate(this%gwfflowja, 'GWFFLOWJA', this%memoryPath)
179  !
180  ! -- deallocate scalars
181  call mem_deallocate(this%flows_from_file)
182  call mem_deallocate(this%iflowsupdated)
183  call mem_deallocate(this%igwfspdis)
184  call mem_deallocate(this%igwfstrgss)
185  call mem_deallocate(this%igwfstrgsy)
186  call mem_deallocate(this%iubud)
187  call mem_deallocate(this%iuhds)
188  call mem_deallocate(this%iumvr)
189  call mem_deallocate(this%iugrb)
190  call mem_deallocate(this%nflowpack)
191  call mem_deallocate(this%idryinactive)
192  !
193  ! -- deallocate parent
194  call this%NumericalPackageType%da()

◆ fmi_df()

subroutine flowmodelinterfacemodule::fmi_df ( class(flowmodelinterfacetype this,
class(disbasetype), intent(in), pointer  dis,
integer(i4b), intent(in)  idryinactive 
)
private

Definition at line 82 of file FlowModelInterface.f90.

83  ! -- modules
84  ! -- dummy
85  class(FlowModelInterfaceType) :: this
86  class(DisBaseType), pointer, intent(in) :: dis
87  integer(I4B), intent(in) :: idryinactive
88  ! -- formats
89  character(len=*), parameter :: fmtfmi = &
90  "(1x,/1x,'FMI -- FLOW MODEL INTERFACE, VERSION 2, 8/17/2023', &
91  &' INPUT READ FROM MEMPATH: ', A, //)"
92  character(len=*), parameter :: fmtfmi0 = &
93  "(1x,/1x,'FMI -- FLOW MODEL INTERFACE,'&
94  &' VERSION 2, 8/17/2023')"
95  !
96  ! --print a message identifying the FMI package.
97  if (this%iout > 0) then
98  if (this%inunit /= 0) then
99  write (this%iout, fmtfmi) this%input_mempath
100  else
101  write (this%iout, fmtfmi0)
102  if (this%flows_from_file) then
103  write (this%iout, '(a)') ' FLOWS ARE ASSUMED TO BE ZERO.'
104  else
105  write (this%iout, '(a)') ' FLOWS PROVIDED BY A GWF MODEL IN THIS &
106  &SIMULATION'
107  end if
108  end if
109  end if
110  !
111  ! -- Store pointers
112  this%dis => dis
113  !
114  ! -- Read fmi options
115  if (this%inunit /= 0) then
116  call this%source_options()
117  end if
118  !
119  ! -- Read packagedata options
120  if (this%inunit /= 0 .and. this%flows_from_file) then
121  call this%source_packagedata()
122  call this%initialize_gwfterms_from_bfr()
123  end if
124  !
125  ! -- If GWF-Model exchange is active, setup flow terms
126  if (.not. this%flows_from_file) then
127  call this%initialize_gwfterms_from_gwfbndlist()
128  end if
129  !
130  ! -- Set flag that stops dry flows from being deactivated in a GWE
131  ! transport model since conduction will still be simulated.
132  ! 0: GWE (skip deactivation step); 1: GWT (default: use existing code)
133  this%idryinactive = idryinactive

◆ get_package_index()

subroutine flowmodelinterfacemodule::get_package_index ( class(flowmodelinterfacetype this,
character(len=*), intent(in)  name,
integer(i4b), intent(inout)  idx 
)
private

Definition at line 1038 of file FlowModelInterface.f90.

1039  use bndmodule, only: bndtype, getbndfromlist
1040  class(FlowModelInterfaceType) :: this
1041  character(len=*), intent(in) :: name
1042  integer(I4B), intent(inout) :: idx
1043  ! -- local
1044  integer(I4B) :: ip
1045  !
1046  ! -- Look through all the packages and return the index with name
1047  idx = 0
1048  do ip = 1, size(this%flowpacknamearray)
1049  if (this%flowpacknamearray(ip) == name) then
1050  idx = ip
1051  exit
1052  end if
1053  end do
1054  if (idx == 0) then
1055  call store_error('Error in get_package_index. Could not find '//name, &
1056  terminate=.true.)
1057  end if
This module contains the base boundary package.
class(bndtype) function, pointer, public getbndfromlist(list, idx)
Get boundary from package list.
@ brief BndType
Here is the call graph for this function:

◆ initialize_bfr()

subroutine flowmodelinterfacemodule::initialize_bfr ( class(flowmodelinterfacetype this)

Definition at line 554 of file FlowModelInterface.f90.

555  class(FlowModelInterfaceType) :: this
556  integer(I4B) :: ncrbud
557  call this%bfr%initialize(this%iubud, this%iout, ncrbud)
558  ! todo: need to run through the budget terms
559  ! and do some checking

◆ initialize_gwfterms_from_bfr()

subroutine flowmodelinterfacemodule::initialize_gwfterms_from_bfr ( class(flowmodelinterfacetype this)
private

initialize terms and figure out how many different terms and packages are contained within the file

Definition at line 833 of file FlowModelInterface.f90.

834  ! -- modules
836  ! -- dummy
837  class(FlowModelInterfaceType) :: this
838  ! -- local
839  integer(I4B) :: nflowpack
840  integer(I4B) :: i, ip
841  integer(I4B) :: naux
842  logical :: found_flowja
843  logical :: found_dataspdis
844  logical :: found_datasat
845  logical :: found_stoss
846  logical :: found_stosy
847  integer(I4B), dimension(:), allocatable :: imap
848  !
849  ! -- Calculate the number of gwf flow packages
850  allocate (imap(this%bfr%nbudterms))
851  imap(:) = 0
852  nflowpack = 0
853  found_flowja = .false.
854  found_dataspdis = .false.
855  found_datasat = .false.
856  found_stoss = .false.
857  found_stosy = .false.
858  do i = 1, this%bfr%nbudterms
859  select case (trim(adjustl(this%bfr%budtxtarray(i))))
860  case ('FLOW-JA-FACE')
861  found_flowja = .true.
862  case ('DATA-SPDIS')
863  found_dataspdis = .true.
864  this%igwfspdis = 1
865  case ('DATA-SAT')
866  found_datasat = .true.
867  case ('STO-SS')
868  found_stoss = .true.
869  this%igwfstrgss = 1
870  case ('STO-SY')
871  found_stosy = .true.
872  this%igwfstrgsy = 1
873  case default
874  nflowpack = nflowpack + 1
875  imap(i) = 1
876  end select
877  end do
878  !
879  ! -- allocate gwfpackage arrays
880  call this%allocate_gwfpackages(nflowpack)
881  !
882  ! -- Copy the package name and aux names from budget file reader
883  ! to the gwfpackages derived-type variable
884  ip = 1
885  do i = 1, this%bfr%nbudterms
886  if (imap(i) == 0) cycle
887  call this%gwfpackages(ip)%set_name(this%bfr%dstpackagenamearray(i), &
888  this%bfr%budtxtarray(i))
889  naux = this%bfr%nauxarray(i)
890  call this%gwfpackages(ip)%set_auxname(naux, this%bfr%auxtxtarray(1:naux, i))
891  ip = ip + 1
892  end do
893  !
894  ! -- Copy just the package names for the boundary packages into
895  ! the flowpacknamearray
896  ip = 1
897  do i = 1, size(imap)
898  if (imap(i) == 1) then
899  this%flowpacknamearray(ip) = this%bfr%dstpackagenamearray(i)
900  ip = ip + 1
901  end if
902  end do
903  !
904  ! -- Error if specific discharge, saturation or flowja not found
905  if (.not. found_dataspdis) then
906  write (errmsg, '(4x,a)') 'SPECIFIC DISCHARGE NOT FOUND IN &
907  &BUDGET FILE. SAVE_SPECIFIC_DISCHARGE AND &
908  &SAVE_FLOWS MUST BE ACTIVATED IN THE NPF PACKAGE.'
909  call store_error(errmsg)
910  end if
911  if (.not. found_datasat) then
912  write (errmsg, '(4x,a)') 'SATURATION NOT FOUND IN &
913  &BUDGET FILE. SAVE_SATURATION AND &
914  &SAVE_FLOWS MUST BE ACTIVATED IN THE NPF PACKAGE.'
915  call store_error(errmsg)
916  end if
917  if (.not. found_flowja) then
918  write (errmsg, '(4x,a)') 'FLOWJA NOT FOUND IN &
919  &BUDGET FILE. SAVE_FLOWS MUST &
920  &BE ACTIVATED IN THE NPF PACKAGE.'
921  call store_error(errmsg)
922  end if
923  if (count_errors() > 0) then
924  call store_error_filename(this%input_fname)
925  end if
Here is the call graph for this function:

◆ initialize_gwfterms_from_gwfbndlist()

subroutine flowmodelinterfacemodule::initialize_gwfterms_from_gwfbndlist ( class(flowmodelinterfacetype this)

Definition at line 929 of file FlowModelInterface.f90.

930  ! -- modules
931  use bndmodule, only: bndtype, getbndfromlist
932  ! -- dummy
933  class(FlowModelInterfaceType) :: this
934  ! -- local
935  integer(I4B) :: ngwfpack
936  integer(I4B) :: ngwfterms
937  integer(I4B) :: ip
938  integer(I4B) :: imover
939  integer(I4B) :: ntomvr
940  integer(I4B) :: iterm
941  character(len=LENPACKAGENAME) :: budtxt
942  class(BndType), pointer :: packobj => null()
943  !
944  ! -- determine size of gwf terms
945  ngwfpack = this%gwfbndlist%Count()
946  !
947  ! -- Count number of to-mvr terms, but do not include advanced packages
948  ! as those mover terms are not losses from the cell, but rather flows
949  ! within the advanced package
950  ntomvr = 0
951  do ip = 1, ngwfpack
952  packobj => getbndfromlist(this%gwfbndlist, ip)
953  imover = packobj%imover
954  if (packobj%isadvpak /= 0) imover = 0
955  if (imover /= 0) then
956  ntomvr = ntomvr + 1
957  end if
958  end do
959  !
960  ! -- Allocate arrays in fmi of size ngwfterms, which is the number of
961  ! packages plus the number of packages with mover terms.
962  ngwfterms = ngwfpack + ntomvr
963  call this%allocate_gwfpackages(ngwfterms)
964  !
965  ! -- Assign values in the fmi package
966  iterm = 1
967  do ip = 1, ngwfpack
968  !
969  ! -- set and store names
970  packobj => getbndfromlist(this%gwfbndlist, ip)
971  budtxt = adjustl(packobj%text)
972  call this%gwfpackages(iterm)%set_name(packobj%packName, budtxt)
973  this%flowpacknamearray(iterm) = packobj%packName
974  iterm = iterm + 1
975  !
976  ! -- if this package has a mover associated with it, then add another
977  ! term that corresponds to the mover flows
978  imover = packobj%imover
979  if (packobj%isadvpak /= 0) imover = 0
980  if (imover /= 0) then
981  budtxt = trim(adjustl(packobj%text))//'-TO-MVR'
982  call this%gwfpackages(iterm)%set_name(packobj%packName, budtxt)
983  this%flowpacknamearray(iterm) = packobj%packName
984  this%igwfmvrterm(iterm) = 1
985  iterm = iterm + 1
986  end if
987  end do
Here is the call graph for this function:

◆ initialize_hfr()

subroutine flowmodelinterfacemodule::initialize_hfr ( class(flowmodelinterfacetype this)
private

Definition at line 722 of file FlowModelInterface.f90.

723  class(FlowModelInterfaceType) :: this
724  call this%hfr%initialize(this%iuhds, this%iout)
725  ! todo: need to run through the head terms
726  ! and do some checking

◆ read_grid()

subroutine flowmodelinterfacemodule::read_grid ( class(flowmodelinterfacetype this)

Definition at line 414 of file FlowModelInterface.f90.

415  ! -- modules
416  use dismodule, only: distype
417  use disvmodule, only: disvtype
418  use disumodule, only: disutype
419  use dis2dmodule, only: dis2dtype
420  use disv2dmodule, only: disv2dtype
421  use disv1dmodule, only: disv1dtype
422  ! -- dummy
423  class(FlowModelInterfaceType) :: this
424  ! -- local
425  integer(I4B) :: user_nodes
426  integer(I4B), allocatable :: idomain1d(:), idomain2d(:, :), idomain3d(:, :, :)
427  ! -- formats
428  character(len=*), parameter :: fmtdiserr = &
429  "('Error in ',a,': Models do not have the same discretization. &
430  &GWF model has ', i0, ' user nodes, this model has ', i0, '. &
431  &Ensure discretization packages, including IDOMAIN, are identical.')"
432  character(len=*), parameter :: fmtidomerr = &
433  "('Error in ',a,': models do not have the same discretization. &
434  &Models have different IDOMAIN arrays. &
435  &Ensure discretization packages, including IDOMAIN, are identical.')"
436 
437  call this%gfr%initialize(this%iugrb)
438 
439  ! check grid equivalence
440  select case (this%gfr%grid_type)
441  case ('DIS')
442  select type (dis => this%dis)
443  type is (distype)
444  user_nodes = this%gfr%read_int("NCELLS")
445  if (user_nodes /= this%dis%nodesuser) then
446  write (errmsg, fmtdiserr) &
447  trim(this%text), user_nodes, this%dis%nodesuser
448  call store_error(errmsg, terminate=.true.)
449  end if
450  idomain1d = this%gfr%read_int_1d("IDOMAIN")
451  idomain3d = reshape(idomain1d, [ &
452  this%gfr%read_int("NCOL"), &
453  this%gfr%read_int("NROW"), &
454  this%gfr%read_int("NLAY") &
455  ])
456  if (.not. all(dis%idomain == idomain3d)) then
457  write (errmsg, fmtidomerr) trim(this%text)
458  call store_error(errmsg, terminate=.true.)
459  end if
460  end select
461  case ('DISV')
462  select type (dis => this%dis)
463  type is (disvtype)
464  user_nodes = this%gfr%read_int("NCELLS")
465  if (user_nodes /= this%dis%nodesuser) then
466  write (errmsg, fmtdiserr) &
467  trim(this%text), user_nodes, this%dis%nodesuser
468  call store_error(errmsg, terminate=.true.)
469  end if
470  idomain1d = this%gfr%read_int_1d("IDOMAIN")
471  idomain2d = reshape(idomain1d, [ &
472  this%gfr%read_int("NCPL"), &
473  this%gfr%read_int("NLAY") &
474  ])
475  if (.not. all(dis%idomain == idomain2d)) then
476  write (errmsg, fmtidomerr) trim(this%text)
477  call store_error(errmsg, terminate=.true.)
478  end if
479  end select
480  case ('DISU')
481  select type (dis => this%dis)
482  type is (disutype)
483  user_nodes = this%gfr%read_int("NODES")
484  if (user_nodes /= this%dis%nodesuser) then
485  write (errmsg, fmtdiserr) &
486  trim(this%text), user_nodes, this%dis%nodesuser
487  call store_error(errmsg, terminate=.true.)
488  end if
489  idomain1d = this%gfr%read_int_1d("IDOMAIN")
490  if (.not. all(dis%idomain == idomain1d)) then
491  write (errmsg, fmtidomerr) trim(this%text)
492  call store_error(errmsg, terminate=.true.)
493  end if
494  end select
495  case ('DIS2D')
496  select type (dis => this%dis)
497  type is (dis2dtype)
498  user_nodes = this%gfr%read_int("NCELLS")
499  if (user_nodes /= this%dis%nodesuser) then
500  write (errmsg, fmtdiserr) &
501  trim(this%text), user_nodes, this%dis%nodesuser
502  call store_error(errmsg, terminate=.true.)
503  end if
504  idomain1d = this%gfr%read_int_1d("IDOMAIN")
505  idomain2d = reshape(idomain1d, [ &
506  this%gfr%read_int("NCOL"), &
507  this%gfr%read_int("NROW") &
508  ])
509  if (.not. all(dis%idomain == idomain2d)) then
510  write (errmsg, fmtidomerr) trim(this%text)
511  call store_error(errmsg, terminate=.true.)
512  end if
513  end select
514  case ('DISV2D')
515  select type (dis => this%dis)
516  type is (disv2dtype)
517  user_nodes = this%gfr%read_int("NODES")
518  if (user_nodes /= this%dis%nodesuser) then
519  write (errmsg, fmtdiserr) &
520  trim(this%text), user_nodes, this%dis%nodesuser
521  call store_error(errmsg, terminate=.true.)
522  end if
523  idomain1d = this%gfr%read_int_1d("IDOMAIN")
524  if (.not. all(dis%idomain == idomain1d)) then
525  write (errmsg, fmtidomerr) trim(this%text)
526  call store_error(errmsg, terminate=.true.)
527  end if
528  end select
529  case ('DISV1D')
530  select type (dis => this%dis)
531  type is (disv1dtype)
532  user_nodes = this%gfr%read_int("NCELLS")
533  if (user_nodes /= this%dis%nodesuser) then
534  write (errmsg, fmtdiserr) &
535  trim(this%text), user_nodes, this%dis%nodesuser
536  call store_error(errmsg, terminate=.true.)
537  end if
538  idomain1d = this%gfr%read_int_1d("IDOMAIN")
539  if (.not. all(dis%idomain == idomain1d)) then
540  write (errmsg, fmtidomerr) trim(this%text)
541  call store_error(errmsg, terminate=.true.)
542  end if
543  end select
544  end select
545 
546  if (allocated(idomain3d)) deallocate (idomain3d)
547  if (allocated(idomain2d)) deallocate (idomain2d)
548  if (allocated(idomain1d)) deallocate (idomain1d)
549 
550  call this%gfr%finalize()
Definition: Dis.f90:1
Structured grid discretization.
Definition: Dis2d.f90:23
Structured grid discretization.
Definition: Dis.f90:23
Unstructured grid discretization.
Definition: Disu.f90:28
Vertex grid discretization.
Definition: Disv2d.f90:24
Vertex grid discretization.
Definition: Disv.f90:24
Here is the call graph for this function:

◆ source_options()

subroutine flowmodelinterfacemodule::source_options ( class(flowmodelinterfacetype this)

Definition at line 300 of file FlowModelInterface.f90.

301  ! -- modules
303  ! -- dummy
304  class(FlowModelInterfaceType) :: this
305  ! -- local
306  logical(LGP) :: found_ipakcb
307  character(len=*), parameter :: fmtisvflow = &
308  "(4x,'CELL-BY-CELL FLOW INFORMATION WILL BE SAVED TO BINARY FILE &
309  &WHENEVER ICBCFL IS NOT ZERO AND FLOW IMBALANCE CORRECTION ACTIVE.')"
310 
311  ! -- source package input
312  call mem_set_value(this%ipakcb, 'SAVE_FLOWS', this%input_mempath, &
313  found_ipakcb)
314 
315  write (this%iout, '(1x,a)') 'PROCESSING FMI OPTIONS'
316 
317  if (found_ipakcb) then
318  this%ipakcb = -1
319  write (this%iout, fmtisvflow)
320  end if
321 
322  write (this%iout, '(1x,a)') 'END OF FMI OPTIONS'

◆ source_packagedata()

subroutine flowmodelinterfacemodule::source_packagedata ( class(flowmodelinterfacetype this)

Definition at line 327 of file FlowModelInterface.f90.

328  ! -- modules
332  use openspecmodule, only: access, form
335  ! -- dummy
336  class(FlowModelInterfaceType) :: this
337  ! -- local
338  type(CharacterStringType), dimension(:), contiguous, &
339  pointer :: flowtypes
340  type(CharacterStringType), dimension(:), contiguous, &
341  pointer :: fileops
342  type(CharacterStringType), dimension(:), contiguous, &
343  pointer :: fnames
344  character(len=LINELENGTH) :: flowtype, fileop, fname
345  integer(I4B) :: inunit, n
346  logical(LGP) :: exist
347 
348  call mem_setptr(flowtypes, 'FLOWTYPE', this%input_mempath)
349  call mem_setptr(fileops, 'FILEIN', this%input_mempath)
350  call mem_setptr(fnames, 'FNAME', this%input_mempath)
351 
352  write (this%iout, '(1x,a)') 'PROCESSING FMI PACKAGEDATA'
353 
354  do n = 1, size(flowtypes)
355  flowtype = flowtypes(n)
356  fileop = fileops(n)
357  fname = fnames(n)
358 
359  inquire (file=trim(fname), exist=exist)
360  if (.not. exist) then
361  call store_error('Could not find file '//trim(fname))
362  cycle
363  end if
364 
365  if (fileop /= 'FILEIN') then
366  call store_error('Unexpected packagedata input keyword read: "' &
367  //trim(fileop)//'".')
368  cycle
369  end if
370 
371  select case (flowtype)
372  case ('GWFBUDGET')
373  inunit = getunit()
374  call openfile(inunit, this%iout, fname, 'DATA(BINARY)', form, &
375  access, 'UNKNOWN')
376  this%iubud = inunit
377  call this%initialize_bfr()
378  case ('GWFHEAD')
379  inunit = getunit()
380  call openfile(inunit, this%iout, fname, 'DATA(BINARY)', form, &
381  access, 'UNKNOWN')
382  this%iuhds = inunit
383  call this%initialize_hfr()
384  case ('GWFMOVER')
385  inunit = getunit()
386  call openfile(inunit, this%iout, fname, 'DATA(BINARY)', form, &
387  access, 'UNKNOWN')
388  this%iumvr = inunit
389  call budgetobject_cr_bfr(this%mvrbudobj, 'MVT', this%iumvr, &
390  this%iout)
391  call this%mvrbudobj%fill_from_bfr(this%dis, this%iout)
392  case ('GWFGRID')
393  inunit = getunit()
394  call openfile(inunit, this%iout, fname, 'DATA(BINARY)', &
395  form, access, 'UNKNOWN')
396  this%iugrb = inunit
397  call this%read_grid()
398  case default
399  write (errmsg, '(a,3(1x,a))') &
400  'UNKNOWN', trim(adjustl(this%text)), 'PACKAGEDATA:', trim(flowtype)
401  call store_error(errmsg)
402  end select
403  end do
404 
405  write (this%iout, '(1x,a)') 'END OF FMI PACKAGEDATA'
406 
407  if (count_errors() > 0) then
408  call store_error_filename(this%input_fname)
409  end if
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
integer(i4b), parameter lenpackagename
maximum length of the package name
Definition: Constants.f90:23
real(dp), parameter dem6
real constant 1e-6
Definition: Constants.f90:109
subroutine, public urdaux(naux, inunit, iout, lloc, istart, istop, auxname, line, text)
Read auxiliary variables from an input line.
integer(i4b) function, public getunit()
Get a free unit number.
subroutine, public openfile(iu, iout, fname, ftype, fmtarg_opt, accarg_opt, filstat_opt, mode_opt)
Open a file.
Definition: InputOutput.f90:30
character(len=20) access
Definition: OpenSpec.f90:7
character(len=20) form
Definition: OpenSpec.f90:7
This class is used to store a single deferred-length character string. It was designed to work in an ...
Definition: CharString.f90:23
Here is the call graph for this function: