MODFLOW 6  version 6.6.0.dev0
USGS Modular Hydrologic Model
timeseriesmodule Module Reference

Data Types

type  timeseriestype
 
type  timeseriesfiletype
 
type  timeseriescontainertype
 

Functions/Subroutines

subroutine, public constructtimeseriesfile (newTimeSeriesFile)
 Construct time series file. More...
 
type(timeseriesfiletype) function, pointer castastimeseriesfiletype (obj)
 Cast an unlimited polymorphic object as class(TimeSeriesFileType) More...
 
type(timeseriesfiletype) function, pointer, public castastimeseriesfileclass (obj)
 Cast an unlimited polymorphic object as class(TimeSeriesFileType) More...
 
subroutine, public addtimeseriesfiletolist (list, tsfile)
 Add time series file to list. More...
 
type(timeseriesfiletype) function, pointer, public gettimeseriesfilefromlist (list, idx)
 Get time series from list. More...
 
logical function, public sametimeseries (ts1, ts2)
 Compare two time series; if they are identical, return true. More...
 
real(dp) function getvalue (this, time0, time1, extendToEndOfSimulation)
 Get time series value. More...
 
subroutine initialize_time_series (this, tsfile, name, autoDeallocate)
 Initialize time series. More...
 
subroutine get_surrounding_records (this, time, tsrecEarlier, tsrecLater)
 Get surrounding records. More...
 
subroutine get_surrounding_nodes (this, time, nodeEarlier, nodeLater)
 Get surrounding nodes. More...
 
logical function read_next_record (this)
 Read next record. More...
 
real(dp) function get_value_at_time (this, time, extendToEndOfSimulation)
 Get value for a time. More...
 
real(dp) function get_integrated_value (this, time0, time1, extendToEndOfSimulation)
 Get integrated value. More...
 
real(dp) function get_average_value (this, time0, time1, extendToEndOfSimulation)
 Get average value. More...
 
subroutine get_latest_preceding_node (this, time, tslNode)
 Get latest preceding node. More...
 
subroutine ts_da (this)
 Deallocate. More...
 
subroutine addtimeseriesrecord (this, tsr)
 Add ts record. More...
 
type(timeseriesrecordtype) function, pointer getcurrenttimeseriesrecord (this)
 Get current ts record. More...
 
type(timeseriesrecordtype) function, pointer getprevioustimeseriesrecord (this)
 Get previous ts record. More...
 
type(timeseriesrecordtype) function, pointer getnexttimeseriesrecord (this)
 Get next ts record. More...
 
type(timeseriesrecordtype) function, pointer gettimeseriesrecord (this, time, epsi)
 Get ts record. More...
 
subroutine reset (this)
 Reset. More...
 
subroutine inserttsr (this, tsr)
 Insert a time series record. More...
 
double precision function findlatesttime (this, readToEnd)
 Find latest time. More...
 
subroutine clear (this, destroy)
 Clear the list of time series records. More...
 
integer(i4b) function count (this)
 Count number of time series. More...
 
type(timeseriestype) function, pointer gettimeseries (this, indx)
 Get time series. More...
 
subroutine initializetsfile (this, filename, iout, autoDeallocate)
 Open time-series tsfile file and read options and first record, which may contain data to define multiple time series. More...
 
logical function read_tsfile_line (this)
 Read time series file line. More...
 
subroutine tsf_da (this)
 Deallocate memory. More...
 

Function/Subroutine Documentation

◆ addtimeseriesfiletolist()

subroutine, public timeseriesmodule::addtimeseriesfiletolist ( type(listtype), intent(inout)  list,
class(timeseriesfiletype), intent(inout), pointer  tsfile 
)

Definition at line 139 of file TimeSeries.f90.

140  ! -- dummy
141  type(ListType), intent(inout) :: list
142  class(TimeSeriesFileType), pointer, intent(inout) :: tsfile
143  ! -- local
144  class(*), pointer :: obj => null()
145  !
146  obj => tsfile
147  call list%Add(obj)
Here is the caller graph for this function:

◆ addtimeseriesrecord()

