16 integer(I4B),
pointer :: npoints
17 real(dp),
pointer,
dimension(:),
contiguous :: station => null()
18 real(dp),
pointer,
dimension(:),
contiguous :: height => null()
19 real(dp),
pointer,
dimension(:),
contiguous :: roughfraction => null()
20 logical(LGP),
pointer,
dimension(:),
contiguous :: valid => null()
24 integer(I4B),
pointer :: iout => null()
25 integer(I4B),
pointer :: iprpak => null()
26 integer(I4B),
pointer :: nreaches => null()
27 integer(I4B),
pointer :: invalid => null()
28 character(len=LINELENGTH),
dimension(:),
allocatable :: filenames
29 integer(I4B),
pointer,
dimension(:),
contiguous :: npoints => null()
31 pointer,
dimension(:),
contiguous :: cross_sections => null()
60 integer(I4B),
pointer,
intent(in) :: iout
61 integer(I4B),
pointer,
intent(in) :: iprpak
62 integer(I4B),
pointer,
intent(in) :: nreaches
65 if (
associated(this))
then
77 this%nreaches => nreaches
86 subroutine initialize(this, ncrossptstot, ncrosspts, iacross, &
87 station, height, roughfraction)
90 integer(I4B),
intent(in) :: ncrossptstot
91 integer(I4B),
dimension(this%nreaches),
intent(in) :: ncrosspts
92 integer(I4B),
dimension(this%nreaches + 1),
intent(in) :: iacross
93 real(DP),
dimension(ncrossptstot),
intent(in) :: station
94 real(DP),
dimension(ncrossptstot),
intent(in) :: height
95 real(DP),
dimension(ncrossptstot),
intent(in) :: roughfraction
99 integer(I4B) :: npoints
105 allocate (this%invalid)
111 allocate (this%filenames(this%nreaches))
112 allocate (this%npoints(this%nreaches))
113 allocate (this%cross_sections(this%nreaches))
114 do n = 1, this%nreaches
115 npoints = ncrosspts(n)
116 allocate (this%cross_sections(n)%npoints)
117 allocate (this%cross_sections(n)%station(npoints))
118 allocate (this%cross_sections(n)%height(npoints))
119 allocate (this%cross_sections(n)%roughfraction(npoints))
120 allocate (this%cross_sections(n)%valid(npoints))
124 do n = 1, this%nreaches
125 this%filenames(n) =
'NONE'
126 this%cross_sections(n)%npoints = ncrosspts(n)
127 this%npoints(n) = ncrosspts(n)
129 i1 = iacross(n + 1) - 1
132 this%cross_sections(n)%station(ipos) = station(i)
133 this%cross_sections(n)%height(ipos) = height(i)
134 this%cross_sections(n)%roughfraction(ipos) = roughfraction(i)
135 this%cross_sections(n)%valid(ipos) = .true.
153 integer(I4B),
intent(in) :: irch
154 real(DP),
intent(in) :: width
155 character(len=*),
intent(in) :: filename
157 character(len=LINELENGTH) :: tag
158 character(len=LINELENGTH) :: keyword
161 logical :: endOfBlock
175 write (tag,
"('Reach',1x,i0,1x,'(',a, ')')") &
176 irch, trim(adjustl(filename))
180 call openfile(iu, this%iout, filename,
'SFR TABLE')
181 call parser%Initialize(iu, this%iout)
184 call parser%GetBlock(
'DIMENSIONS', isfound, ierr, supportopenclose=.true.)
190 if (this%iprpak /= 0)
then
191 write (this%iout,
'(/1x,a)') &
192 'PROCESSING '//trim(adjustl(tag))//
' DIMENSIONS'
195 call parser%GetNextLine(endofblock)
197 call parser%GetStringCaps(keyword)
198 select case (keyword)
200 n = parser%GetInteger()
202 write (
errmsg,
'(a)')
'Table NROW must be > 0'
206 j = parser%GetInteger()
209 write (
errmsg,
'(a,1x,i0)')
'Table NCOL must be >= ', jmin
214 'UNKNOWN '//trim(adjustl(tag))//
' DIMENSIONS keyword: ', &
219 if (this%iprpak /= 0)
then
220 write (this%iout,
'(1x,a)') &
221 'END OF '//trim(adjustl(tag))//
' DIMENSIONS'
224 call store_error(
'Required DIMENSIONS block not found.')
230 'NROW not specified in the table DIMENSIONS block'
235 'NCOL not specified in the table DIMENSIONS block'
244 this%filenames(irch) = filename
245 this%npoints(irch) = n
248 deallocate (this%cross_sections(irch)%npoints)
249 deallocate (this%cross_sections(irch)%station)
250 deallocate (this%cross_sections(irch)%height)
251 deallocate (this%cross_sections(irch)%roughfraction)
252 deallocate (this%cross_sections(irch)%valid)
255 allocate (this%cross_sections(irch)%npoints)
256 allocate (this%cross_sections(irch)%station(n))
257 allocate (this%cross_sections(irch)%height(n))
258 allocate (this%cross_sections(irch)%roughfraction(n))
259 allocate (this%cross_sections(irch)%valid(n))
262 this%cross_sections(irch)%npoints = n
265 call parser%GetBlock(
'TABLE', isfound, ierr, supportopenclose=.true.)
271 if (this%iprpak /= 0)
then
272 write (this%iout,
'(/1x,a)') &
273 'PROCESSING '//trim(adjustl(tag))//
' TABLE'
277 call parser%GetNextLine(endofblock)
280 if (ipos > this%npoints(irch))
then
283 this%cross_sections(irch)%station(ipos) = parser%GetDouble() * width
284 this%cross_sections(irch)%height(ipos) = parser%GetDouble()
286 this%cross_sections(irch)%roughfraction(ipos) = parser%GetDouble()
288 this%cross_sections(irch)%roughfraction(ipos) = done
290 this%cross_sections(irch)%valid(ipos) = .true.
293 if (this%iprpak /= 0)
then
294 write (this%iout,
'(1x,a)') &
295 'END OF '//trim(adjustl(tag))//
' TABLE'
298 call store_error(
'Required TABLE block not found.')
302 if (ipos /= this%npoints(irch))
then
303 write (
errmsg,
'(a,1x,i0,1x,a,1x,i0,1x,a)') &
304 'NROW set to', this%npoints(irch),
'but', ipos,
'rows were read'
315 call this%validate(irch)
331 integer(I4B),
intent(in) :: irch
333 logical(LGP) :: station_error
334 logical(LGP) :: height_error
335 logical(LGP) :: height_zero_error
336 logical(LGP) :: roughness_error
337 character(len=LINELENGTH) :: filename
344 real(DP) :: roughfraction
349 real(DP),
dimension(:),
allocatable :: heights
350 real(DP),
dimension(:),
allocatable :: unique_heights
351 real(DP),
dimension(3) :: factor
354 station_error = .false.
355 height_error = .false.
356 height_zero_error = .true.
357 roughness_error = .false.
358 npts = this%npoints(irch)
362 station = this%cross_sections(irch)%station(n)
363 if (station < dzero)
then
364 station_error = .true.
366 height = this%cross_sections(irch)%height(n)
367 if (height < dzero)
then
368 height_error = .true.
369 else if (height == dzero)
then
370 height_zero_error = .false.
372 roughfraction = this%cross_sections(irch)%roughfraction(n)
373 if (roughfraction <= dzero)
then
374 roughness_error = .true.
376 if (station_error .and. height_error .and. &
377 roughness_error)
then
383 if (station_error .or. height_error .or. &
384 height_zero_error .or. roughness_error)
then
385 filename = this%filenames(irch)
386 if (station_error)
then
387 write (
errmsg,
'(3a,1x,i0,1x,a)') &
388 "All xfraction data in '", trim(adjustl(filename)), &
389 "' for reach", irch,
'must be greater than or equal to zero.'
392 if (height_error)
then
393 write (
errmsg,
'(3a,1x,i0,1x,a)') &
394 "All height data in '", trim(adjustl(filename)), &
395 "' for reach", irch,
'must be greater than or equal to zero.'
398 if (height_zero_error)
then
399 write (
errmsg,
'(3a,1x,i0,1x,a)') &
400 "At least one height data value in '", trim(adjustl(filename)), &
401 "' for reach", irch,
'must be equal to zero.'
404 if (roughness_error)
then
405 write (
errmsg,
'(3a,1x,i0,1x,a)') &
406 "All manfraction data in '", trim(adjustl(filename)), &
407 "' for reach", irch,
'must be greater than zero.'
413 allocate (heights(npts))
415 heights(n) = this%cross_sections(irch)%height(n)
423 do n = 1,
size(unique_heights)
424 if (unique_heights(n) <= dzero) cycle
427 height = unique_heights(n) + real(i, dp) *
dem6
429 this%cross_sections(irch)%height, height)
431 this%cross_sections(irch)%height, height)
437 dc0 = (factor(2) - factor(1)) /
dem6
438 dc1 = (factor(3) - factor(2)) /
dem6
441 if (dc0 < dzero .or. dc1 < dzero)
then
442 this%invalid = this%invalid + 1
443 height = unique_heights(n)
445 if (this%cross_sections(irch)%height(i) == height)
then
446 this%cross_sections(irch)%valid(i) = .false.
454 deallocate (unique_heights)
463 subroutine output(this, widths, roughs, kstp, kper)
468 real(DP),
dimension(this%nreaches),
intent(in) :: widths
469 real(DP),
dimension(this%nreaches),
intent(in) :: roughs
470 integer(I4B),
intent(in),
optional :: kstp
471 integer(I4B),
intent(in),
optional :: kper
473 character(len=LINELENGTH) :: title
474 character(len=LINELENGTH) :: text
475 character(len=LINELENGTH) :: filename
476 character(len=10) :: cvalid
477 logical(LGP) :: transient
478 integer(I4B) :: kkstp
479 integer(I4B) :: kkper
482 integer(I4B) :: ntabcols
483 integer(I4B) :: ninvalid_reaches
485 real(DP) :: xfraction
488 integer(I4B),
dimension(this%nreaches) :: reach_fail
495 if (
present(kstp))
then
498 if (
present(kper))
then
503 if (kkstp > 0 .and. kkper > 0)
then
510 do irch = 1, this%nreaches
511 filename = this%filenames(irch)
515 if (trim(adjustl(filename)) /=
'NONE')
then
516 do n = 1, this%npoints(irch)
517 if (.not. this%cross_sections(irch)%valid(n))
then
518 reach_fail(irch) = reach_fail(irch) + 1
525 do irch = 1, this%nreaches
526 filename = this%filenames(irch)
529 if (trim(adjustl(filename)) /=
'NONE')
then
533 if (this%iprpak > 0 .or. reach_fail(irch) > 0)
then
536 if (reach_fail(irch) > 0)
then
543 write (title,
'(a,1x,i0,1x,3a)') &
544 'CROSS_SECTION DATA FOR REACH', irch,
"FROM TAB6 FILE ('", &
545 trim(adjustl(filename)),
"')"
546 call table_cr(this%inputtab, trim(adjustl(filename)), title)
547 call this%inputtab%table_df(this%npoints(irch), ntabcols, &
548 this%iout, finalize=.false., &
551 call this%inputtab%set_kstpkper(kkstp, kkper)
554 call this%inputtab%initialize_column(text, 20, alignment=
tableft)
556 call this%inputtab%initialize_column(text, 20, alignment=
tableft)
558 call this%inputtab%initialize_column(text, 20, alignment=
tableft)
560 call this%inputtab%initialize_column(text, 20, alignment=
tableft)
561 text =
"MANNING'S ROUGHNESS COEFFICIENT"
562 call this%inputtab%initialize_column(text, 20, alignment=
tableft)
563 if (reach_fail(irch) > 0)
then
564 text =
'NEEDS ADJUSTMENT'
565 call this%inputtab%initialize_column(text, 10, alignment=
tableft)
573 do n = 1, this%npoints(irch)
574 xfraction = this%cross_sections(irch)%station(n) / width
575 r = this%cross_sections(irch)%roughfraction(n) * rough
576 call this%inputtab%add_term(xfraction)
577 call this%inputtab%add_term(this%cross_sections(irch)%station(n))
578 call this%inputtab%add_term(this%cross_sections(irch)%height(n))
579 call this%inputtab%add_term(&
580 &this%cross_sections(irch)%roughfraction(n))
581 call this%inputtab%add_term(r)
582 if (reach_fail(irch) > 0)
then
583 if (this%cross_sections(irch)%valid(n))
then
588 call this%inputtab%add_term(cvalid)
593 call this%inputtab%finalize_table()
599 if (this%invalid > 0)
then
601 do irch = 1, this%nreaches
602 if (reach_fail(irch) > 0)
then
603 ninvalid_reaches = ninvalid_reaches + 1
606 write (
warnmsg,
'(a,1x,i0,7(1x,a))') &
607 'Cross-section data for', ninvalid_reaches, &
608 'reaches include one or more points that result in a', &
609 'non-unique depth-conveyance relation. This occurs when', &
610 'there are horizontal sections at non-zero heights', &
611 '(for example, flat overbank sections). This can usually', &
612 'be resolved by adding a small slope to these flat', &
613 'sections. See the cross-section tables in the model', &
614 'listing file for more information.'
629 integer(I4B) :: nptstot
634 do n = 1, this%nreaches
635 nptstot = nptstot + this%npoints(n)
644 subroutine pack(this, ncrossptstot, ncrosspts, iacross, &
645 station, height, roughfraction)
648 integer(I4B),
intent(in) :: ncrossptstot
649 integer(I4B),
dimension(this%nreaches),
intent(inout) :: ncrosspts
650 integer(I4B),
dimension(this%nreaches + 1),
intent(inout) :: iacross
651 real(DP),
dimension(ncrossptstot),
intent(inout) :: station
652 real(DP),
dimension(ncrossptstot),
intent(inout) :: height
653 real(DP),
dimension(ncrossptstot),
intent(inout) :: roughfraction
657 integer(I4B) :: npoints
663 do n = 1, this%nreaches
664 npoints = this%npoints(n)
665 ncrosspts(n) = npoints
667 station(ipos) = this%cross_sections(n)%station(i)
668 height(ipos) = this%cross_sections(n)%height(i)
669 roughfraction(ipos) = this%cross_sections(n)%roughfraction(i)
672 iacross(n + 1) = ipos
688 deallocate (this%npoints)
689 nullify (this%npoints)
690 do n = 1, this%nreaches
691 deallocate (this%cross_sections(n)%npoints)
692 nullify (this%cross_sections(n)%npoints)
693 deallocate (this%cross_sections(n)%station)
694 nullify (this%cross_sections(n)%station)
695 deallocate (this%cross_sections(n)%height)
696 nullify (this%cross_sections(n)%height)
697 deallocate (this%cross_sections(n)%roughfraction)
698 nullify (this%cross_sections(n)%roughfraction)
699 deallocate (this%cross_sections(n)%valid)
700 nullify (this%cross_sections(n)%valid)
702 deallocate (this%cross_sections)
703 nullify (this%cross_sections)
706 if (
associated(this%inputtab))
then
707 call this%inputtab%table_da()
708 deallocate (this%inputtab)
709 nullify (this%inputtab)
713 deallocate (this%invalid)
714 nullify (this%invalid)
718 nullify (this%iprpak)
719 nullify (this%nreaches)
This module contains block parser methods.
This module contains simulation constants.
integer(i4b), parameter linelength
maximum length of a standard line
@ tableft
left justified table column
real(dp), parameter dtwothirds
real constant 2/3
real(dp), parameter dem6
real constant 1e-6
real(dp), parameter dzero
real constant zero
integer(i4b), parameter iuoc
open/close file unit number
real(dp), parameter done
real constant 1
This module contains stateless sfr subroutines and functions.
real(dp) function, public get_cross_section_area(npts, stations, heights, d)
Calculate the cross-sectional area for a reach.
real(dp) function, public get_hydraulic_radius(npts, stations, heights, d)
Calculate the hydraulic radius for a reach.
This module defines variable data types.
integer(i4b) function get_ncrossptstot(this)
Get the total number of cross-section points.
subroutine, public cross_section_cr(this, iout, iprpak, nreaches)
Create a cross-section object.
subroutine output(this, widths, roughs, kstp, kper)
Write cross-section tables.
subroutine destroy(this)
Deallocate the cross-section object.
subroutine read_table(this, irch, width, filename)
Read a cross-section table.
subroutine validate(this, irch)
Validate cross-section tables.
subroutine initialize(this, ncrossptstot, ncrosspts, iacross, station, height, roughfraction)
Initialize a cross-section object.
subroutine pack(this, ncrossptstot, ncrosspts, iacross, station, height, roughfraction)
Pack the cross-section object.
This module contains simulation methods.
subroutine, public store_warning(msg, substring)
Store warning message.
subroutine, public store_error(msg, terminate)
Store an error message.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
character(len=maxcharlen) warnmsg
warning message string
subroutine, public table_cr(this, name, title)