24 integer(I4B),
public :: iout = 0
26 character(len=LENMODELNAME) :: modelname
28 type(
listtype),
pointer,
private :: boundtaslinks => null()
29 character(len=LINELENGTH),
allocatable,
dimension(:) :: tasfiles
31 character(len=LENTIMESERIESNAME),
allocatable,
dimension(:) :: tasnames
61 character(len=*),
intent(in) :: modelname
62 integer(I4B),
intent(in) :: iout
64 if (
present(dis))
then
68 this%modelname = modelname
70 allocate (this%boundTasLinks)
71 allocate (this%tasfiles(0))
84 integer(I4B) :: nfiles
89 nfiles =
size(this%tasfiles)
90 allocate (this%taslist(nfiles))
91 allocate (this%tasnames(nfiles))
95 tasptr => this%taslist(i)
96 call tasptr%tas_init(this%tasfiles(i), this%modelname, &
97 this%iout, this%tasnames(i))
114 integer(I4B) :: i, j, nlinks, nvals, isize1, isize2, inunit
115 real(DP) :: begintime, endtime
117 character(len=*),
parameter :: fmt5 = &
118 "(/,'Time-array-series controlled arrays in stress period ', &
119 &i0, ', time step ', i0, ':')"
120 10
format(
'"', a,
'" package: ', a,
' array obtained from time-array series "', &
125 endtime = begintime +
delt
130 if (
associated(this%boundTasLinks))
then
131 nlinks = this%boundTasLinks%Count()
134 if (taslink%Iprpak == 1 .and. i == 1)
then
137 if (taslink%UseDefaultProc)
then
138 timearrayseries => taslink%timeArraySeries
139 nvals =
size(taslink%BndArray)
142 call timearrayseries%GetAverageValues(nvals, taslink%BndArray, &
146 if (taslink%ConvertFlux)
then
147 call this%tasmgr_convert_flux(taslink)
152 if (taslink%Iprpak == 1)
then
153 write (this%iout, 10) trim(taslink%PackageName), &
154 trim(taslink%Text), &
155 trim(taslink%timeArraySeries%Name)
158 if (i == nlinks)
then
159 write (this%iout,
'()')
167 if (taslink%UseDefaultProc)
then
168 if (
associated(taslink%RMultArray))
then
169 isize1 =
size(taslink%BndArray)
170 isize2 =
size(taslink%RMultArray)
171 if (isize1 == isize2 .and. isize1 == nvals)
then
173 taslink%BndArray(j) = taslink%BndArray(j) * taslink%RMultArray(j)
176 errmsg =
'Size mismatch between boundary and multiplier arrays'// &
177 ' using time-array series: '// &
178 trim(taslink%TimeArraySeries%Name)
180 inunit = taslink%TimeArraySeries%GetInunit()
203 n = this%boundTasLinks%Count()
210 do i = 1,
size(this%taslist)
211 call this%taslist(i)%da()
215 call this%boundTasLinks%Clear(.true.)
216 deallocate (this%boundTasLinks)
217 deallocate (this%tasfiles)
220 deallocate (this%taslist)
221 deallocate (this%tasnames)
225 this%boundTasLinks => null()
238 character(len=*),
intent(in) :: fname
243 indx =
size(this%tasfiles)
244 this%tasfiles(indx) = fname
258 character(len=*),
intent(in) :: pkgName
260 integer(I4B) :: i, j, nlinks
264 nlinks = this%boundTasLinks%Count()
267 if (
associated(taslink))
then
268 do j = 1,
size(taslink%BndArray)
269 taslink%BndArray(j) =
dzero
275 if (
associated(this%boundTasLinks))
then
277 nlinks = this%boundTasLinks%Count()
280 if (
associated(taslink))
then
282 call this%boundTasLinks%RemoveNode(i, .true.)
294 tasName, text, convertFlux, nodelist, inunit)
297 character(len=*),
intent(in) :: pkgName
298 real(DP),
dimension(:),
pointer :: bndArray
299 integer(I4B),
intent(in) :: iprpak
300 character(len=*),
intent(in) :: tasName
301 character(len=*),
intent(in) :: text
302 logical,
intent(in) :: convertFlux
303 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(in) :: nodelist
304 integer(I4B),
intent(in) :: inunit
306 integer(I4B) :: i, nfiles, iloc
307 character(LINELENGTH) :: ermsg
312 nfiles =
size(this%tasnames)
315 if (this%tasnames(i) == tasname)
then
321 ermsg =
'Error: Time-array series "'//trim(tasname)//
'" not found.'
325 tasptr => this%taslist(iloc)
330 pkgname, bndarray, iprpak, &
332 newtaslink%ConvertFlux = convertflux
333 newtaslink%nodelist => nodelist
336 call this%tasmgr_add_link(newtaslink)
347 integer(I4B),
intent(in) :: indx
353 if (
associated(this%boundTasLinks))
then
369 if (
associated(this%boundtaslinks))
then
389 integer(I4B) :: i, n, noder
392 if (.not. (
associated(this%dis) .and. &
393 associated(taslink%nodelist)))
then
394 errmsg =
'Programming error. Cannot convert flux. Verify that '&
395 &
'a valid DIS instance and nodelist were provided.'
400 n =
size(taslink%BndArray)
402 noder = taslink%nodelist(i)
404 area = this%dis%get_area(noder)
405 taslink%BndArray(i) = taslink%BndArray(i) * area
This module contains simulation constants.
integer(i4b), parameter linelength
maximum length of a standard line
integer(i4b), parameter lenmodelname
maximum length of the model name
integer(i4b), parameter lenhugeline
maximum length of a huge line
integer(i4b), parameter lentimeseriesname
maximum length of a time series name
real(dp), parameter dzero
real constant zero
This module defines variable data types.
This module contains simulation methods.
subroutine, public store_error(msg, terminate)
Store an error message.
subroutine, public store_error_unit(iunit, terminate)
Store the file unit number.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
real(dp), pointer, public totimc
simulation time at start of time step
integer(i4b), pointer, public kstp
current time step number
integer(i4b), pointer, public kper
current stress period number
real(dp), pointer, public delt
length of the current time step
type(timearrayserieslinktype) function, pointer, public gettimearrayserieslinkfromlist(list, idx)
Get time-array series from a list and return.
subroutine, public addtimearrayserieslinktolist(list, tasLink)
Add time-array series to list.
subroutine, public constructtimearrayserieslink(newTasLink, timeArraySeries, pkgName, bndArray, iprpak, text)
Construct a time series of arrays that are linked.
subroutine tasmanager_df(this)
Define the time-array series manager.
subroutine tasmgr_add_link(this, tasLink)
Add a time arrays series link.
subroutine tasmgr_ad(this)
Time step (or subtime step) advance.
subroutine reset(this, pkgName)
Zero out arrays that are represented with time series.
subroutine maketaslink(this, pkgName, bndArray, iprpak, tasName, text, convertFlux, nodelist, inunit)
Make link from time-array series to package array.
subroutine tasmgr_convert_flux(this, tasLink)
Convert the array from a flux to a flow rate by multiplying by the cell area.
subroutine, public tasmanager_cr(this, dis, modelname, iout)
Create the time-array series manager.
integer(i4b) function countlinks(this)
Count number of links in the boundtaslinks list.
subroutine add_tasfile(this, fname)
Add a time-array series file.
type(timearrayserieslinktype) function, pointer getlink(this, indx)
Get link from the boundtaslinks list.
subroutine tasmgr_da(this)
Deallocate.
A generic heterogeneous doubly-linked list.