subroutine timeseriesmodule::addtimeseriesrecord ( class(timeseriestype this,
type(timeseriesrecordtype), intent(inout), pointer  tsr 
)
private

Definition at line 808 of file TimeSeries.f90.

809  ! -- dummy
810  class(TimeSeriesType) :: this
811  type(TimeSeriesRecordType), pointer, intent(inout) :: tsr
812  ! -- local
813  class(*), pointer :: obj => null()
814  !
815  obj => tsr
816  call this%list%Add(obj)

◆ castastimeseriesfileclass()

type(timeseriesfiletype) function, pointer, public timeseriesmodule::castastimeseriesfileclass ( class(*), intent(inout), pointer  obj)

Definition at line 122 of file TimeSeries.f90.

123  ! -- dummy
124  class(*), pointer, intent(inout) :: obj
125  ! -- return
126  type(TimeSeriesFileType), pointer :: res
127  !
128  res => null()
129  if (.not. associated(obj)) return
130  !
131  select type (obj)
132  class is (timeseriesfiletype)
133  res => obj
134  end select
Here is the caller graph for this function:

◆ castastimeseriesfiletype()

type(timeseriesfiletype) function, pointer timeseriesmodule::castastimeseriesfiletype ( class(*), intent(inout), pointer  obj)
private

Definition at line 105 of file TimeSeries.f90.

106  ! -- dummy
107  class(*), pointer, intent(inout) :: obj
108  ! -- return
109  type(TimeSeriesFileType), pointer :: res
110  !
111  res => null()
112  if (.not. associated(obj)) return
113  !
114  select type (obj)
115  type is (timeseriesfiletype)
116  res => obj
117  end select
Here is the caller graph for this function:

◆ clear()

subroutine timeseriesmodule::clear ( class(timeseriestype), intent(inout)  this,
logical, intent(in), optional  destroy 
)
private

Definition at line 1018 of file TimeSeries.f90.

1019  ! -- dummy
1020  class(TimeSeriesType), intent(inout) :: this
1021  logical, optional, intent(in) :: destroy
1022  !
1023  call this%list%Clear(destroy)

◆ constructtimeseriesfile()

subroutine, public timeseriesmodule::constructtimeseriesfile ( type(timeseriesfiletype), intent(inout), pointer  newTimeSeriesFile)

Definition at line 95 of file TimeSeries.f90.

96  ! -- dummy
97  type(TimeSeriesFileType), pointer, intent(inout) :: newTimeSeriesFile
98  !
99  allocate (newtimeseriesfile)
100  allocate (newtimeseriesfile%parser)
Here is the caller graph for this function:

◆ count()

integer(i4b) function timeseriesmodule::count ( class(timeseriesfiletype this)
private

Definition at line 1030 of file TimeSeries.f90.

1031  ! -- return
1032  integer(I4B) :: Count
1033  ! -- dummy
1034  class(TimeSeriesFileType) :: this
1035  !
1036  if (associated(this%timeSeries)) then
1037  count = size(this%timeSeries)
1038  else
1039  count = 0
1040  end if

◆ findlatesttime()

double precision function timeseriesmodule::findlatesttime ( class(timeseriestype), intent(inout)  this,
logical, intent(in), optional  readToEnd 
)
private

Definition at line 991 of file TimeSeries.f90.

992  ! -- dummy
993  class(TimeSeriesType), intent(inout) :: this
994  logical, intent(in), optional :: readToEnd
995  ! -- local
996  integer :: nrecords
997  type(TimeSeriesRecordType), pointer :: tsr
998  class(*), pointer :: obj => null()
999  ! -- return
1000  double precision :: endtime
1001  !
1002  ! -- If the caller requested the very last time in the series (readToEnd is true), check that we have first read all records
1003  if (present(readtoend)) then
1004  if (readtoend) then
1005  do while (this%read_next_record())
1006  end do
1007  end if
1008  end if
1009  !
1010  nrecords = this%list%Count()
1011  obj => this%list%GetItem(nrecords)
1012  tsr => castastimeseriesrecordtype(obj)
1013  endtime = tsr%tsrTime
Here is the call graph for this function:

◆ get_average_value()

real(dp) function timeseriesmodule::get_average_value ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time0,
real(dp), intent(in)  time1,
logical, intent(in)  extendToEndOfSimulation 
)
private

Return a time-weighted average value for a specified time span. Units: (ts-value-unit)

Definition at line 697 of file TimeSeries.f90.

698  ! -- return
699  real(DP) :: get_average_value
700  ! -- dummy
701  class(TimeSeriesType), intent(inout) :: this
702  real(DP), intent(in) :: time0
703  real(DP), intent(in) :: time1
704  logical, intent(in) :: extendToEndOfSimulation
705  ! -- local
706  real(DP) :: timediff, value, valueIntegrated
707  !
708  timediff = time1 - time0
709  if (timediff > 0) then
710  valueintegrated = this%get_integrated_value(time0, time1, &
711  extendtoendofsimulation)
712  if (this%iMethod == linearend) then
713  value = valueintegrated
714  else
715  value = valueintegrated / timediff
716  end if
717  else
718  ! -- time0 and time1 are the same
719  value = this%get_value_at_time(time0, extendtoendofsimulation)
720  end if
721  get_average_value = value

◆ get_integrated_value()

real(dp) function timeseriesmodule::get_integrated_value ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time0,
real(dp), intent(in)  time1,
logical, intent(in)  extendToEndOfSimulation 
)
private

Return an integrated value for a specified time span. Units: (ts-value-unit)*time

Definition at line 561 of file TimeSeries.f90.

562  ! -- return
563  real(DP) :: get_integrated_value
564  ! -- dummy
565  class(TimeSeriesType), intent(inout) :: this
566  real(DP), intent(in) :: time0
567  real(DP), intent(in) :: time1
568  logical, intent(in) :: extendToEndOfSimulation
569  ! -- local
570  real(DP) :: area, currTime, nextTime, ratio0, ratio1, t0, t01, t1, &
571  timediff, value, value0, value1, valuediff, currVal, nextVal
572  logical :: ldone, lprocess
573  type(ListNodeType), pointer :: tslNodePreceding => null()
574  type(ListNodeType), pointer :: currNode => null(), nextnode => null()
575  type(TimeSeriesRecordType), pointer :: currRecord => null()
576  type(TimeSeriesRecordType), pointer :: nextRecord => null()
577  class(*), pointer :: currObj => null(), nextobj => null()
578  ! -- formats
579 10 format('Error encountered while performing integration', &
580  ' for time series "', a, '" for time interval: ', g12.5, ' to ', g12.5)
581  !
582  value = dzero
583  ldone = .false.
584  t1 = -done
585  call this%get_latest_preceding_node(time0, tslnodepreceding)
586  if (associated(tslnodepreceding)) then
587  currnode => tslnodepreceding
588  do while (.not. ldone)
589  currobj => currnode%GetItem()
590  currrecord => castastimeseriesrecordtype(currobj)
591  currtime = currrecord%tsrTime
592  if (is_close(currtime, time1)) then
593  ! Current node time = time1 so should be ldone
594  ldone = .true.
595  elseif (currtime < time1) then
596  if (.not. associated(currnode%nextNode)) then
597  ! -- try to read the next record
598  if (.not. this%read_next_record()) then
599  if (.not. extendtoendofsimulation) then
600  write (errmsg, 10) trim(this%Name), time0, time1
601  call store_error(errmsg, terminate=.true.)
602  end if
603  end if
604  end if
605  !
606  currval = currrecord%tsrValue
607  lprocess = .false.
608  if (associated(currnode%nextNode)) then
609  nextnode => currnode%nextNode
610  nextobj => nextnode%GetItem()
611  nextrecord => castastimeseriesrecordtype(nextobj)
612  nexttime = nextrecord%tsrTime
613  nextval = nextrecord%tsrValue
614  lprocess = .true.
615  elseif (extendtoendofsimulation) then
616  ! -- Last time series value extends forever, so integrate the final value over all simulation time after the end of the series
617  nexttime = time1
618  nextval = currval
619  lprocess = .true.
620  end if
621  !
622  if (lprocess) then
623  ! -- determine lower and upper limits of time span of interest
624  ! within current interval
625  if (currtime > time0 .or. is_close(currtime, time0)) then
626  t0 = currtime
627  else
628  t0 = time0
629  end if
630  if (nexttime < time1 .or. is_close(nexttime, time1)) then
631  t1 = nexttime
632  else
633  t1 = time1
634  end if
635  ! -- find area of rectangle or trapezoid delimited by t0 and t1
636  t01 = t1 - t0
637  select case (this%iMethod)
638  case (stepwise)
639  ! -- compute area of a rectangle
640  value0 = currval
641  area = value0 * t01
642  case (linear, linearend)
643  ! -- compute area of a trapezoid
644  timediff = nexttime - currtime
645  ratio0 = (t0 - currtime) / timediff
646  ratio1 = (t1 - currtime) / timediff
647  valuediff = nextval - currval
648  value0 = currval + ratio0 * valuediff
649  value1 = currval + ratio1 * valuediff
650  if (this%iMethod == linear) then
651  area = 0.5d0 * t01 * (value0 + value1)
652  elseif (this%iMethod == linearend) then
653  area = dzero
654  value = value1
655  end if
656  end select
657  ! -- add area to integrated value
658  value = value + area
659  end if
660  end if
661  !
662  ! -- Are we done yet?
663  if (t1 > time1) then
664  ldone = .true.
665  elseif (is_close(t1, time1)) then
666  ldone = .true.
667  else
668  ! -- We are not done yet
669  if (.not. associated(currnode%nextNode)) then
670  ! -- Not done and no more data, so try to read the next record
671  if (.not. this%read_next_record()) then
672  write (errmsg, 10) trim(this%Name), time0, time1
673  call store_error(errmsg, terminate=.true.)
674  end if
675  elseif (associated(currnode%nextNode)) then
676  currnode => currnode%nextNode
677  end if
678  end if
679  end do
680  end if
681  !
682  get_integrated_value = value
683  if (this%autoDeallocate) then
684  if (associated(tslnodepreceding)) then
685  if (associated(tslnodepreceding%prevNode)) then
686  call this%list%DeallocateBackward(tslnodepreceding%prevNode)
687  end if
688  end if
689  end if
Here is the call graph for this function:

