16 use,
intrinsic :: iso_fortran_env, only: iostat_end
25 character(len=LENTIMESERIESNAME),
public :: name =
''
27 integer(I4B),
private :: inunit = 0
28 integer(I4B),
private :: iout = 0
30 real(dp),
private :: sfac =
done
31 character(len=LINELENGTH),
private :: datafile =
''
32 logical,
private :: autodeallocate = .true.
33 type(
listtype),
pointer,
private :: list => null()
34 character(len=LENMODELNAME) :: modelname
62 character(len=*),
intent(in) :: filename
66 10
format(
'Error: Time-array-series file "', a,
'" does not exist.')
70 allocate (newtas%list)
73 inquire (file=filename, exist=lex)
75 write (
errmsg, 10) trim(filename)
78 newtas%datafile = filename
88 subroutine tas_init(this, fname, modelname, iout, tasname, autoDeallocate)
91 character(len=*),
intent(in) :: fname
92 character(len=*),
intent(in) :: modelname
93 integer(I4B),
intent(in) :: iout
94 character(len=*),
intent(inout) :: tasname
95 logical,
optional,
intent(in) :: autoDeallocate
97 integer(I4B) :: istatus
99 integer(I4B) :: inunit
100 character(len=40) :: keyword, keyvalue
101 logical :: found, continueread, endOfBlock
104 if (
present(autodeallocate)) this%autoDeallocate = autodeallocate
105 this%dataFile = fname
109 this%modelname = modelname
115 call openfile(inunit, 0, fname,
'TAS6')
118 call this%parser%Initialize(this%inunit, this%iout)
121 continueread = .false.
125 call this%parser%GetBlock(
'ATTRIBUTES', found, ierr, &
126 supportopenclose=.true.)
127 if (.not. found)
then
128 errmsg =
'Error: Attributes block not found in file: '// &
131 call this%parser%StoreErrorUnit()
137 call this%parser%GetNextLine(endofblock)
141 call this%parser%GetStringCaps(keyword)
144 call this%parser%GetStringCaps(keyvalue)
145 select case (keyword)
150 select case (keyvalue)
156 errmsg =
'Unknown interpolation method: "'//trim(keyvalue)//
'"'
158 call this%parser%StoreErrorUnit()
160 case (
'AUTODEALLOCATE')
161 this%autoDeallocate = (keyvalue ==
'TRUE')
163 read (keyvalue, *, iostat=istatus) this%sfac
164 if (istatus /= 0)
then
165 errmsg =
'Error reading numeric SFAC value from "'//trim(keyvalue) &
168 call this%parser%StoreErrorUnit()
171 errmsg =
'Unknown option found in ATTRIBUTES block: "'// &
174 call this%parser%StoreErrorUnit()
179 if (this%Name ==
'')
then
180 errmsg =
'Name not specified for time array series in file: '// &
183 call this%parser%StoreErrorUnit()
186 errmsg =
'Interpolation method not specified for time'// &
187 ' array series in file: '//trim(this%dataFile)
189 call this%parser%StoreErrorUnit()
194 errmsg =
'Error(s) encountered initializing time array series from file: ' &
195 //trim(this%dataFile)
197 call this%parser%StoreErrorUnit()
201 if (.not. this%read_next_array())
then
202 errmsg =
'Error encountered reading time-array data from file: '// &
205 call this%parser%StoreErrorUnit()
218 integer(I4B),
intent(in) :: nvals
219 real(DP),
dimension(nvals),
intent(inout) :: values
220 real(DP),
intent(in) :: time0
221 real(DP),
intent(in) :: time1
226 timediff = time1 - time0
227 if (timediff > 0)
then
228 call this%get_integrated_values(nvals, values, time0, time1)
230 values(i) = values(i) / timediff
234 call this%get_values_at_time(nvals, values, time0)
262 real(DP),
intent(in) :: time
266 real(DP) :: time0, time1
270 type(
timearraytype),
pointer :: ta => null(), ta0 => null(), ta1 => null()
271 class(*),
pointer :: obj
276 if (
associated(this%list%firstNode))
then
277 currnode => this%list%firstNode
283 if (
associated(currnode))
then
284 if (
associated(currnode%nextNode))
then
285 obj => currnode%nextNode%GetItem()
287 if (ta%taTime <= time)
then
288 currnode => currnode%nextNode
294 if (.not. this%read_next_array())
exit
301 if (
associated(currnode))
then
305 obj => node0%GetItem()
308 do while (time0 > time)
309 if (
associated(node0%prevNode))
then
310 node0 => node0%prevNode
311 obj => node0%GetItem()
321 obj => node1%GetItem()
324 do while (time1 < time)
325 if (
associated(node1%nextNode))
then
326 node1 => node1%nextNode
327 obj => node1%GetItem()
332 if (.not. this%read_next_array())
then
341 if (time0 <= time) taearlier => ta0
342 if (time1 >= time) talater => ta1
358 integer(I4B) :: i, ierr, istart, istat, istop, lloc, nrow, ncol, &
360 logical :: lopen, isfound
362 character(len=LENMEMPATH) :: mempath
363 integer(I4B),
dimension(:),
contiguous,
pointer :: mshape
373 mempath =
create_mem_path(component=this%modelname, subcomponent=
'DIS')
379 if (
size(mshape) == 2)
then
380 nodesperlayer = mshape(2)
383 else if (
size(mshape) == 3)
then
384 nodesperlayer = mshape(2) * mshape(3)
388 errmsg =
'Time array series is not supported for selected &
389 &discretization type.'
391 call this%parser%StoreErrorUnit()
395 inquire (unit=this%inunit, opened=lopen)
400 call this%parser%GetBlock(
'TIME', isfound, ierr, &
401 supportopenclose=.false.)
403 ta%taTime = this%parser%GetDouble()
405 call readarray(this%parser%iuactive, ta%taArray, this%Name, &
406 size(mshape), ncol, nrow, 1, nodesperlayer, &
410 do i = 1, nodesperlayer
411 ta%taArray(i) = ta%taArray(i) * this%sfac
419 call this%parser%terminateblock()
433 integer(I4B),
intent(in) :: nvals
434 real(DP),
dimension(nvals),
intent(inout) :: values
435 real(DP),
intent(in) :: time
437 integer(I4B) :: i, ierr
438 real(DP) :: ratio, time0, time1, timediff, timediffi, val0, val1, &
443 10
format(
'Error getting array at time ', g10.3, &
444 ' for time-array series "', a,
'"')
447 call this%get_surrounding_records(time, taearlier, talater)
448 if (
associated(taearlier))
then
449 if (
associated(talater))
then
453 values = taearlier%taArray
454 elseif (this%iMethod ==
linear)
then
456 time0 = taearlier%taTime
457 time1 = talater%tatime
458 timediff = time1 - time0
459 timediffi = time - time0
460 if (timediff > 0)
then
461 ratio = timediffi / timediff
468 val0 = taearlier%taArray(i)
469 val1 = talater%taArray(i)
470 valdiff = val1 - val0
471 values(i) = val0 + (ratio * valdiff)
477 if (
is_close(taearlier%taTime, time))
then
478 values = taearlier%taArray
483 values = taearlier%taArray
490 if (
associated(talater))
then
491 if (
is_close(talater%taTime, time))
then
492 values = talater%taArray
505 write (
errmsg, 10) time, trim(this%Name)
521 integer(I4B),
intent(in) :: nvals
522 real(DP),
dimension(nvals),
intent(inout) :: values
523 real(DP),
intent(in) :: time0
524 real(DP),
intent(in) :: time1
527 real(DP) :: area, currTime, nextTime, ratio0, ratio1, t0, &
528 t01, t1, timediff,
value, value0, value1, valuediff
531 type(
listnodetype),
pointer :: currNode => null(), nextnode => null()
532 type(
timearraytype),
pointer :: currRecord => null(), nextrecord => null()
533 class(*),
pointer :: currObj => null(), nextobj => null()
535 10
format(
'Error encountered while performing integration', &
536 ' for time-array series "', a,
'" for time interval: ', &
537 g12.5,
' to ', g12.5)
543 call this%get_latest_preceding_node(time0, precnode)
544 if (
associated(precnode))
then
546 do while (.not. ldone)
547 currobj => currnode%GetItem()
549 currtime = currrecord%taTime
550 if (currtime < time1)
then
551 if (.not.
associated(currnode%nextNode))
then
553 if (.not. this%read_next_array())
then
554 write (
errmsg, 10) trim(this%Name), time0, time1
559 if (
associated(currnode%nextNode))
then
560 nextnode => currnode%nextNode
561 nextobj => nextnode%GetItem()
563 nexttime = nextrecord%taTime
566 if (currtime >= time0)
then
571 if (nexttime <= time1)
then
579 select case (this%iMethod)
583 value0 = currrecord%taArray(i)
586 values(i) = values(i) + area
591 timediff = nexttime - currtime
592 ratio0 = (t0 - currtime) / timediff
593 ratio1 = (t1 - currtime) / timediff
594 valuediff = nextrecord%taArray(i) - currrecord%taArray(i)
595 value0 = currrecord%taArray(i) + ratio0 * valuediff
596 value1 = currrecord%taArray(i) + ratio1 * valuediff
597 area = 0.5d0 * t01 * (value0 + value1)
599 values(i) = values(i) + area
603 write (
errmsg, 10) trim(this%Name), time0, time1
605 call store_error(
'(Probable programming error)', terminate=.true.)
613 if (t1 >= time1)
then
616 if (.not.
associated(currnode%nextNode))
then
618 if (.not. this%read_next_array())
then
619 write (
errmsg, 10) trim(this%Name), time0, time1
621 call this%parser%StoreErrorUnit()
624 if (
associated(currnode%nextNode))
then
625 currnode => currnode%nextNode
627 write (
errmsg, 10) trim(this%Name), time0, time1
629 call store_error(
'(Probable programming error)', terminate=.true.)
635 if (this%autoDeallocate)
then
636 if (
associated(precnode))
then
637 if (
associated(precnode%prevNode))
then
638 call this%DeallocateBackward(precnode%prevNode)
659 class(*),
pointer :: obj => null()
661 if (
associated(fromnode))
then
663 if (
associated(fromnode%nextNode))
then
664 this%list%firstNode => fromnode%nextNode
666 this%list%firstNode => null()
670 do while (
associated(current))
671 prev => current%prevNode
672 obj => current%GetItem()
677 call this%list%RemoveNode(current, .true.)
693 real(DP),
intent(in) :: time
701 class(*),
pointer :: obj => null()
704 if (
associated(this%list%firstNode))
then
705 currnode => this%list%firstNode
708 &get_latest_preceding_node', &
716 if (
associated(currnode))
then
717 if (
associated(currnode%nextNode))
then
718 obj => currnode%nextNode%GetItem()
720 if (ta%taTime < time .or.
is_close(ta%taTime, time))
then
721 currnode => currnode%nextNode
727 if (.not. this%read_next_array())
exit
734 if (
associated(currnode))
then
738 obj => node0%GetItem()
741 do while (time0 > time)
742 if (
associated(node0%prevNode))
then
743 node0 => node0%prevNode
744 obj => node0%GetItem()
753 if (time0 <= time) tslnode => node0
769 n = this%list%Count()
776 call this%list%Clear(.true.)
777 deallocate (this%list)
789 class(*),
pointer,
intent(inout) :: obj
794 if (.not.
associated(obj))
return
809 type(
listtype),
intent(inout) :: list
810 integer,
intent(in) :: indx
814 class(*),
pointer :: obj
816 obj => list%GetItem(indx)
This module contains block parser methods.
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 lentimeseriesname
maximum length of a time series name
@ undefined
undefined interpolation
@ linear
linear interpolation
@ stepwise
stepwise interpolation
real(dp), parameter dzero
real constant zero
integer(i4b), parameter lenmempath
maximum length of the memory path
real(dp), parameter done
real constant 1
This module defines variable data types.
pure logical function, public is_close(a, b, rtol, atol, symmetric)
Check if a real value is approximately equal to another.
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
This module contains simulation methods.
subroutine, public store_error(msg, terminate)
Store an error message.
integer(i4b) function, public count_errors()
Return number of errors.
subroutine, public store_error_unit(iunit, terminate)
Store the file unit number.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
type(timearraytype) function, pointer, public castastimearraytype(obj)
Cast an unlimited polymorphic object as TimeArrayType.
type(timearraytype) function, pointer, public gettimearrayfromlist(list, indx)
Retrieve a time array from a list.
subroutine, public constructtimearray(newTa, modelname)
Construct time array.
subroutine, public addtimearraytolist(list, timearray)
Add a time array to a to list.
subroutine, public constructtimearrayseries(newTas, filename)
Allocate a new TimeArraySeriesType object.
subroutine tas_init(this, fname, modelname, iout, tasname, autoDeallocate)
Initialize the time array series.
subroutine get_latest_preceding_node(this, time, tslNode)
Return pointer to ListNodeType object for the node representing the latest preceding time in the time...
subroutine getaveragevalues(this, nvals, values, time0, time1)
Populate an array time-weighted average value for a specified time span.
subroutine tas_da(this)
Deallocate memory.
subroutine get_values_at_time(this, nvals, values, time)
Return an array of values for a specified time, same units as time-series values.
subroutine deallocatebackward(this, fromNode)
Deallocate fromNode and all previous nodes in list; reassign firstNode.
subroutine get_surrounding_records(this, time, taEarlier, taLater)
Get surrounding records.
type(timearrayseriestype) function, pointer, public gettimearrayseriesfromlist(list, indx)
Get time array from list.
subroutine get_integrated_values(this, nvals, values, time0, time1)
Populates an array with integrated values for a specified time span.
type(timearrayseriestype) function, pointer, public castastimearrayseriestype(obj)
Cast an unlimited polymorphic object as class(TimeArraySeriesType)
integer(i4b) function getinunit(this)
Return unit number.
logical function read_next_array(this)
Read next time array from input file and append to list.
A generic heterogeneous doubly-linked list.