Skip to content

Commit

Permalink
merge develop into transport (#322)
Browse files Browse the repository at this point in the history 
* fix(uzf): fix indexing error in UZF (#274)

* Error introduced as part of recent UZF refactoring
* Closes #273

* fix(memory): some variables not deallocated (#278)

Implemented new check in develop mode so code bombs with error if memory manager variable not deallocated

* refactor(budobj): new budget object for advanced packages (#279)

* single code base for writing binary budget files for advanced packages
* single code base for creating and writing budget tables to list file for advanced packages
* implemented for MAW, UZF, LAK, SFR, and MVR
* closes #277
* update mf6exes from 2.0 to 3.0
* will allow generalized transport calculations for advanced packages

* refactor(advanced packages): read static data as part of df() (#281)

* refactor(advanced packages): modify advanced packages to read all static data as part of df()
* modify setup_budobj to include the connectivity so that it is available to other models

* fix(budterm): initialize nlist to zero (#283)

* initialize nlist to zero (#284)

* fix(lak): initialize chdratin and chdratout to zero (#286)

* fix(lak): added sign checks for user specified lak flow terms (#288)

* User-specified values for RAINFALL, EVAPORATION, RUNOFF, INFLOW, and WITHDRAWAL must be positive.  The program worked if these values were negative, but that doesn't necessarily make sense and is probably an input error.
* Updated definition file to reflect these changes
* Corrected minor typo in lake definition file
* Updated release notes to reflect this change
* Close #287

* feat(sfr): add storage term to sfr budget (#293)

* Also includes a reach volume term written as an aux variable.  This is needed for transport.
* initialize str so non-ascii characters don't show up in output files
* updated notes for these changes

* * fix(csub): Fix CSUB binary budget data saved as IMETH=6 datatype

closes #290

* refactor(sfr): Refactor SFR Package to remove use of Geometry objects (#296)

* refactor(sfr): Refactor SFR Package to remove use of Geometry objects

* doc(pak-ts): update description of package timeseries variables. 

* ci(yml): update yml to clone shallow copy of flopy and pymake repos

* doc(release): Update release notes

* Closes #276, and #289

* refactor(xt3d): accumulate flowja instead of set (#306)

This change is required for transport, which accumulates terms in flowja

* fix(dis): corrected connection vector error in DIS package (#308)

Code was incorrectly calculating cell center elevations when nozee was .false.  This error would have affected XT3D simulations with unconfined flow and non-zero values specified for ANGLE2 in the NPF Package.  It also would affect dispersive transport with XT3D for unconfined conditions.

* Update sln-ims-example.dat

* refactor(maw): refactor MAW conductance calculation for issue 305 (#310)

Add traps to catch 1) skin factors that are <= 0 when using the SKIN
conductance equation and 2) and negative saturated conductances values.

Closes #305

* feat(tableobj): Add a generic table object for lst file output (#303)

Full implementation for SFR package. Partial implementation for LAK, MAW, and UZF packages.

* feat(tableobj): update MAW package to use tableobj for data output (#315)

Also fix some budget reporting issues in the MAW package

* fix(auxmult): auxmult fix when auxmult and bound are in time series (#316)

* change order of time series interpolation so aux is done first in case it is an auxmult column
* partially addresses #314
* updated release notes

* refactor(BoundaryPackage): Add use of TableObject for print_flows option (#317)

* feat(lnf): update n-point geometry package data

* Revert "feat(lnf): update n-point geometry package data"

This reverts commit 3460cfb.

* fix(lak): revise the way outlet to-mvr flows are stored in budobj

* fix(lak): revise the way outlet to-mvr flows are stored in budobj (#321)

* fix(lak): revise the way outlet to-mvr flows are stored in budobj

* updated release notes

* added flow table to cnc and ssm

Co-authored-by: Hughes, J.D. <jdhughes@usgs.gov>
  • Loading branch information
@langevin-usgs @jdhughes-usgs
langevin-usgs and jdhughes-usgs authored Feb 27, 2020
1 parent 2336f32 commit f55cf79
Show file tree
Hide file tree
Showing 20 changed files with 587 additions and 149 deletions.
1 change: 1 addition & 0 deletions autotest/test_gwf_mvr01.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -179,6 +179,7 @@ def get_model(idx, dir):
lak = flopy.mf6.ModflowGwflak(gwf, mover=True, nlakes=nlakes,
noutlets=noutlets,
print_stage=True,
print_flows=True,
packagedata=packagedata,
connectiondata=connectiondata,
outlets=outlets,
Expand Down
3 changes: 3 additions & 0 deletions autotest/test_gwt_dsp01.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,7 @@ def get_model(idx, dir):
chd = flopy.mf6.ModflowGwfchd(gwf, maxbound=len(c),
stress_period_data=c,
save_flows=False,
print_flows=True,
pname='CHD-1')

# output control
Expand Down Expand Up @@ -166,6 +167,7 @@ def get_model(idx, dir):
cncs = {0: [[(0, 0, 0), 1.0]]}
cnc = flopy.mf6.ModflowGwtcnc(gwt, maxbound=len(cncs),
stress_period_data=cncs,
print_flows=True,
save_flows=False,
pname='CNC-1')

Expand All @@ -174,6 +176,7 @@ def get_model(idx, dir):

# sources
ssm = flopy.mf6.ModflowGwtssm(gwt, sources=[[]],
print_flows=True,
filename='{}.ssm'.format(gwtname))

# output control
Expand Down
1 change: 1 addition & 0 deletions doc/ReleaseNotes/ReleaseNotes.tex
View file
Original file line number Diff line number Diff line change
Expand Up @@ -138,6 +138,7 @@ \section{History}
\item The LAK Package would accept negative user-input values for RAINFALL, EVAPORATION, RUNOFF, INFLOW, and WITHDRAWAL even though the user guide mentioned that negative values are not allowed for these flow terms. Error checks were added to ensure these values are specified as positive.
\item Added a storage term to the SFR Package binary output file. This term is always zero with the present implementation. An auxiliary variable, called VOLUME, is also written with the storage budget term. This term contains the calculated water volume in the reach.
\item Added additional error trapping in the MAW Package to catch divide by zero errors when calculating the saturated conductance for wells using the SKIN CONDEQN in connections where the cell transmissivity (the product of geometric mean of the horizontal hydraulic conductivity and cell thickness) and well transmissivity (the product of HK\_SKIN and screen thickness) is equal to one. Also add error trapping for well connections using the 1) SKIN CONDEQN where the contrast between the cell and well transmissivities are less than one and 2) SKIN and MEAN CONDEQN where the calculated connection saturated conductance is less than zero.
\item For the Lake Package, the outlet number was written as ID1 and ID2 for the TO-MVR record in the binary budget file. This has been changed so that the lake number of the connected outlet is written to ID1 and ID2. This change was implemented so that lake budgets can be calculated using the information in the lake budget file.
\end{itemize}

\underline{SOLUTION}
Expand Down
2 changes: 1 addition & 1 deletion doc/mf6io/gwf/binaryoutput.tex
View file
Original file line number Diff line number Diff line change
Expand Up @@ -483,7 +483,7 @@ \subsubsection{LAK, MAW, SFR, and UZF Packages}
\texttt{CONSTANT} & 6 & 1 / \texttt{nlakes} & Calculated flow to maintain constant stage for lake. The lake number is written to (\texttt{ID1}) and (\texttt{ID2}). \\
\texttt{EXT-OUTFLOW} & 6 & 1 / \texttt{nlakes} & Calculated outflow to external boundaries (is nonzero for lakes with outlets not connected to another lake). The lake number is written to (\texttt{ID1}) and (\texttt{ID2}). \\
\texttt{FROM-MVR} & 6 & 1 / \texttt{nlakes} & Calculated flow to lake from the MVR Package. Only saved if MVR Package is used in the LAK Package. The lake number is written to (\texttt{ID1}) and (\texttt{ID2}). \\
\texttt{TO-MVR} & 6 & 1 / \texttt{noutlets} & Calculated flow from a lake outlet to the MVR Package. Only saved if MVR Package is used in the LAK Package. The lake outlet number is written to (\texttt{ID1}) and (\texttt{ID2}). \\
\texttt{TO-MVR} & 6 & 1 / \texttt{noutlets} & Calculated flow from a lake outlet to the MVR Package. Only saved if MVR Package is used in the LAK Package. The lake number \texttt{LAKEIN} for the connected outlet is written to (\texttt{ID1}) and (\texttt{ID2}). \\
\texttt{AUXILIARY} & 6 & \texttt{naux}+1 / \texttt{nlakes} & Auxiliary variables, if specified in the LAK Package, are saved to this flow term. The first entry of the \texttt{DATA2D} column has a value of zero. The lake number is written to (\texttt{ID1}) and (\texttt{ID2}).
\label{table:binarylak}
\end{longtable}
Expand Down
146 changes: 119 additions & 27 deletions src/Exchange/GwfGwfExchange.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,8 @@ module GwfGwfExchangeModule
use ArrayHandlersModule, only: ExpandArray
use BaseModelModule, only: GetBaseModelFromList
use BaseExchangeModule, only: BaseExchangeType, AddBaseExchangeToList
use ConstantsModule, only: LENBOUNDNAME, NAMEDBOUNDFLAG
use ConstantsModule, only: LENBOUNDNAME, NAMEDBOUNDFLAG, LINELENGTH, &
TABCENTER, TABLEFT
use ListsModule, only: basemodellist
use NumericalExchangeModule, only: NumericalExchangeType
use NumericalModelModule, only: NumericalModelType
Expand All @@ -16,6 +17,7 @@ module GwfGwfExchangeModule
use SimModule, only: count_errors, store_error, &
store_error_unit, ustop
use BlockParserModule, only: BlockParserType
use TableModule, only: TableType, table_cr

implicit none

Expand Down Expand Up @@ -46,7 +48,13 @@ module GwfGwfExchangeModule
type(ObsType), pointer :: obs => null() ! observation object
character(len=LENBOUNDNAME), dimension(:), &
pointer, contiguous :: boundname => null() ! boundnames
!
! -- table objects
type(TableType), pointer :: outputtab1 => null()
type(TableType), pointer :: outputtab2 => null()

contains

procedure :: exg_df => gwf_gwf_df
procedure :: exg_ac => gwf_gwf_ac
procedure :: exg_mc => gwf_gwf_mc
Expand Down Expand Up @@ -830,7 +838,7 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
! ------------------------------------------------------------------------------
! -- modules
use ConstantsModule, only: DZERO, LENBUDTXT, LENPACKAGENAME
use TdisModule, only: kstp, kper
!use TdisModule, only: kstp, kper
! -- dummy
class(GwfExchangeType) :: this
integer(I4B), intent(inout) :: icnvg
Expand All @@ -841,23 +849,53 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
character(len=LENPACKAGENAME+4) :: packname1
character(len=LENPACKAGENAME+4) :: packname2
character(len=LENBUDTXT), dimension(1) :: budtxt
character(len=20) :: nodestr
integer(I4B) :: ntabrows
integer(I4B) :: nodeu
real(DP), dimension(2, 1) :: budterm
integer(I4B) :: i, n1, n2, n1u, n2u
integer(I4B) :: ibinun1, ibinun2
integer(I4B) :: ibdlbl
integer(I4B) :: icbcfl, ibudfl
real(DP) :: ratin, ratout, rrate, deltaqgnc
! -- formats
character(len=*), parameter :: fmttkk = &
"(1X,/1X,A,' PERIOD ',I0,' STEP ',I0)"
! ------------------------------------------------------------------------------
!
! -- initialize local variables
budtxt(1) = ' FLOW-JA-FACE'
packname1 = 'EXG '//this%name
packname1 = adjustr(packname1)
packname2 = 'EXG '//this%name
packname2 = adjustr(packname2)
!
! -- update output tables
if (this%iprflow /= 0) then
!
! -- update titles
if (this%gwfmodel1%oc%oc_save('BUDGET')) then
call this%outputtab1%set_title(packname1)
end if
if (this%gwfmodel2%oc%oc_save('BUDGET')) then
call this%outputtab2%set_title(packname2)
end if
!
! -- update maxbound of tables
ntabrows = 0
do i = 1, this%nexg
n1 = this%nodem1(i)
n2 = this%nodem2(i)
!
! -- If both cells are active then calculate flow rate
if (this%gwfmodel1%ibound(n1) /= 0 .and. &
this%gwfmodel2%ibound(n2) /= 0) then
ntabrows = ntabrows + 1
end if
end do
if (ntabrows > 0) then
call this%outputtab1%set_maxbound(ntabrows)
call this%outputtab2%set_maxbound(ntabrows)
end if
end if
!
! -- Print and write budget terms for model 1
!
! -- Set binary unit numbers for saving flows
Expand Down Expand Up @@ -886,7 +924,6 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
! Initialize accumulators
ratin = DZERO
ratout = DZERO
ibdlbl = 0
!
! -- Loop through all exchanges
do i = 1, this%nexg
Expand Down Expand Up @@ -917,12 +954,13 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
! -- Print the individual rates to model list files if requested
if(this%iprflow /= 0) then
if(this%gwfmodel1%oc%oc_save('BUDGET')) then
if(ibdlbl == 0) write(this%gwfmodel1%iout,fmttkk) packname1, &
kper, kstp
call this%gwfmodel1%dis%print_list_entry(i, n1, rrate, &
this%gwfmodel1%iout, bname)
endif
ibdlbl = 1
!
! -- set nodestr and write outputtab table
nodeu = this%gwfmodel1%dis%get_nodeuser(n1)
call this%gwfmodel1%dis%nodeu_to_string(nodeu, nodestr)
call this%outputtab1%print_list_entry(i, trim(adjustl(nodestr)), &
rrate, bname)
end if
endif
if(rrate < DZERO) then
ratout = ratout - rrate
Expand Down Expand Up @@ -974,7 +1012,6 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
! Initialize accumulators
ratin = DZERO
ratout = DZERO
ibdlbl = 0
!
! -- Loop through all exchanges
do i = 1, this%nexg
Expand Down Expand Up @@ -1005,12 +1042,13 @@ subroutine gwf_gwf_bd(this, icnvg, isuppress_output, isolnid)
! -- Print the individual rates to model list files if requested
if(this%iprflow /= 0) then
if(this%gwfmodel2%oc%oc_save('BUDGET')) then
if(ibdlbl == 0) write(this%gwfmodel2%iout,fmttkk) packname2, &
kper, kstp
call this%gwfmodel2%dis%print_list_entry(i, n2, -rrate, &
this%gwfmodel2%iout, bname)
endif
ibdlbl = 1
!
! -- set nodestr and write outputtab table
nodeu = this%gwfmodel2%dis%get_nodeuser(n2)
call this%gwfmodel2%dis%nodeu_to_string(nodeu, nodestr)
call this%outputtab2%print_list_entry(i, trim(adjustl(nodestr)), &
-rrate, bname)
end if
endif
if(rrate < DZERO) then
ratout = ratout - rrate
Expand Down Expand Up @@ -1743,6 +1781,26 @@ subroutine gwf_gwf_da(this)
call this%obs%obs_da()
deallocate(this%obs)
!
! -- arrays
call mem_deallocate(this%ihc)
call mem_deallocate(this%cl1)
call mem_deallocate(this%cl2)
call mem_deallocate(this%hwva)
call mem_deallocate(this%condsat)
deallocate(this%boundname)
!
! -- output table objects
if (associated(this%outputtab1)) then
call this%outputtab1%table_da()
deallocate(this%outputtab1)
nullify(this%outputtab1)
end if
if (associated(this%outputtab2)) then
call this%outputtab2%table_da()
deallocate(this%outputtab2)
nullify(this%outputtab2)
end if
!
! -- scalars
call mem_deallocate(this%icellavg)
call mem_deallocate(this%ivarcv)
Expand All @@ -1756,14 +1814,6 @@ subroutine gwf_gwf_da(this)
call mem_deallocate(this%inamedbound)
call mem_deallocate(this%satomega)
!
! -- arrays
call mem_deallocate(this%ihc)
call mem_deallocate(this%cl1)
call mem_deallocate(this%cl2)
call mem_deallocate(this%hwva)
call mem_deallocate(this%condsat)
deallocate(this%boundname)
!
! -- return
return
end subroutine gwf_gwf_da
Expand All @@ -1781,7 +1831,9 @@ subroutine allocate_arrays(this)
! -- dummy
class(GwfExchangeType) :: this
! -- local
character(len=LINELENGTH) :: text
character(len=LENORIGIN) :: origin
integer(I4B) :: ntabcol
! ------------------------------------------------------------------------------
!
! -- create the origin name
Expand All @@ -1804,6 +1856,46 @@ subroutine allocate_arrays(this)
endif
this%boundname(:) = ''
!
! -- allocate and initialize the output table
if (this%iprflow /= 0) then
!
! -- dimension table
ntabcol = 3
if (this%inamedbound > 0) then
ntabcol = ntabcol + 1
end if
!
! -- initialize the output table objects
! outouttab1
call table_cr(this%outputtab1, this%name, ' ')
call this%outputtab1%table_df(this%nexg, ntabcol, this%gwfmodel1%iout, &
transient=.TRUE.)
text = 'NUMBER'
call this%outputtab1%initialize_column(text, 10, alignment=TABCENTER)
text = 'CELLID'
call this%outputtab1%initialize_column(text, 20, alignment=TABLEFT)
text = 'RATE'
call this%outputtab1%initialize_column(text, 15, alignment=TABCENTER)
if (this%inamedbound > 0) then
text = 'NAME'
call this%outputtab1%initialize_column(text, 20, alignment=TABLEFT)
end if
! outouttab2
call table_cr(this%outputtab2, this%name, ' ')
call this%outputtab2%table_df(this%nexg, ntabcol, this%gwfmodel2%iout, &
transient=.TRUE.)
text = 'NUMBER'
call this%outputtab2%initialize_column(text, 10, alignment=TABCENTER)
text = 'CELLID'
call this%outputtab2%initialize_column(text, 20, alignment=TABLEFT)
text = 'RATE'
call this%outputtab2%initialize_column(text, 15, alignment=TABCENTER)
if (this%inamedbound > 0) then
text = 'NAME'
call this%outputtab2%initialize_column(text, 20, alignment=TABLEFT)
end if
end if
!
! -- return
return
end subroutine allocate_arrays
Expand Down
23 changes: 15 additions & 8 deletions src/Model/GroundWaterFlow/gwf3chd8.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -241,7 +241,7 @@ subroutine chd_bd(this, x, idvfl, icbcfl, ibudfl, icbcun, iprobs, &
! SPECIFICATIONS:
! ------------------------------------------------------------------------------
! -- modules
use TdisModule, only: kstp, kper, delt
use TdisModule, only: delt
use ConstantsModule, only: LENBOUNDNAME
use BudgetModule, only: BudgetType
! -- dummy
Expand All @@ -257,19 +257,21 @@ subroutine chd_bd(this, x, idvfl, icbcfl, ibudfl, icbcun, iprobs, &
integer(I4B), dimension(:), optional, intent(in) :: imap
integer(I4B), optional, intent(in) :: iadv
! -- local
character(len=20) :: nodestr
integer(I4B) :: nodeu
integer(I4B) :: i, node, ibinun, n2
real(DP) :: rrate, chin, chout, q
integer(I4B) :: ibdlbl, naux, ipos
integer(I4B) :: naux, ipos
! -- for observations
character(len=LENBOUNDNAME) :: bname
! -- formats
character(len=*), parameter :: fmttkk = &
"(1X,/1X,A,' PERIOD ',I0,' STEP ',I0)"
! ------------------------------------------------------------------------------
!
! -- initialize local variables
chin = DZERO
chout = DZERO
ibdlbl = 0
!
! -- Set unit number for binary output
if(this%ipakcb < 0) then
Expand All @@ -291,6 +293,11 @@ subroutine chd_bd(this, x, idvfl, icbcfl, ibudfl, icbcun, iprobs, &
!
! -- If no boundaries, skip flow calculations.
if(this%nbound > 0) then
!
! -- reset size of table
if (this%iprflow /= 0) then
call this%outputtab%set_maxbound(this%nbound)
end if
!
! -- Loop through each boundary calculating flow.
do i = 1, this%nbound
Expand Down Expand Up @@ -327,11 +334,11 @@ subroutine chd_bd(this, x, idvfl, icbcfl, ibudfl, icbcun, iprobs, &
! -- Print the individual rates if requested(this%iprflow<0)
if (ibudfl /= 0) then
if(this%iprflow /= 0) then
if(ibdlbl == 0) write(this%iout,fmttkk) &
this%text // ' (' // trim(this%name) // ')', kper, kstp
call this%dis%print_list_entry(i, node, rrate, this%iout, &
bname)
ibdlbl=1
!
! -- set nodestr and write outputtab table
nodeu = this%dis%get_nodeuser(node)
call this%dis%nodeu_to_string(nodeu, nodestr)
call this%outputtab%print_list_entry(i, nodestr, rrate, bname)
end if
end if
!
Expand Down
7 changes: 4 additions & 3 deletions src/Model/GroundWaterFlow/gwf3lak8.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4178,7 +4178,7 @@ subroutine lak_ot(this, kstp, kper, iout, ihedfl, ibudfl)
this%lakename(n), n, q, ALIGNMENT=TABLEFT)
end if
call UWWORD(line, iloc, 6, TABINTEGER, text, n, q, SEP=' ')
call UWWORD(line, iloc, 11, TABREAL, text, n, this %xnewpak(n))
call UWWORD(line, iloc, 11, TABREAL, text, n, this%xnewpak(n))
write(iout, '(1X,A)') line(1:iloc)
end do
end if
Expand Down Expand Up @@ -5803,7 +5803,7 @@ subroutine lak_setup_budobj(this)
this%name_model, &
this%name, &
maxlist, .false., .false., &
naux)
naux, ordered_id1=.false.)
end if
!
! --
Expand Down Expand Up @@ -5992,11 +5992,12 @@ subroutine lak_fill_budobj(this)
idx = idx + 1
call this%budobj%budterm(idx)%reset(this%noutlets)
do n = 1, this%noutlets
n1 = this%lakein(n)
q = this%pakmvrobj%get_qtomvr(n)
if (q > DZERO) then
q = -q
end if
call this%budobj%budterm(idx)%update_term(n, n, q)
call this%budobj%budterm(idx)%update_term(n1, n1, q)
end do

end if
Expand Down
Loading

0 comments on commit f55cf79

Please sign in to comment.