◆ get_latest_preceding_node()

subroutine timeseriesmodule::get_latest_preceding_node ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time,
type(listnodetype), intent(inout), pointer  tslNode 
)
private

Return pointer to ListNodeType object for the node representing the latest preceding time in the time series

Definition at line 729 of file TimeSeries.f90.

730  ! -- dummy
731  class(TimeSeriesType), intent(inout) :: this
732  real(DP), intent(in) :: time
733  type(ListNodeType), pointer, intent(inout) :: tslNode
734  ! -- local
735  real(DP) :: time0
736  type(ListNodeType), pointer :: currNode => null()
737  type(ListNodeType), pointer :: tsNode0 => null()
738  type(TimeSeriesRecordType), pointer :: tsr => null()
739  type(TimeSeriesRecordType), pointer :: tsrec0 => null()
740  class(*), pointer :: obj => null()
741  !
742  tslnode => null()
743  if (associated(this%list%firstNode)) then
744  currnode => this%list%firstNode
745  else
746  call store_error('probable programming error in &
747  &get_latest_preceding_node', &
748  terminate=.true.)
749  end if
750  !
751  ! -- If the next node is earlier than time of interest, advance along
752  ! linked list until the next node is later than time of interest.
753  do
754  if (associated(currnode)) then
755  if (associated(currnode%nextNode)) then
756  obj => currnode%nextNode%GetItem()
757  tsr => castastimeseriesrecordtype(obj)
758  if (tsr%tsrTime < time .or. is_close(tsr%tsrTime, time)) then
759  currnode => currnode%nextNode
760  else
761  exit
762  end if
763  else
764  ! -- read another record
765  if (.not. this%read_next_record()) exit
766  end if
767  else
768  exit
769  end if
770  end do
771  !
772  if (associated(currnode)) then
773  !
774  ! -- find earlier record
775  tsnode0 => currnode
776  obj => tsnode0%GetItem()
777  tsrec0 => castastimeseriesrecordtype(obj)
778  time0 = tsrec0%tsrTime
779  do while (time0 > time)
780  if (associated(tsnode0%prevNode)) then
781  tsnode0 => tsnode0%prevNode
782  obj => tsnode0%GetItem()
783  tsrec0 => castastimeseriesrecordtype(obj)
784  time0 = tsrec0%tsrTime
785  else
786  exit
787  end if
788  end do
789  end if
790  !
791  if (time0 < time .or. is_close(time0, time)) tslnode => tsnode0
Here is the call graph for this function:

◆ get_surrounding_nodes()

subroutine timeseriesmodule::get_surrounding_nodes ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time,
type(listnodetype), intent(inout), pointer  nodeEarlier,
type(listnodetype), intent(inout), pointer  nodeLater 
)
private

This subroutine is for working with time series already entirely stored in memory – it does not read data from a file.

Definition at line 362 of file TimeSeries.f90.

363  ! -- dummy
364  class(TimeSeriesType), intent(inout) :: this
365  real(DP), intent(in) :: time
366  type(ListNodeType), pointer, intent(inout) :: nodeEarlier
367  type(ListNodeType), pointer, intent(inout) :: nodeLater
368  ! -- local
369  real(DP) :: time0, time1
370  type(ListNodeType), pointer :: currNode => null()
371  type(ListNodeType), pointer :: tsNode0 => null()
372  type(ListNodeType), pointer :: tsNode1 => null()
373  type(TimeSeriesRecordType), pointer :: tsr => null(), tsrec0 => null()
374  type(TimeSeriesRecordType), pointer :: tsrec1 => null()
375  type(TimeSeriesRecordType), pointer :: tsrecEarlier
376  type(TimeSeriesRecordType), pointer :: tsrecLater
377  class(*), pointer :: obj => null()
378  !
379  tsrecearlier => null()
380  tsreclater => null()
381  nodeearlier => null()
382  nodelater => null()
383  !
384  if (associated(this%list%firstNode)) then
385  currnode => this%list%firstNode
386  end if
387  !
388  ! -- If the next node is earlier than time of interest, advance along
389  ! linked list until the next node is later than time of interest.
390  do
391  if (associated(currnode)) then
392  if (associated(currnode%nextNode)) then
393  obj => currnode%nextNode%GetItem()
394  tsr => castastimeseriesrecordtype(obj)
395  if (tsr%tsrTime < time .and. .not. is_close(tsr%tsrTime, time)) then
396  currnode => currnode%nextNode
397  else
398  exit
399  end if
400  else
401  exit
402  end if
403  else
404  exit
405  end if
406  end do
407  !
408  if (associated(currnode)) then
409  !
410  ! -- find earlier record
411  tsnode0 => currnode
412  obj => tsnode0%GetItem()
413  tsrec0 => castastimeseriesrecordtype(obj)
414  time0 = tsrec0%tsrTime
415  do while (time0 > time)
416  if (associated(tsnode0%prevNode)) then
417  tsnode0 => tsnode0%prevNode
418  obj => tsnode0%GetItem()
419  tsrec0 => castastimeseriesrecordtype(obj)
420  time0 = tsrec0%tsrTime
421  else
422  exit
423  end if
424  end do
425  !
426  ! -- find later record
427  tsnode1 => currnode
428  obj => tsnode1%GetItem()
429  tsrec1 => castastimeseriesrecordtype(obj)
430  time1 = tsrec1%tsrTime
431  do while (time1 < time .and. .not. is_close(time1, time))
432  if (associated(tsnode1%nextNode)) then
433  tsnode1 => tsnode1%nextNode
434  obj => tsnode1%GetItem()
435  tsrec1 => castastimeseriesrecordtype(obj)
436  time1 = tsrec1%tsrTime
437  else
438  exit
439  end if
440  end do
441  !
442  end if
443  !
444  if (time0 < time .or. is_close(time0, time)) then
445  tsrecearlier => tsrec0
446  nodeearlier => tsnode0
447  end if
448  if (time1 > time .or. is_close(time1, time)) then
449  tsreclater => tsrec1
450  nodelater => tsnode1
451  end if
Here is the call graph for this function:

