31 character(len=LENMEMPATH) :: memorypath
32 character(len=LENMEMPATH) :: input_mempath =
''
33 character(len=LENMODELNAME),
pointer :: name_model => null()
34 character(len=LINELENGTH),
pointer :: input_fname => null()
35 integer(I4B),
pointer :: inunit => null()
36 integer(I4B),
pointer :: iout => null()
37 integer(I4B),
pointer :: nodes => null()
38 integer(I4B),
pointer :: nodesuser => null()
39 integer(I4B),
pointer :: nja => null()
40 integer(I4B),
pointer :: njas => null()
41 integer(I4B),
pointer :: lenuni => null()
42 integer(I4B),
pointer :: ndim => null()
43 integer(I4B),
pointer :: icondir => null()
44 integer(I4B),
pointer :: nogrb => null()
45 real(dp),
dimension(:),
pointer,
contiguous :: xc => null()
46 real(dp),
dimension(:),
pointer,
contiguous :: yc => null()
47 real(dp),
pointer :: yorigin => null()
48 real(dp),
pointer :: xorigin => null()
49 real(dp),
pointer :: angrot => null()
50 integer(I4B),
dimension(:),
pointer,
contiguous :: mshape => null()
51 real(dp),
dimension(:),
pointer,
contiguous :: top => null()
52 real(dp),
dimension(:),
pointer,
contiguous :: bot => null()
53 real(dp),
dimension(:),
pointer,
contiguous :: area => null()
56 real(dp),
dimension(:),
pointer,
contiguous :: dbuff => null()
57 integer(I4B),
dimension(:),
pointer,
contiguous :: ibuff => null()
58 integer(I4B),
dimension(:),
pointer,
contiguous :: nodereduced => null()
59 integer(I4B),
dimension(:),
pointer,
contiguous :: nodeuser => null()
129 call store_error(
'Programmer error: dis_df must be overridden', &
139 integer(I4B),
intent(in) :: moffset
142 integer(I4B) :: i, j, ipos, iglo, jglo
145 do ipos = this%con%ia(i), this%con%ia(i + 1) - 1
146 j = this%con%ja(ipos)
149 call sparse%addconnection(iglo, jglo, 1)
155 subroutine dis_mc(this, moffset, idxglo, matrix_sln)
158 integer(I4B),
intent(in) :: moffset
159 integer(I4B),
dimension(:),
intent(inout) :: idxglo
162 integer(I4B) :: i, j, ipos, iglo, jglo
166 do ipos = this%con%ia(i), this%con%ia(i + 1) - 1
167 j = this%con%ja(ipos)
169 idxglo(ipos) = matrix_sln%get_position(iglo, jglo)
178 integer(I4B),
dimension(:),
intent(in) :: icelltype
180 integer(I4B),
dimension(:),
allocatable :: ict
181 integer(I4B) :: nu, nr
184 allocate (ict(this%nodesuser))
185 do nu = 1, this%nodesuser
186 nr = this%get_nodenumber(nu, 0)
188 ict(nu) = icelltype(nr)
194 if (this%nogrb == 0)
call this%write_grb(ict)
200 integer(I4B),
dimension(:),
intent(in) :: icelltype
201 call store_error(
'Programmer error: write_grb must be overridden', &
213 deallocate (this%name_model)
214 deallocate (this%input_fname)
242 call this%con%con_da()
243 deallocate (this%con)
249 integer(I4B),
intent(in) :: nodeu
250 character(len=*),
intent(inout) :: str
252 call store_error(
'Programmer error: nodeu_to_string must be overridden', &
259 integer(I4B),
intent(in) :: nodeu
260 integer(I4B),
dimension(:),
intent(inout) :: arr
262 call store_error(
'Programmer error: nodeu_to_array must be overridden', &
269 integer(I4B),
intent(in) :: noder
270 integer(I4B) :: nodenumber
272 if (this%nodes < this%nodesuser)
then
273 nodenumber = this%nodeuser(noder)
281 integer(I4B),
intent(in) :: nodeu
282 integer(I4B),
intent(in) :: icheck
283 integer(I4B) :: nodenumber
286 call store_error(
'Programmer error: get_nodenumber_idx1 must be overridden', &
292 integer(I4B),
intent(in) :: k, j
293 integer(I4B),
intent(in) :: icheck
294 integer(I4B) :: nodenumber
297 call store_error(
'Programmer error: get_nodenumber_idx2 must be overridden', &
303 integer(I4B),
intent(in) :: k, i, j
304 integer(I4B),
intent(in) :: icheck
305 integer(I4B) :: nodenumber
308 call store_error(
'Programmer error: get_nodenumber_idx3 must be overridden', &
317 integer(I4B),
intent(in) :: noden
318 integer(I4B),
intent(in) :: nodem
319 integer(I4B),
intent(in) :: ihc
320 real(DP),
intent(inout) :: xcomp
321 real(DP),
intent(inout) :: ycomp
322 real(DP),
intent(inout) :: zcomp
323 integer(I4B),
intent(in) :: ipos
325 call store_error(
'Programmer error: connection_normal must be overridden', &
333 xcomp, ycomp, zcomp, conlen)
335 integer(I4B),
intent(in) :: noden
336 integer(I4B),
intent(in) :: nodem
337 logical,
intent(in) :: nozee
338 real(DP),
intent(in) :: satn
339 real(DP),
intent(in) :: satm
340 integer(I4B),
intent(in) :: ihc
341 real(DP),
intent(inout) :: xcomp
342 real(DP),
intent(inout) :: ycomp
343 real(DP),
intent(inout) :: zcomp
344 real(DP),
intent(inout) :: conlen
346 call store_error(
'Programmer error: connection_vector must be overridden', &
352 real(dp),
intent(in) :: x
353 real(dp),
intent(in) :: y
354 real(dp),
intent(in) :: xorigin
355 real(dp),
intent(in) :: yorigin
356 real(dp),
intent(in) :: angrot
357 real(dp),
intent(out) :: xglo
358 real(dp),
intent(out) :: yglo
367 if (ang /=
dzero)
then
368 xglo = x * cos(ang) - y * sin(ang)
369 yglo = x * sin(ang) + y * cos(ang)
373 xglo = xglo + xorigin
374 yglo = yglo + yorigin
380 character(len=*),
intent(out) :: dis_type
382 dis_type =
"Not implemented"
383 call store_error(
'Programmer error: get_dis_type must be overridden', &
391 integer(I4B) :: dis_enum
394 call store_error(
'Programmer error: get_dis_enum must be overridden', &
402 character(len=*),
intent(in) :: name_model
403 character(len=*),
intent(in) :: input_mempath
404 logical(LGP) :: found
410 allocate (this%name_model)
411 allocate (this%input_fname)
413 call mem_allocate(this%inunit,
'INUNIT', this%memoryPath)
416 call mem_allocate(this%nodesuser,
'NODESUSER', this%memoryPath)
418 call mem_allocate(this%icondir,
'ICONDIR', this%memoryPath)
420 call mem_allocate(this%xorigin,
'XORIGIN', this%memoryPath)
421 call mem_allocate(this%yorigin,
'YORIGIN', this%memoryPath)
422 call mem_allocate(this%angrot,
'ANGROT', this%memoryPath)
425 call mem_allocate(this%lenuni,
'LENUNI', this%memoryPath)
428 this%name_model = name_model
429 this%input_mempath = input_mempath
430 this%input_fname =
''
447 this%input_mempath, found)
456 call mem_allocate(this%mshape, this%ndim,
'MSHAPE', this%memoryPath)
457 call mem_allocate(this%xc, this%nodes,
'XC', this%memoryPath)
458 call mem_allocate(this%yc, this%nodes,
'YC', this%memoryPath)
459 call mem_allocate(this%top, this%nodes,
'TOP', this%memoryPath)
460 call mem_allocate(this%bot, this%nodes,
'BOT', this%memoryPath)
461 call mem_allocate(this%area, this%nodes,
'AREA', this%memoryPath)
464 this%mshape(1) = this%nodes
467 if (this%nodes < this%nodesuser)
then
468 isize = this%nodesuser
474 call mem_allocate(this%dbuff, isize,
'DBUFF', this%name_model)
475 call mem_allocate(this%ibuff, isize,
'IBUFF', this%name_model)
485 flag_string, allow_zero)
result(nodeu)
488 integer(I4B),
intent(inout) :: lloc
489 integer(I4B),
intent(inout) :: istart
490 integer(I4B),
intent(inout) :: istop
491 integer(I4B),
intent(in) :: in
492 integer(I4B),
intent(in) :: iout
493 character(len=*),
intent(inout) :: line
494 logical,
optional,
intent(in) :: flag_string
495 logical,
optional,
intent(in) :: allow_zero
496 integer(I4B) :: nodeu
499 call store_error(
'Programmer error: nodeu_from_string must be overridden', &
512 allow_zero)
result(nodeu)
515 character(len=*),
intent(inout) :: cellid
516 integer(I4B),
intent(in) :: inunit
517 integer(I4B),
intent(in) :: iout
518 logical,
optional,
intent(in) :: flag_string
519 logical,
optional,
intent(in) :: allow_zero
520 integer(I4B) :: nodeu
523 call store_error(
'Programmer error: nodeu_from_cellid must be overridden', &
535 flag_string)
result(noder)
538 integer(I4B),
intent(inout) :: lloc
539 integer(I4B),
intent(inout) :: istart
540 integer(I4B),
intent(inout) :: istop
541 integer(I4B),
intent(in) :: in
542 integer(I4B),
intent(in) :: iout
543 character(len=*),
intent(inout) :: line
544 logical,
optional,
intent(in) :: flag_string
545 integer(I4B) :: noder
547 integer(I4B) :: nodeu
548 character(len=LINELENGTH) :: nodestr
549 logical :: flag_string_local
551 if (
present(flag_string))
then
552 flag_string_local = flag_string
554 flag_string_local = .false.
556 nodeu = this%nodeu_from_string(lloc, istart, istop, in, iout, line, &
561 noder = this%get_nodenumber(nodeu, 0)
565 if (noder <= 0 .and. .not. flag_string_local)
then
566 call this%nodeu_to_string(nodeu, nodestr)
568 ' Cell is outside active grid domain: '// &
569 trim(adjustl(nodestr))
583 allow_zero)
result(noder)
585 integer(I4B) :: noder
588 character(len=*),
intent(inout) :: cellid
589 integer(I4B),
intent(in) :: inunit
590 integer(I4B),
intent(in) :: iout
591 logical,
optional,
intent(in) :: flag_string
592 logical,
optional,
intent(in) :: allow_zero
594 integer(I4B) :: nodeu
595 logical :: allowzerolocal
596 character(len=LINELENGTH) :: nodestr
597 logical :: flag_string_local
599 if (
present(flag_string))
then
600 flag_string_local = flag_string
602 flag_string_local = .false.
604 if (
present(allow_zero))
then
605 allowzerolocal = allow_zero
607 allowzerolocal = .false.
610 nodeu = this%nodeu_from_cellid(cellid, inunit, iout, flag_string_local, &
615 noder = this%get_nodenumber(nodeu, 0)
619 if (noder <= 0 .and. .not. flag_string_local)
then
620 call this%nodeu_to_string(nodeu, nodestr)
622 ' Cell is outside active grid domain: '// &
623 trim(adjustl(nodestr))
632 call store_error(
'Programmer error: supports_layers must be overridden', &
642 call store_error(
'Programmer error: get_ncpl must be overridden', &
652 integer(I4B),
intent(in) :: n
653 real(dp),
intent(in) :: x
664 thick = (tp - bt) * sat
672 integer(I4B),
intent(in) :: ic
673 real(DP),
allocatable,
intent(out) :: polyverts(:, :)
674 logical(LGP),
intent(in),
optional :: closed
676 errmsg =
'Programmer error: get_polyverts must be overridden'
685 character(len=*),
intent(inout) :: line
686 integer(I4B),
intent(inout) :: lloc
687 integer(I4B),
intent(inout) :: istart
688 integer(I4B),
intent(inout) :: istop
689 integer(I4B),
intent(in) :: in
690 integer(I4B),
intent(in) :: iout
691 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(inout) :: iarray
692 character(len=*),
intent(in) :: aname
694 errmsg =
'Programmer error: read_int_array must be overridden'
703 character(len=*),
intent(inout) :: line
704 integer(I4B),
intent(inout) :: lloc
705 integer(I4B),
intent(inout) :: istart
706 integer(I4B),
intent(inout) :: istop
707 integer(I4B),
intent(in) :: in
708 integer(I4B),
intent(in) :: iout
709 real(DP),
dimension(:),
pointer,
contiguous,
intent(inout) :: darray
710 character(len=*),
intent(in) :: aname
712 errmsg =
'Programmer error: read_dbl_array must be overridden'
720 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(in) :: ibuff1
721 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(inout) :: ibuff2
723 integer(I4B) :: nodeu
724 integer(I4B) :: noder
726 do nodeu = 1, this%nodesuser
727 noder = this%get_nodenumber(nodeu, 0)
728 if (noder <= 0) cycle
729 ibuff2(noder) = ibuff1(nodeu)
737 real(DP),
dimension(:),
pointer,
contiguous,
intent(in) :: buff1
738 real(DP),
dimension(:),
pointer,
contiguous,
intent(inout) :: buff2
740 integer(I4B) :: nodeu
741 integer(I4B) :: noder
743 do nodeu = 1, this%nodesuser
744 noder = this%get_nodenumber(nodeu, 0)
745 if (noder <= 0) cycle
746 buff2(noder) = buff1(nodeu)
757 subroutine read_list(this, line_reader, in, iout, iprpak, nlist, &
758 inamedbound, iauxmultcol, nodelist, rlist, auxvar, &
759 auxname, boundname, label, pkgname, tsManager, iscloc, &
772 integer(I4B),
intent(in) :: in
773 integer(I4B),
intent(in) :: iout
774 integer(I4B),
intent(in) :: iprpak
775 integer(I4B),
intent(inout) :: nlist
776 integer(I4B),
intent(in) :: inamedbound
777 integer(I4B),
intent(in) :: iauxmultcol
778 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(inout) :: nodelist
779 real(DP),
dimension(:, :),
pointer,
contiguous,
intent(inout) :: rlist
780 real(DP),
dimension(:, :),
pointer,
contiguous,
intent(inout) :: auxvar
781 character(len=LENAUXNAME),
dimension(:),
intent(inout) :: auxname
782 character(len=LENBOUNDNAME),
dimension(:),
pointer,
contiguous, &
783 intent(inout) :: boundname
784 character(len=*),
intent(in) :: label
785 character(len=*),
intent(in) :: pkgName
787 integer(I4B),
intent(in) :: iscloc
788 integer(I4B),
intent(in),
optional :: indxconvertflux
791 integer(I4B) :: nodeu, noder
792 character(len=LINELENGTH) :: nodestr
793 integer(I4B) :: ii, jj
794 real(DP),
pointer :: bndElem => null()
800 call lstrdobj%read_list(line_reader, in, iout, nlist, inamedbound, &
801 this%mshape, nodelist, rlist, auxvar, auxname, &
806 if (lstrdobj%ntxtrlist > 0)
then
807 do l = 1, lstrdobj%ntxtrlist
808 ii = lstrdobj%idxtxtrow(l)
809 jj = lstrdobj%idxtxtcol(l)
811 bndelem => rlist(jj, ii)
813 pkgname,
'BND', tsmanager, iprpak, &
815 if (
associated(tslinkbnd))
then
820 if (iauxmultcol > 0 .and. jj == iscloc)
then
821 tslinkbnd%RMultiplier => auxvar(iauxmultcol, ii)
825 if (lstrdobj%inamedbound == 1)
then
826 tslinkbnd%BndName = lstrdobj%boundname(tslinkbnd%IRow)
831 if (
present(indxconvertflux))
then
832 if (indxconvertflux == jj)
then
833 tslinkbnd%convertflux = .true.
835 noder = this%get_nodenumber(nodeu, 0)
836 tslinkbnd%CellArea = this%get_area(noder)
845 if (lstrdobj%ntxtauxvar > 0)
then
846 do l = 1, lstrdobj%ntxtauxvar
847 ii = lstrdobj%idxtxtauxrow(l)
848 jj = lstrdobj%idxtxtauxcol(l)
850 bndelem => auxvar(jj, ii)
852 pkgname,
'AUX', tsmanager, iprpak, &
854 if (lstrdobj%inamedbound == 1)
then
855 if (
associated(tslinkaux))
then
856 tslinkaux%BndName = lstrdobj%boundname(tslinkaux%IRow)
863 if (iauxmultcol > 0)
then
865 rlist(iscloc, l) = rlist(iscloc, l) * auxvar(iauxmultcol, l)
870 if (iprpak /= 0)
then
871 call lstrdobj%write_list()
877 if (this%nodes < this%nodesuser)
then
880 noder = this%get_nodenumber(nodeu, 0)
882 call this%nodeu_to_string(nodeu, nodestr)
884 ' Cell is outside active grid domain: '// &
885 trim(adjustl(nodestr))
906 icolbnd, aname, inunit, iout)
909 integer(I4B),
intent(in) :: ncolbnd
910 integer(I4B),
intent(in) :: maxbnd
911 integer(I4B),
dimension(maxbnd) :: nodelist
912 real(DP),
dimension(ncolbnd, maxbnd),
intent(inout) :: darray
913 integer(I4B),
intent(in) :: icolbnd
914 character(len=*),
intent(in) :: aname
915 integer(I4B),
intent(in) :: inunit
916 integer(I4B),
intent(in) :: iout
918 errmsg =
'Programmer error: read_layer_array must be overridden'
927 cdatafmp, nvaluesp, nwidthp, editdesc, dinact)
930 real(DP),
dimension(:),
pointer,
contiguous,
intent(inout) :: darray
931 integer(I4B),
intent(in) :: iout
932 integer(I4B),
intent(in) :: iprint
933 integer(I4B),
intent(in) :: idataun
934 character(len=*),
intent(in) :: aname
935 character(len=*),
intent(in) :: cdatafmp
936 integer(I4B),
intent(in) :: nvaluesp
937 integer(I4B),
intent(in) :: nwidthp
938 character(len=*),
intent(in) :: editdesc
939 real(DP),
intent(in) :: dinact
941 errmsg =
'Programmer error: record_array must be overridden'
949 real(DP),
dimension(:),
intent(in) :: flowja
950 integer(I4B),
intent(in) :: ibinun
951 integer(I4B),
intent(in) :: iout
953 character(len=16),
dimension(1) :: text
955 data text(1)/
' FLOW-JA-FACE'/
958 call ubdsv1(
kstp,
kper, text(1), ibinun, flowja,
size(flowja), 1, 1, &
966 integer(I4B),
intent(in) :: noder
967 character(len=*),
intent(inout) :: str
969 integer(I4B) :: nodeu
971 nodeu = this%get_nodeuser(noder)
972 call this%nodeu_to_string(nodeu, str)
979 integer(I4B),
intent(in) :: noder
980 integer(I4B),
dimension(:),
intent(inout) :: arr
982 integer(I4B) :: nodeu
984 nodeu = this%get_nodeuser(noder)
985 call this%nodeu_to_array(nodeu, arr)
990 dstmodel, dstpackage, naux, auxtxt, &
993 character(len=16),
intent(in) :: text
994 character(len=16),
intent(in) :: textmodel
995 character(len=16),
intent(in) :: textpackage
996 character(len=16),
intent(in) :: dstmodel
997 character(len=16),
intent(in) :: dstpackage
998 integer(I4B),
intent(in) :: naux
999 character(len=16),
dimension(:),
intent(in) :: auxtxt
1000 integer(I4B),
intent(in) :: ibdchn
1001 integer(I4B),
intent(in) :: nlist
1002 integer(I4B),
intent(in) :: iout
1004 errmsg =
'Programmer error: record_srcdst_list_header must be overridden'
1010 naux, aux, olconv, olconv2)
1013 integer(I4B),
intent(in) :: ibdchn
1014 integer(I4B),
intent(in) :: noder
1015 integer(I4B),
intent(in) :: noder2
1016 real(DP),
intent(in) :: q
1017 integer(I4B),
intent(in) :: naux
1018 real(DP),
dimension(naux),
intent(in) :: aux
1019 logical,
optional,
intent(in) :: olconv
1020 logical,
optional,
intent(in) :: olconv2
1024 integer(I4B) :: nodeu
1025 integer(I4B) :: nodeu2
1028 if (
present(olconv))
then
1034 nodeu = this%get_nodeuser(noder)
1038 if (
present(olconv2))
then
1044 nodeu2 = this%get_nodeuser(noder2)
1048 call ubdsvd(ibdchn, nodeu, nodeu2, q, naux, aux)
1057 integer(I4B),
intent(in) :: maxbnd
1058 integer(I4B),
dimension(:),
pointer,
contiguous :: darray
1059 integer(I4B),
dimension(maxbnd),
intent(inout) :: nodelist
1060 integer(I4B),
intent(inout) :: nbound
1061 character(len=*),
intent(in) :: aname
1063 errmsg =
'Programmer error: nlarray_to_nodelist must be overridden'
1071 integer(I4B),
intent(inout) :: n
1072 integer(I4B),
dimension(:),
intent(in) :: ibound
1074 integer(I4B) :: m, ii, iis
1075 logical done, bottomcell
1080 do while (.not. done)
1082 cloop:
do ii = this%con%ia(n) + 1, this%con%ia(n + 1) - 1
1084 iis = this%con%jas(ii)
1085 if (this%con%ihc(iis) == 0 .and. m > n)
then
1088 bottomcell = .false.
1091 if (ibound(m) /= 0)
then
1101 if (bottomcell) done = .true.
1108 integer(I4B),
intent(in) :: node
1111 area = this%area(node)
1125 real(dp) :: area_factor
1128 integer(I4B),
intent(in) :: node
1129 integer(I4B),
intent(in) :: idx_conn
1131 real(dp) :: area_node
1132 real(dp) :: area_conn
1135 area_node = this%area(node)
1136 area_conn = this%con%hwva(idx_conn)
1139 area_factor = area_conn / area_node
1152 integer(I4B),
intent(in) :: n
1153 integer(I4B),
intent(in) :: m
1154 integer(I4B),
intent(in) :: idx_conn
1155 real(DP),
intent(out) :: width_n
1156 real(DP),
intent(out) :: width_m
1158 integer(I4B) :: isympos
1161 isympos = this%con%jas(idx_conn)
1162 width_n = this%con%hwva(isympos)
1174 select case (this%get_dis_enum())
1187 select case (this%get_dis_enum())
1200 select case (this%get_dis_enum())
real(dp) function get_area(this, node)
Return the cell area for the given node.
subroutine dis_df(this)
Define the discretization.
subroutine noder_to_string(this, noder, str)
Convert reduced node number to string (nodenumber), (k,j) or (k,i,j)
integer(i4b) function get_nodenumber_idx2(this, k, j, icheck)
real(dp) function get_area_factor(this, node, idx_conn)
@ brief Calculate the area factor for the cell connection
subroutine read_layer_array(this, nodelist, darray, ncolbnd, maxbnd, icolbnd, aname, inunit, iout)
Read a 2d double array into col icolbnd of darray.
integer(i4b) function nodeu_from_string(this, lloc, istart, istop, in, iout, line, flag_string, allow_zero)
Convert a string to a user nodenumber.
subroutine allocate_scalars(this, name_model, input_mempath)
Allocate and initialize scalar variables.
subroutine dis_mc(this, moffset, idxglo, matrix_sln)
Map cell connections in the numerical solution coefficient matrix.
logical(lgp) function is_1d(this)
@Brief return true if grid is one dimensional
subroutine highest_active(this, n, ibound)
Find the first highest active cell beneath cell n.
subroutine record_array(this, darray, iout, iprint, idataun, aname, cdatafmp, nvaluesp, nwidthp, editdesc, dinact)
Record a double precision array.
subroutine record_srcdst_list_header(this, text, textmodel, textpackage, dstmodel, dstpackage, naux, auxtxt, ibdchn, nlist, iout)
Record list header for imeth=6.
integer(i4b) function get_nodenumber_idx3(this, k, i, j, icheck)
subroutine record_srcdst_list_entry(this, ibdchn, noder, noder2, q, naux, aux, olconv, olconv2)
Record list header.
subroutine write_grb(this, icelltype)
Write a binary grid file.
subroutine, public dis_transform_xy(x, y, xorigin, yorigin, angrot, xglo, yglo)
Get global (x, y) coordinates from cell-local coordinates.
subroutine read_int_array(this, line, lloc, istart, istop, iout, in, iarray, aname)
Read an integer array.
logical(lgp) function is_3d(this)
@Brief return true if grid is three dimensional
logical function supports_layers(this)
Indicates whether the grid discretization supports layers.
subroutine connection_normal(this, noden, nodem, ihc, xcomp, ycomp, zcomp, ipos)
Get normal vector components between the cell and a given neighbor. The normal points outward from th...
subroutine nodeu_to_string(this, nodeu, str)
Convert a user nodenumber to a string (nodenumber), (k,j), or (k,i,j)
integer(i4b) function noder_from_string(this, lloc, istart, istop, in, iout, line, flag_string)
Convert a string to a reduced nodenumber.
subroutine fill_int_array(this, ibuff1, ibuff2)
Fill an integer array.
subroutine get_polyverts(this, ic, polyverts, closed)
Get a 2D array of polygon vertices, listed in clockwise order beginning with the lower left corner.
subroutine dis_da(this)
@brier Deallocate variables
subroutine fill_dbl_array(this, buff1, buff2)
Fill a double precision array.
subroutine dis_ac(this, moffset, sparse)
Add connections to sparse cell connectivity matrix.
subroutine record_connection_array(this, flowja, ibinun, iout)
Record a connection-based double precision array.
subroutine nodeu_to_array(this, nodeu, arr)
Convert a user nodenumber to an array (nodenumber), (k,j), or (k,i,j)
integer(i4b) function get_nodenumber_idx1(this, nodeu, icheck)
integer(i4b) function nodeu_from_cellid(this, cellid, inunit, iout, flag_string, allow_zero)
Convert a cellid string to a user nodenumber.
logical(lgp) function is_2d(this)
@Brief return true if grid is two dimensional
subroutine noder_to_array(this, noder, arr)
Convert reduced node number to array (nodenumber), (k,j) or (k,i,j)
subroutine get_dis_type(this, dis_type)
Get the discretization type (DIS, DISV, or DISU)
integer(i4b) function get_nodeuser(this, noder)
Convert a reduced nodenumber to a user node number.
subroutine get_flow_width(this, n, m, idx_conn, width_n, width_m)
@ brief Calculate the flow width between two cells
subroutine read_dbl_array(this, line, lloc, istart, istop, iout, in, darray, aname)
Read a double precision array.
subroutine dis_ar(this, icelltype)
Allocate and setup variables, and write binary grid file.
subroutine allocate_arrays(this)
Allocate and initialize arrays.
integer(i4b) function get_dis_enum(this)
Get the discretization type enumeration.
integer(i4b) function noder_from_cellid(this, cellid, inunit, iout, flag_string, allow_zero)
Convert cellid string to reduced nodenumber.
subroutine nlarray_to_nodelist(this, darray, nodelist, maxbnd, nbound, aname)
Convert an integer array to nodelist.
integer(i4b) function get_ncpl(this)
Return number of cells per layer. This is nodes for a DISU grid, as there are no layers.
real(dp) function get_cell_volume(this, n, x)
Return volume of cell n based on x value passed.
subroutine read_list(this, line_reader, in, iout, iprpak, nlist, inamedbound, iauxmultcol, nodelist, rlist, auxvar, auxname, boundname, label, pkgname, tsManager, iscloc, indxconvertflux)
Read a list using the list reader.
subroutine connection_vector(this, noden, nodem, nozee, satn, satm, ihc, xcomp, ycomp, zcomp, conlen)
Get unit vector components between the cell and a given neighbor. Saturation must be provided to comp...
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
@ disu
DISV6 discretization.
@ dis
DIS6 discretization.
@ dis1d
DIS1D6 discretization.
@ disv2d
DISV2D6 discretization.
@ disv1d
DISV1D6 discretization.
@ dis2d
DIS2D6 discretization.
@ disv
DISU6 discretization.
@ disu1d
DISU1D6 discretization.
@ disundef
undefined discretization
@ disu2d
DISU2D6 discretization.
integer(i4b), parameter lenauxname
maximum length of a aux variable
integer(i4b), parameter lenboundname
maximum length of a bound name
real(dp), parameter dpio180
real constant
real(dp), parameter dzero
real constant zero
integer(i4b), parameter lenmempath
maximum length of the memory path
This module defines variable data types.
Generic List Reader Module.
This module contains the LongLineReaderType.
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
real(dp) function squadraticsaturation(top, bot, x, eps)
@ brief sQuadraticSaturation
real(dp), pointer, public pertim
time relative to start of stress period
real(dp), pointer, public totim
time relative to start of simulation
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
subroutine, public read_value_or_time_series(textInput, ii, jj, bndElem, pkgName, auxOrBnd, tsManager, iprpak, tsLink)
Call this subroutine if the time-series link is available or needed.