◆ get_surrounding_records()

subroutine timeseriesmodule::get_surrounding_records ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time,
type(timeseriesrecordtype), intent(inout), pointer  tsrecEarlier,
type(timeseriesrecordtype), intent(inout), pointer  tsrecLater 
)
private

Definition at line 270 of file TimeSeries.f90.

271  ! -- dummy
272  class(TimeSeriesType), intent(inout) :: this
273  real(DP), intent(in) :: time
274  type(TimeSeriesRecordType), pointer, intent(inout) :: tsrecEarlier
275  type(TimeSeriesRecordType), pointer, intent(inout) :: tsrecLater
276  ! -- local
277  real(DP) :: time0, time1
278  type(ListNodeType), pointer :: currNode => null()
279  type(ListNodeType), pointer :: tsNode0 => null()
280  type(ListNodeType), pointer :: tsNode1 => null()
281  type(TimeSeriesRecordType), pointer :: tsr => null(), tsrec0 => null()
282  type(TimeSeriesRecordType), pointer :: tsrec1 => null()
283  class(*), pointer :: obj => null()
284  !
285  tsrecearlier => null()
286  tsreclater => null()
287  !
288  if (associated(this%list%firstNode)) then
289  currnode => this%list%firstNode
290  end if
291  !
292  ! -- If the next node is earlier than time of interest, advance along
293  ! linked list until the next node is later than time of interest.
294  do
295  if (associated(currnode)) then
296  if (associated(currnode%nextNode)) then
297  obj => currnode%nextNode%GetItem()
298  tsr => castastimeseriesrecordtype(obj)
299  if (tsr%tsrTime < time .and. .not. is_close(tsr%tsrTime, time)) then
300  currnode => currnode%nextNode
301  else
302  exit
303  end if
304  else
305  ! -- read another record
306  if (.not. this%read_next_record()) exit
307  end if
308  else
309  exit
310  end if
311  end do
312  !
313  if (associated(currnode)) then
314  !
315  ! -- find earlier record
316  tsnode0 => currnode
317  obj => tsnode0%GetItem()
318  tsrec0 => castastimeseriesrecordtype(obj)
319  time0 = tsrec0%tsrTime
320  do while (time0 > time)
321  if (associated(tsnode0%prevNode)) then
322  tsnode0 => tsnode0%prevNode
323  obj => tsnode0%GetItem()
324  tsrec0 => castastimeseriesrecordtype(obj)
325  time0 = tsrec0%tsrTime
326  else
327  exit
328  end if
329  end do
330  !
331  ! -- find later record
332  tsnode1 => currnode
333  obj => tsnode1%GetItem()
334  tsrec1 => castastimeseriesrecordtype(obj)
335  time1 = tsrec1%tsrTime
336  do while (time1 < time .and. .not. is_close(time1, time))
337  if (associated(tsnode1%nextNode)) then
338  tsnode1 => tsnode1%nextNode
339  obj => tsnode1%GetItem()
340  tsrec1 => castastimeseriesrecordtype(obj)
341  time1 = tsrec1%tsrTime
342  else
343  ! -- get next record
344  if (.not. this%read_next_record()) then
345  ! -- end of file reached, so exit loop
346  exit
347  end if
348  end if
349  end do
350  !
351  end if
352  !
353  if (time0 < time .or. is_close(time0, time)) tsrecearlier => tsrec0
354  if (time1 > time .or. is_close(time1, time)) tsreclater => tsrec1
Here is the call graph for this function:

◆ get_value_at_time()

real(dp) function timeseriesmodule::get_value_at_time ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time,
logical, intent(in)  extendToEndOfSimulation 
)
private

Return a value for a specified time, same units as time-series values

Definition at line 478 of file TimeSeries.f90.

479  ! -- return
480  real(DP) :: get_value_at_time
481  ! -- dummy
482  class(TimeSeriesType), intent(inout) :: this
483  real(DP), intent(in) :: time ! time of interest
484  logical, intent(in) :: extendToEndOfSimulation
485  ! -- local
486  integer(I4B) :: ierr
487  real(DP) :: ratio, time0, time1, timediff, timediffi, val0, val1, &
488  valdiff
489  type(TimeSeriesRecordType), pointer :: tsrEarlier => null()
490  type(TimeSeriesRecordType), pointer :: tsrLater => null()
491  ! -- formats
492 10 format('Error getting value at time ', g10.3, ' for time series "', a, '"')
493  !
494  ierr = 0
495  call this%get_surrounding_records(time, tsrearlier, tsrlater)
496  if (associated(tsrearlier)) then
497  if (associated(tsrlater)) then
498  ! -- values are available for both earlier and later times
499  if (this%iMethod == stepwise) then
500  get_value_at_time = tsrearlier%tsrValue
501  elseif (this%iMethod == linear .or. this%iMethod == linearend) then
502  ! -- For get_value_at_time, result is the same for either
503  ! linear method.
504  ! -- Perform linear interpolation.
505  time0 = tsrearlier%tsrTime
506  time1 = tsrlater%tsrtime
507  timediff = time1 - time0
508  timediffi = time - time0
509  if (timediff > 0) then
510  ratio = timediffi / timediff
511  else
512  ! -- should not happen if TS does not contain duplicate times
513  ratio = 0.5d0
514  end if
515  val0 = tsrearlier%tsrValue
516  val1 = tsrlater%tsrValue
517  valdiff = val1 - val0
518  get_value_at_time = val0 + (ratio * valdiff)
519  else
520  ierr = 1
521  end if
522  else
523  if (extendtoendofsimulation .or. is_close(tsrearlier%tsrTime, time)) then
524  get_value_at_time = tsrearlier%tsrValue
525  else
526  ! -- Only earlier time is available, and it is not time of interest;
527  ! however, if method is STEPWISE, use value for earlier time.
528  if (this%iMethod == stepwise) then
529  get_value_at_time = tsrearlier%tsrValue
530  else
531  ierr = 1
532  end if
533  end if
534  end if
535  else
536  if (associated(tsrlater)) then
537  if (is_close(tsrlater%tsrTime, time)) then
538  get_value_at_time = tsrlater%tsrValue
539  else
540  ! -- only later time is available, and it is not time of interest
541  ierr = 1
542  end if
543  else
544  ! -- Neither earlier nor later time is available.
545  ! This should never happen!
546  ierr = 1
547  end if
548  end if
549  !
550  if (ierr > 0) then
551  write (errmsg, 10) time, trim(this%Name)
552  call store_error(errmsg, terminate=.true.)
553  end if
Here is the call graph for this function:

◆ getcurrenttimeseriesrecord()

type(timeseriesrecordtype) function, pointer timeseriesmodule::getcurrenttimeseriesrecord ( class(timeseriestype this)
private

Definition at line 821 of file TimeSeries.f90.

822  ! -- dummy
823  class(TimeSeriesType) :: this
824  ! -- result
825  type(TimeSeriesRecordType), pointer :: res
826  ! -- local
827  class(*), pointer :: obj => null()
828  !
829  obj => null()
830  res => null()
831  obj => this%list%GetItem()
832  if (associated(obj)) then
833  res => castastimeseriesrecordtype(obj)
834  end if
Here is the call graph for this function:

◆ getnexttimeseriesrecord()

type(timeseriesrecordtype) function, pointer timeseriesmodule::getnexttimeseriesrecord ( class(timeseriestype this)
private

Definition at line 857 of file TimeSeries.f90.

858  ! -- dummy
859  class(TimeSeriesType) :: this
860  ! -- result
861  type(TimeSeriesRecordType), pointer :: res
862  ! -- local
863  class(*), pointer :: obj => null()
864  !
865  obj => null()
866  res => null()
867  obj => this%list%GetNextItem()
868  if (associated(obj)) then
869  res => castastimeseriesrecordtype(obj)
870  end if
Here is the call graph for this function:

◆ getprevioustimeseriesrecord()

type(timeseriesrecordtype) function, pointer timeseriesmodule::getprevioustimeseriesrecord ( class(timeseriestype this)
private

Definition at line 839 of file TimeSeries.f90.

840  ! -- dummy
841  class(TimeSeriesType) :: this
842  ! -- result
843  type(TimeSeriesRecordType), pointer :: res
844  ! -- local
845  class(*), pointer :: obj => null()
846  !
847  obj => null()
848  res => null()
849  obj => this%list%GetPreviousItem()
850  if (associated(obj)) then
851  res => castastimeseriesrecordtype(obj)
852  end if
Here is the call graph for this function:

◆ gettimeseries()

type(timeseriestype) function, pointer timeseriesmodule::gettimeseries ( class(timeseriesfiletype this,
integer(i4b), intent(in)  indx 
)
private

Definition at line 1045 of file TimeSeries.f90.

1046  ! -- dummy
1047  class(TimeSeriesFileType) :: this
1048  integer(I4B), intent(in) :: indx
1049  ! -- return
1050  type(TimeSeriesType), pointer :: res
1051  !
1052  res => null()
1053  if (indx > 0 .and. indx <= this%nTimeSeries) then
1054  res => this%timeSeries(indx)
1055  end if

◆ gettimeseriesfilefromlist()

type(timeseriesfiletype) function, pointer, public timeseriesmodule::gettimeseriesfilefromlist ( type(listtype), intent(inout)  list,
integer(i4b), intent(in)  idx 
)

Definition at line 152 of file TimeSeries.f90.

153  ! -- dummy
154  type(ListType), intent(inout) :: list
155  integer(I4B), intent(in) :: idx
156  ! -- return
157  type(TimeSeriesFileType), pointer :: res
158  ! -- local
159  class(*), pointer :: obj => null()
160  !
161  obj => list%GetItem(idx)
162  res => castastimeseriesfiletype(obj)
163  !
164  if (.not. associated(res)) then
165  res => castastimeseriesfileclass(obj)
166  end if
Here is the call graph for this function:
Here is the caller graph for this function:

◆ gettimeseriesrecord()

type(timeseriesrecordtype) function, pointer timeseriesmodule::gettimeseriesrecord ( class(timeseriestype this,
double precision, intent(in)  time,
double precision, intent(in)  epsi 
)
private

Definition at line 875 of file TimeSeries.f90.

876  ! -- dummy
877  class(TimeSeriesType) :: this
878  double precision, intent(in) :: time
879  double precision, intent(in) :: epsi
880  ! -- result
881  type(TimeSeriesRecordType), pointer :: res
882  ! -- local
883  type(TimeSeriesRecordType), pointer :: tsr
884  !
885  call this%list%Reset()
886  res => null()
887  do
888  tsr => this%GetNextTimeSeriesRecord()
889  if (associated(tsr)) then
890  if (is_close(tsr%tsrTime, time)) then
891  res => tsr
892  exit
893  end if
894  if (tsr%tsrTime > time) exit
895  else
896  exit
897  end if
898  end do
Here is the call graph for this function:

◆ getvalue()

real(dp) function timeseriesmodule::getvalue ( class(timeseriestype), intent(inout)  this,
real(dp), intent(in)  time0,
real(dp), intent(in)  time1,
logical, intent(in), optional  extendToEndOfSimulation 
)
private

If iMethod is STEPWISE or LINEAR: Return a time-weighted average value for a specified time span. If iMethod is LINEAREND: Return value at time1. Time0 argument is ignored. Units: (ts-value-unit)

Definition at line 209 of file TimeSeries.f90.

210  ! -- return
211  real(DP) :: GetValue
212  ! -- dummy
213  class(TimeSeriesType), intent(inout) :: this
214  real(DP), intent(in) :: time0
215  real(DP), intent(in) :: time1
216  logical, intent(in), optional :: extendToEndOfSimulation
217  ! -- local
218  logical :: extend
219  !
220  if (present(extendtoendofsimulation)) then
221  extend = extendtoendofsimulation
222  else
223  extend = .false.
224  end if
225  !
226  select case (this%iMethod)
227  case (stepwise, linear)
228  getvalue = this%get_average_value(time0, time1, extend)
229  case (linearend)
230  getvalue = this%get_value_at_time(time1, extend)
231  end select

◆ initialize_time_series()

subroutine timeseriesmodule::initialize_time_series ( class(timeseriestype), intent(inout)  this,
class(timeseriesfiletype), target  tsfile,
character(len=*), intent(in)  name,
logical, intent(in), optional  autoDeallocate 
)
private

Open time-series file and read options and first time-series record.

Definition at line 238 of file TimeSeries.f90.

239  ! -- dummy
240  class(TimeSeriesType), intent(inout) :: this
241  class(TimeSeriesFileType), target :: tsfile
242  character(len=*), intent(in) :: name
243  logical, intent(in), optional :: autoDeallocate
244  ! -- local
245  character(len=LENTIMESERIESNAME) :: tsNameTemp
246  !
247  ! -- Assign the time-series tsfile, name, and autoDeallocate
248  this%tsfile => tsfile
249  ! Store time-series name as all caps
250  tsnametemp = name
251  call upcase(tsnametemp)
252  this%Name = tsnametemp
253  !
254  this%iMethod = undefined
255  !
256  if (present(autodeallocate)) this%autoDeallocate = autodeallocate
257  !
258  ! -- allocate the list
259  allocate (this%list)
260  !
261  ! -- ensure that NAME has been specified
262  if (this%Name == '') then
263  errmsg = 'Name not specified for time series.'
264  call store_error(errmsg, terminate=.true.)
265  end if
Here is the call graph for this function:

◆ initializetsfile()

subroutine timeseriesmodule::initializetsfile ( class(timeseriesfiletype), intent(inout), target  this,
character(len=*), intent(in)  filename,
integer(i4b), intent(in)  iout,
logical, intent(in), optional  autoDeallocate 
)
private

Definition at line 1061 of file TimeSeries.f90.

1062  ! -- dummy
1063  class(TimeSeriesFileType), target, intent(inout) :: this
1064  character(len=*), intent(in) :: filename
1065  integer(I4B), intent(in) :: iout
1066  logical, optional, intent(in) :: autoDeallocate
1067  ! -- local
1068  integer(I4B) :: iMethod, istatus, j, nwords
1069  integer(I4B) :: ierr, inunit
1070  logical :: autoDeallocateLocal = .true.
1071  logical :: continueread, found, endOfBlock
1072  logical :: methodWasSet
1073  real(DP) :: sfaclocal
1074  character(len=40) :: keyword, keyvalue
1075  character(len=:), allocatable :: line
1076  character(len=LENTIMESERIESNAME), allocatable, dimension(:) :: words
1077  !
1078  ! -- Initialize some variables
1079  if (present(autodeallocate)) autodeallocatelocal = autodeallocate
1080  imethod = undefined
1081  methodwasset = .false.
1082  !
1083  ! -- Assign members
1084  this%iout = iout
1085  this%datafile = filename
1086  !
1087  ! -- Open the time-series tsfile input file
1088  this%inunit = getunit()
1089  inunit = this%inunit
1090  call openfile(inunit, 0, filename, 'TS6')
1091  !
1092  ! -- Initialize block parser
1093  call this%parser%Initialize(this%inunit, this%iout)
1094  !
1095  ! -- Read the ATTRIBUTES block and count time series
1096  continueread = .false.
1097  ierr = 0
1098  !
1099  ! -- get BEGIN line of ATTRIBUTES block
1100  call this%parser%GetBlock('ATTRIBUTES', found, ierr, &
1101  supportopenclose=.true.)
1102  if (ierr /= 0) then
1103  ! end of file
1104  errmsg = 'End-of-file encountered while searching for'// &
1105  ' ATTRIBUTES in time-series '// &
1106  'input file "'//trim(this%datafile)//'"'
1107  call store_error(errmsg)
1108  call this%parser%StoreErrorUnit()
1109  elseif (.not. found) then
1110  errmsg = 'ATTRIBUTES block not found in time-series '// &
1111  'tsfile input file "'//trim(this%datafile)//'"'
1112  call store_error(errmsg)
1113  call this%parser%StoreErrorUnit()
1114  end if
1115  !
1116  ! -- parse ATTRIBUTES entries
1117  do
1118  ! -- read a line from input
1119  call this%parser%GetNextLine(endofblock)
1120  if (endofblock) exit
1121  !
1122  ! -- get the keyword
1123  call this%parser%GetStringCaps(keyword)
1124  !
1125  ! support either NAME or NAMES as equivalent keywords
1126  if (keyword == 'NAMES') keyword = 'NAME'
1127  !
1128  if (keyword /= 'NAME' .and. keyword /= 'METHODS' .and. &
1129  keyword /= 'SFACS') then
1130  ! -- get the word following the keyword (the key value)
1131  call this%parser%GetStringCaps(keyvalue)
1132  end if
1133  !
1134  select case (keyword)
1135  case ('NAME')
1136 ! line = line(istart:linelen)
1137  call this%parser%GetRemainingLine(line)
1138  call parseline(line, nwords, words, this%parser%iuactive)
1139  this%nTimeSeries = nwords
1140  ! -- Allocate the timeSeries array and initialize each
1141  ! time series.
1142  allocate (this%timeSeries(this%nTimeSeries))
1143  do j = 1, this%nTimeSeries
1144  call this%timeSeries(j)%initialize_time_series(this, words(j), &
1145  autodeallocatelocal)
1146  end do
1147  case ('METHOD')
1148  methodwasset = .true.
1149  if (this%nTimeSeries == 0) then
1150  errmsg = 'Error: NAME attribute not provided before METHOD in file: ' &
1151  //trim(filename)
1152  call store_error(errmsg)
1153  call this%parser%StoreErrorUnit()
1154  end if
1155  select case (keyvalue)
1156  case ('STEPWISE')
1157  imethod = stepwise
1158  case ('LINEAR')
1159  imethod = linear
1160  case ('LINEAREND')
1161  imethod = linearend
1162  case default
1163  errmsg = 'Unknown interpolation method: "'//trim(keyvalue)//'"'
1164  call store_error(errmsg)
1165  end select
1166  do j = 1, this%nTimeSeries
1167  this%timeSeries(j)%iMethod = imethod
1168  end do
1169  case ('METHODS')
1170  methodwasset = .true.
1171  if (this%nTimeSeries == 0) then
1172  errmsg = 'Error: NAME attribute not provided before METHODS in file: ' &
1173  //trim(filename)
1174  call store_error(errmsg)
1175  call this%parser%StoreErrorUnit()
1176  end if
1177  call this%parser%GetRemainingLine(line)
1178  call parseline(line, nwords, words, this%parser%iuactive)
1179  if (nwords < this%nTimeSeries) then
1180  errmsg = 'METHODS attribute does not list a method for'// &
1181  ' all time series.'
1182  call store_error(errmsg)
1183  call this%parser%StoreErrorUnit()
1184  end if
1185  do j = 1, this%nTimeSeries
1186  call upcase(words(j))
1187  select case (words(j))
1188  case ('STEPWISE')
1189  imethod = stepwise
1190  case ('LINEAR')
1191  imethod = linear
1192  case ('LINEAREND')
1193  imethod = linearend
1194  case default
1195  errmsg = 'Unknown interpolation method: "'//trim(words(j))//'"'
1196  call store_error(errmsg)
1197  end select
1198  this%timeSeries(j)%iMethod = imethod
1199  end do
1200  case ('SFAC')
1201  if (this%nTimeSeries == 0) then
1202  errmsg = 'NAME attribute not provided before SFAC in file: ' &
1203  //trim(filename)
1204  call store_error(errmsg)
1205  call this%parser%StoreErrorUnit()
1206  end if
1207  read (keyvalue, *, iostat=istatus) sfaclocal
1208  if (istatus /= 0) then
1209  errmsg = 'Error reading numeric value from: "'//trim(keyvalue)//'"'
1210  call store_error(errmsg)
1211  end if
1212  do j = 1, this%nTimeSeries
1213  this%timeSeries(j)%sfac = sfaclocal
1214  end do
1215  case ('SFACS')
1216  if (this%nTimeSeries == 0) then
1217  errmsg = 'NAME attribute not provided before SFACS in file: ' &
1218  //trim(filename)
1219  call store_error(errmsg)
1220  call this%parser%StoreErrorUnit()
1221  end if
1222  do j = 1, this%nTimeSeries
1223  sfaclocal = this%parser%GetDouble()
1224  this%timeSeries(j)%sfac = sfaclocal
1225  end do
1226  case ('AUTODEALLOCATE')
1227  do j = 1, this%nTimeSeries
1228  this%timeSeries(j)%autoDeallocate = (keyvalue == 'TRUE')
1229  end do
1230  case default
1231  errmsg = 'Unknown option found in ATTRIBUTES block: "'// &
1232  trim(keyword)//'"'
1233  call store_error(errmsg)
1234  call this%parser%StoreErrorUnit()
1235  end select
1236  end do
1237  !
1238  ! -- Get TIMESERIES block
1239  call this%parser%GetBlock('TIMESERIES', found, ierr, &
1240  supportopenclose=.true.)
1241  !
1242  ! -- Read the first line of time-series data
1243  if (.not. this%read_tsfile_line()) then
1244  errmsg = 'Error: No time-series data contained in file: '// &
1245  trim(this%datafile)
1246  call store_error(errmsg)
1247  end if
1248  !
1249  ! -- Ensure method was set
1250  if (.not. methodwasset) then
1251  errmsg = 'Interpolation method was not set. METHOD or METHODS &
1252  &must be specified in the ATTRIBUTES block for this time series file.'
1253  call store_error(errmsg)
1254  end if
1255  !
1256  ! -- Clean up and return
1257  if (allocated(words)) deallocate (words)
1258  !
1259  if (count_errors() > 0) then
1260  call this%parser%StoreErrorUnit()
1261  end if
Here is the call graph for this function:

◆ inserttsr()

subroutine timeseriesmodule::inserttsr ( class(timeseriestype), intent(inout)  this,
type(timeseriesrecordtype), intent(inout), pointer  tsr 
)
private

Definition at line 912 of file TimeSeries.f90.

913  ! -- dummy
914  class(TimeSeriesType), intent(inout) :: this
915  type(TimeSeriesRecordType), pointer, intent(inout) :: tsr
916  ! -- local
917  double precision :: badtime, time, time0, time1
918  type(TimeSeriesRecordType), pointer :: tsrEarlier, tsrLater
919  type(ListNodeType), pointer :: nodeEarlier, nodeLater
920  class(*), pointer :: obj => null()
921  !
922  badtime = -9.0d30
923  time0 = badtime
924  time1 = badtime
925  time = tsr%tsrTime
926  call this%get_surrounding_nodes(time, nodeearlier, nodelater)
927  !
928  if (associated(nodeearlier)) then
929  obj => nodeearlier%GetItem()
930  tsrearlier => castastimeseriesrecordtype(obj)
931  if (associated(tsrearlier)) then
932  time0 = tsrearlier%tsrTime
933  end if
934  end if
935  !
936  if (associated(nodelater)) then
937  obj => nodelater%GetItem()
938  tsrlater => castastimeseriesrecordtype(obj)
939  if (associated(tsrlater)) then
940  time1 = tsrlater%tsrTime
941  end if
942  end if
943  !
944  if (time0 > badtime) then
945  ! Time0 is valid
946  if (time1 > badtime) then
947  ! Both time0 and time1 are valid
948  if (time > time0 .and. time < time1) then
949  ! Insert record between two list nodes
950  obj => tsr
951  call this%list%InsertBefore(obj, nodelater)
952  else
953  ! No need to insert a time series record, but if existing record
954  ! for time of interest has NODATA as tsrValue, replace tsrValue
955  if (time == time0 .and. tsrearlier%tsrValue == dnodata .and. &
956  tsr%tsrValue /= dnodata) then
957  tsrearlier%tsrValue = tsr%tsrValue
958  elseif (time == time1 .and. tsrlater%tsrValue == dnodata .and. &
959  tsr%tsrValue /= dnodata) then
960  tsrlater%tsrValue = tsr%tsrValue
961  end if
962  end if
963  else
964  ! Time0 is valid and time1 is invalid. Just add tsr to the list.
965  call this%AddTimeSeriesRecord(tsr)
966  end if
967  else
968  ! Time0 is invalid, so time1 must be for first node in list
969  if (time1 > badtime) then
970  ! Time 1 is valid
971  if (time < time1) then
972  ! Insert tsr at beginning of list
973  obj => tsr
974  call this%list%InsertBefore(obj, nodelater)
975  elseif (time == time1) then
976  ! No need to insert a time series record, but if existing record
977  ! for time of interest has NODATA as tsrValue, replace tsrValue
978  if (tsrlater%tsrValue == dnodata .and. tsr%tsrValue /= dnodata) then
979  tsrlater%tsrValue = tsr%tsrValue
980  end if
981  end if
982  else
983  ! Both time0 and time1 are invalid. Just add tsr to the list.
984  call this%AddTimeSeriesRecord(tsr)
985  end if
986  end if
Here is the call graph for this function:

◆ read_next_record()

logical function timeseriesmodule::read_next_record ( class(timeseriestype), intent(inout)  this)
private

Read next time-series record from input file

Definition at line 458 of file TimeSeries.f90.

459  ! -- dummy
460  class(TimeSeriesType), intent(inout) :: this
461  !
462  ! -- If we have already encountered the end of the TIMESERIES block, do not try to read any further
463  if (this%tsfile%finishedReading) then
464  read_next_record = .false.
465  return
466  end if
467  !
468  read_next_record = this%tsfile%read_tsfile_line()
469  if (.not. read_next_record) then
470  this%tsfile%finishedReading = .true.
471  end if

◆ read_tsfile_line()

logical function timeseriesmodule::read_tsfile_line ( class(timeseriesfiletype), intent(inout)  this)
private

Definition at line 1265 of file TimeSeries.f90.

1266  ! -- dummy
1267  class(TimeSeriesFileType), intent(inout) :: this
1268  ! -- local
1269  real(DP) :: tsrTime, tsrValue
1270  integer(I4B) :: i
1271  logical :: endOfBlock
1272  type(TimeSeriesRecordType), pointer :: tsRecord => null()
1273  !
1274  read_tsfile_line = .false.
1275  !
1276  ! -- Get an arbitrary length, non-comment, non-blank line
1277  ! from the input file.
1278  call this%parser%GetNextLine(endofblock)
1279  !
1280  ! -- Check if we've reached the end of the TIMESERIES block
1281  if (endofblock) then
1282  return
1283  end if
1284  !
1285  ! -- Get the time
1286  tsrtime = this%parser%GetDouble()
1287  !
1288  ! -- Construct a new record and append a new node to each time series
1289  tsloop: do i = 1, this%nTimeSeries
1290  tsrvalue = this%parser%GetDouble()
1291  if (tsrvalue == dnodata) cycle tsloop
1292  ! -- multiply value by sfac
1293  tsrvalue = tsrvalue * this%timeSeries(i)%sfac
1294  call constructtimeseriesrecord(tsrecord, tsrtime, tsrvalue)
1295  call addtimeseriesrecordtolist(this%timeSeries(i)%list, tsrecord)
1296  end do tsloop
1297  read_tsfile_line = .true.
Here is the call graph for this function:

◆ reset()

subroutine timeseriesmodule::reset ( class(timeseriestype this)
private

Definition at line 903 of file TimeSeries.f90.

904  ! -- dummy
905  class(TimeSeriesType) :: this
906  !
907  call this%list%Reset()

◆ sametimeseries()

logical function, public timeseriesmodule::sametimeseries ( type(timeseriestype), intent(in)  ts1,
type(timeseriestype), intent(in)  ts2 
)

Definition at line 171 of file TimeSeries.f90.

172  ! -- dummy
173  type(TimeSeriesType), intent(in) :: ts1
174  type(TimeSeriesType), intent(in) :: ts2
175  ! -- return
176  logical :: same
177  ! -- local
178  integer :: i, n1, n2
179  type(TimeSeriesRecordType), pointer :: tsr1, tsr2
180  !
181  same = .false.
182  n1 = ts1%list%Count()
183  n2 = ts2%list%Count()
184  if (n1 /= n2) return
185  !
186  call ts1%Reset()
187  call ts2%Reset()
188  !
189  do i = 1, n1
190  tsr1 => ts1%GetNextTimeSeriesRecord()
191  tsr2 => ts2%GetNextTimeSeriesRecord()
192  if (tsr1%tsrTime /= tsr2%tsrTime) return
193  if (tsr1%tsrValue /= tsr2%tsrValue) return
194  end do
195  !
196  same = .true.

◆ ts_da()

subroutine timeseriesmodule::ts_da ( class(timeseriestype), intent(inout)  this)
private

Definition at line 796 of file TimeSeries.f90.

797  ! -- dummy
798  class(TimeSeriesType), intent(inout) :: this
799  !
800  if (associated(this%list)) then
801  call this%list%Clear(.true.)
802  deallocate (this%list)
803  end if

◆ tsf_da()

subroutine timeseriesmodule::tsf_da ( class(timeseriesfiletype), intent(inout)  this)
private

Definition at line 1302 of file TimeSeries.f90.

1303  ! -- dummy
1304  class(TimeSeriesFileType), intent(inout) :: this
1305  ! -- local
1306  integer :: i, n
1307  type(TimeSeriesType), pointer :: ts => null()
1308  !
1309  n = this%Count()
1310  do i = 1, n
1311  ts => this%GetTimeSeries(i)
1312  if (associated(ts)) then
1313  call ts%da()
1314 ! deallocate(ts)
1315  end if
1316  end do
1317  !
1318  deallocate (this%timeSeries)
1319  deallocate (this%parser)