20 character(len=LENMEMPATH) :: memorypath
21 character(len=LENMODELNAME),
pointer :: name_model => null()
22 integer(I4B),
pointer :: nodes => null()
23 integer(I4B),
pointer :: nja => null()
24 integer(I4B),
pointer :: njas => null()
25 integer(I4B),
pointer :: ianglex => null()
26 integer(I4B),
dimension(:),
pointer,
contiguous :: ia => null()
27 integer(I4B),
dimension(:),
pointer,
contiguous :: ja => null()
28 integer(I4B),
dimension(:),
pointer,
contiguous :: mask => null()
29 real(dp),
dimension(:),
pointer,
contiguous :: cl1 => null()
30 real(dp),
dimension(:),
pointer,
contiguous :: cl2 => null()
31 real(dp),
dimension(:),
pointer,
contiguous :: hwva => null()
32 real(dp),
dimension(:),
pointer,
contiguous :: anglex => null()
33 integer(I4B),
dimension(:),
pointer,
contiguous :: isym => null()
34 integer(I4B),
dimension(:),
pointer,
contiguous :: jas => null()
35 integer(I4B),
dimension(:),
pointer,
contiguous :: ihc => null()
36 integer(I4B),
dimension(:),
pointer,
contiguous :: iausr => null()
37 integer(I4B),
dimension(:),
pointer,
contiguous :: jausr => null()
68 deallocate (this%name_model)
77 if (
associated(this%iausr, this%ia))
then
82 if (
associated(this%jausr, this%ja))
then
88 if (
associated(this%mask, this%ja))
then
117 character(len=*),
intent(in) :: name_model
120 allocate (this%name_model)
126 call mem_allocate(this%ianglex,
'IANGLEX', this%memoryPath)
127 this%name_model = name_model
145 call mem_allocate(this%ia, this%nodes + 1,
'IA', this%memoryPath)
146 call mem_allocate(this%ja, this%nja,
'JA', this%memoryPath)
147 call mem_allocate(this%isym, this%nja,
'ISYM', this%memoryPath)
148 call mem_allocate(this%jas, this%nja,
'JAS', this%memoryPath)
149 call mem_allocate(this%hwva, this%njas,
'HWVA', this%memoryPath)
150 call mem_allocate(this%anglex, this%njas,
'ANGLEX', this%memoryPath)
151 call mem_allocate(this%ihc, this%njas,
'IHC', this%memoryPath)
152 call mem_allocate(this%cl1, this%njas,
'CL1', this%memoryPath)
153 call mem_allocate(this%cl2, this%njas,
'CL2', this%memoryPath)
154 call mem_allocate(this%iausr, 1,
'IAUSR', this%memoryPath)
155 call mem_allocate(this%jausr, 1,
'JAUSR', this%memoryPath)
173 integer(I4B),
dimension(:),
intent(in) :: ihctemp
174 real(DP),
dimension(:),
intent(in) :: cl12temp
175 real(DP),
dimension(:),
intent(in) :: hwvatemp
176 real(DP),
dimension(:),
intent(in) :: angldegx
178 integer(I4B) :: ii, n, m
179 integer(I4B),
parameter :: nname = 6
180 character(len=24),
dimension(nname) :: aname(nname)
182 character(len=*),
parameter :: fmtsymerr = &
183 &
"('Error in array: ',a,'.', &
184 &' Array is not symmetric in positions: ',i0,' and ',i0,'.', &
185 &' Values in these positions are: ',1pg15.6,' and ', 1pg15.6, &
186 &' For node ',i0,' connected to node ',i0)"
187 character(len=*),
parameter :: fmtsymerrja = &
188 &
"('Error in array: ',a,'.', &
189 &' Array does not have symmetric counterpart in position ',i0, &
190 &' for cell ',i0,' connected to cell ',i0)"
191 character(len=*),
parameter :: fmtjanmerr = &
192 &
"('Error in array: ',a,'.', &
193 &' First value for cell : ',i0,' must equal ',i0,'.', &
194 &' Found ',i0,' instead.')"
195 character(len=*),
parameter :: fmtjasorterr = &
196 &
"('Error in array: ',a,'.', &
197 &' Entries not sorted for row: ',i0,'.', &
198 &' Offending entries are: ',i0,' and ',i0)"
199 character(len=*),
parameter :: fmtihcerr = &
200 "('IHC must be 0, 1, or 2. Found: ',i0)"
202 data aname(1)/
' IAC'/
204 data aname(3)/
' IHC'/
205 data aname(4)/
' CL12'/
206 data aname(5)/
' HWVA'/
207 data aname(6)/
' ANGLDEGX'/
211 if (this%ja(ii) < 0) this%ja(ii) = -this%ja(ii)
216 m = this%ja(this%ia(n))
218 write (
errmsg, fmtjanmerr) trim(adjustl(aname(2))), n, n, m
221 do ii = this%ia(n) + 1, this%ia(n + 1) - 2
223 if (m > this%ja(ii + 1))
then
224 write (
errmsg, fmtjasorterr) trim(adjustl(aname(2))), n, &
231 call this%parser%StoreErrorUnit()
235 call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
240 do ii = this%ia(n), this%ia(n + 1) - 1
242 if (this%isym(ii) == 0)
then
243 write (
errmsg, fmtsymerrja) trim(adjustl(aname(2))), ii, n, m
249 call this%parser%StoreErrorUnit()
253 call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
257 do ii = this%ia(n) + 1, this%ia(n + 1) - 1
259 if (ihctemp(ii) /= ihctemp(this%isym(ii)))
then
260 write (
errmsg, fmtsymerr) trim(adjustl(aname(3))), ii, this%isym(ii), &
261 ihctemp(ii), ihctemp(this%isym(ii)), n, m
264 this%ihc(this%jas(ii)) = ihctemp(ii)
269 call this%parser%StoreErrorUnit()
274 do ii = this%ia(n) + 1, this%ia(n + 1) - 1
277 this%cl1(this%jas(ii)) = cl12temp(ii)
279 this%cl2(this%jas(ii)) = cl12temp(ii)
291 do ii = this%ia(n) + 1, this%ia(n + 1) - 1
293 if (hwvatemp(ii) /= hwvatemp(this%isym(ii)))
then
294 write (
errmsg, fmtsymerr) trim(adjustl(aname(5))), ii, this%isym(ii), &
295 hwvatemp(ii), hwvatemp(this%isym(ii)), n, m
298 if (ihctemp(ii) < 0 .or. ihctemp(ii) > 2)
then
299 write (
errmsg, fmtihcerr) ihctemp(ii)
302 this%hwva(this%jas(ii)) = hwvatemp(ii)
306 call this%parser%StoreErrorUnit()
310 if (this%ianglex /= 0)
then
312 do ii = this%ia(n) + 1, this%ia(n + 1) - 1
315 this%anglex(this%jas(ii)) = angldegx(ii) *
dpio180
319 do n = 1,
size(this%anglex)
337 character(len=*),
intent(in) :: name_model
338 integer(I4B),
intent(in) :: nodes
339 integer(I4B),
intent(in) :: nja
340 integer(I4B),
intent(in) :: iout
342 character(len=LINELENGTH) :: line
343 character(len=LINELENGTH) :: keyword
344 integer(I4B) :: ii, n, m
346 logical :: isfound, endOfBlock
347 integer(I4B),
parameter :: nname = 2
348 logical,
dimension(nname) :: lname
349 character(len=24),
dimension(nname) :: aname(nname)
351 character(len=*),
parameter :: fmtsymerr = &
352 &
"(/,'Error in array: ',(a),/, &
353 &'Array is not symmetric in positions: ',2i9,/, &
354 &'Values in these positions are: ', 2(1pg15.6))"
355 character(len=*),
parameter :: fmtihcerr = &
356 &
"(/,'IHC must be 0, 1, or 2. Found: ',i0)"
358 data aname(1)/
' IAC'/
362 call this%allocate_scalars(name_model)
365 this%njas = (this%nja - this%nodes) / 2
368 call this%allocate_arrays()
371 call this%parser%GetBlock(
'CONNECTIONDATA', isfound, ierr)
374 write (iout,
'(1x,a)')
'PROCESSING CONNECTIONDATA'
376 call this%parser%GetNextLine(endofblock)
378 call this%parser%GetStringCaps(keyword)
379 select case (keyword)
381 call readarray(this%parser%iuactive, this%ia, aname(1), 1, &
385 call readarray(this%parser%iuactive, this%ja, aname(2), 1, &
390 'Unknown CONNECTIONDATA tag: ', trim(keyword)
392 call this%parser%StoreErrorUnit()
395 write (iout,
'(1x,a)')
'END PROCESSING CONNECTIONDATA'
397 call store_error(
'Required CONNECTIONDATA block not found.')
398 call this%parser%StoreErrorUnit()
403 if (.not. lname(n))
then
405 'Required input was not specified: ', aname(n)
410 call this%parser%StoreErrorUnit()
414 do n = 2, this%nodes + 1
415 this%ia(n) = this%ia(n) + this%ia(n - 1)
417 do n = this%nodes + 1, 2, -1
418 this%ia(n) = this%ia(n - 1) + 1
424 if (this%ja(ii) < 0) this%ja(ii) = -this%ja(ii)
428 call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
429 call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, &
434 do ii = this%ia(n), this%ia(n + 1) - 1
436 if (n /= this%ja(this%isym(ii)))
then
437 write (line, fmtsymerr) aname(2), ii, this%isym(ii)
439 call this%parser%StoreErrorUnit()
445 call this%parser%StoreErrorUnit()
459 real(DP),
dimension(:),
intent(in) :: fleng
461 integer(I4B) :: n, m, ii
465 do ii = this%ia(n) + 1, this%ia(n + 1) - 1
467 this%cl1(this%jas(ii)) = fleng(n) *
dhalf
468 this%cl2(this%jas(ii)) = fleng(m) *
dhalf
480 nrsize, delr, delc, top, bot, nodereduced, &
487 character(len=*),
intent(in) :: name_model
488 integer(I4B),
intent(in) :: nodes
489 integer(I4B),
intent(in) :: ncol
490 integer(I4B),
intent(in) :: nrow
491 integer(I4B),
intent(in) :: nlay
492 integer(I4B),
intent(in) :: nrsize
493 real(DP),
dimension(ncol),
intent(in) :: delr
494 real(DP),
dimension(nrow),
intent(in) :: delc
495 real(DP),
dimension(nodes),
intent(in) :: top
496 real(DP),
dimension(nodes),
intent(in) :: bot
497 integer(I4B),
dimension(:),
target,
intent(in) :: nodereduced
498 integer(I4B),
dimension(:),
intent(in) :: nodeuser
500 integer(I4B),
dimension(:, :, :),
pointer :: nrdcd_ptr => null()
501 integer(I4B),
dimension(:),
allocatable :: rowmaxnnz
503 integer(I4B) :: i, j, k, kk, ierror, isympos, nodesuser
504 integer(I4B) :: nr, mr
507 call this%allocate_scalars(name_model)
514 allocate (rowmaxnnz(this%nodes))
518 call sparse%init(this%nodes, this%nodes, rowmaxnnz)
521 if (nrsize /= 0)
then
522 nrdcd_ptr(1:ncol, 1:nrow, 1:nlay) => nodereduced
532 if (nrsize == 0)
then
533 nr =
get_node(k, i, j, nlay, nrow, ncol)
535 nr = nrdcd_ptr(j, i, k)
540 call sparse%addconnection(nr, nr, 1)
545 if (nrsize == 0)
then
546 mr =
get_node(kk, i, j, nlay, nrow, ncol)
548 mr = nrdcd_ptr(j, i, kk)
553 call sparse%addconnection(nr, mr, 1)
559 if (nrsize == 0)
then
560 mr =
get_node(k, i - 1, j, nlay, nrow, ncol)
562 mr = nrdcd_ptr(j, i - 1, k)
565 call sparse%addconnection(nr, mr, 1)
571 if (nrsize == 0)
then
572 mr =
get_node(k, i, j - 1, nlay, nrow, ncol)
574 mr = nrdcd_ptr(j - 1, i, k)
577 call sparse%addconnection(nr, mr, 1)
583 if (nrsize == 0)
then
584 mr =
get_node(k, i, j + 1, nlay, nrow, ncol)
586 mr = nrdcd_ptr(j + 1, i, k)
589 call sparse%addconnection(nr, mr, 1)
595 if (nrsize == 0)
then
596 mr =
get_node(k, i + 1, j, nlay, nrow, ncol)
598 mr = nrdcd_ptr(j, i + 1, k)
601 call sparse%addconnection(nr, mr, 1)
608 if (nrsize == 0)
then
609 mr =
get_node(kk, i, j, nlay, nrow, ncol)
611 mr = nrdcd_ptr(j, i, kk)
616 call sparse%addconnection(nr, mr, 1)
622 this%nja = sparse%nnz
623 this%njas = (this%nja - this%nodes) / 2
626 call this%allocate_arrays()
629 call sparse%filliaja(this%ia, this%ja, ierror)
630 call sparse%destroy()
633 call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
634 call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
643 if (nrsize == 0)
then
644 nr =
get_node(k, i, j, nlay, nrow, ncol)
646 nr = nrdcd_ptr(j, i, k)
652 if (nrsize == 0)
then
653 mr =
get_node(k, i, j + 1, nlay, nrow, ncol)
655 mr = nrdcd_ptr(j + 1, i, k)
658 this%ihc(isympos) = 1
659 this%cl1(isympos) =
dhalf * delr(j)
660 this%cl2(isympos) =
dhalf * delr(j + 1)
661 this%hwva(isympos) = delc(i)
662 this%anglex(isympos) =
dzero
663 isympos = isympos + 1
669 if (nrsize == 0)
then
670 mr =
get_node(k, i + 1, j, nlay, nrow, ncol)
672 mr = nrdcd_ptr(j, i + 1, k)
675 this%ihc(isympos) = 1
676 this%cl1(isympos) =
dhalf * delc(i)
677 this%cl2(isympos) =
dhalf * delc(i + 1)
678 this%hwva(isympos) = delr(j)
680 isympos = isympos + 1
687 if (nrsize == 0)
then
688 mr =
get_node(kk, i, j, nlay, nrow, ncol)
690 mr = nrdcd_ptr(j, i, kk)
695 this%ihc(isympos) = 0
696 this%cl1(isympos) =
dhalf * (top(nr) - bot(nr))
697 this%cl2(isympos) =
dhalf * (top(mr) - bot(mr))
698 this%hwva(isympos) = delr(j) * delc(i)
699 this%anglex(isympos) =
dzero
700 isympos = isympos + 1
708 deallocate (rowmaxnnz)
712 nodesuser = nlay * nrow * ncol
713 call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
722 nvert, vertex, iavert, javert, cellxy, &
723 top, bot, nodereduced, nodeuser)
731 character(len=*),
intent(in) :: name_model
732 integer(I4B),
intent(in) :: nodes
733 integer(I4B),
intent(in) :: ncpl
734 integer(I4B),
intent(in) :: nlay
735 integer(I4B),
intent(in) :: nrsize
736 integer(I4B),
intent(in) :: nvert
737 real(DP),
dimension(2, nvert),
intent(in) :: vertex
738 integer(I4B),
dimension(:),
intent(in) :: iavert
739 integer(I4B),
dimension(:),
intent(in) :: javert
740 real(DP),
dimension(2, ncpl),
intent(in) :: cellxy
741 real(DP),
dimension(nodes),
intent(in) :: top
742 real(DP),
dimension(nodes),
intent(in) :: bot
743 integer(I4B),
dimension(:),
intent(in) :: nodereduced
744 integer(I4B),
dimension(:),
intent(in) :: nodeuser
746 integer(I4B),
dimension(:),
allocatable :: itemp
748 integer(I4B) :: n, m, ipos, i, j, ierror, nodesuser
752 call this%allocate_scalars(name_model)
759 call cell1%init(nlay, ncpl, nodes, top, bot, iavert, javert, vertex, &
760 cellxy, nodereduced, nodeuser)
761 call cell2%init(nlay, ncpl, nodes, top, bot, iavert, javert, vertex, &
762 cellxy, nodereduced, nodeuser)
767 allocate (itemp(nvert))
771 call vertcellspm%init(nvert, ncpl, itemp)
774 do i = iavert(j), iavert(j + 1) - 1
775 call vertcellspm%addconnection(javert(i), j, 1)
780 call vertexconnect(this%nodes, nrsize, 6, nlay, ncpl, sparse, &
781 vertcellspm, cell1, cell2, nodereduced)
782 this%nja = sparse%nnz
783 this%njas = (this%nja - this%nodes) / 2
786 call this%allocate_arrays()
790 call sparse%filliaja(this%ia, this%ja, ierror)
791 call sparse%destroy()
794 call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
795 call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
799 call cell1%set_nodered(n)
800 do ipos = this%ia(n) + 1, this%ia(n + 1) - 1
803 call cell2%set_nodered(m)
804 call cell1%cprops(cell2, this%hwva(this%jas(ipos)), &
805 this%cl1(this%jas(ipos)), this%cl2(this%jas(ipos)), &
806 this%anglex(this%jas(ipos)), &
807 this%ihc(this%jas(ipos)))
813 nodesuser = nlay * ncpl
814 call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
824 nodereduced, nodeuser, iainp, jainp, &
825 ihcinp, cl12inp, hwvainp, angldegxinp, &
833 character(len=*),
intent(in) :: name_model
834 integer(I4B),
intent(in) :: nodes
835 integer(I4B),
intent(in) :: nodesuser
836 integer(I4B),
intent(in) :: nrsize
837 integer(I4B),
dimension(:),
contiguous,
intent(in) :: nodereduced
838 integer(I4B),
dimension(:),
contiguous,
intent(in) :: nodeuser
839 integer(I4B),
dimension(:),
contiguous,
intent(in) :: iainp
840 integer(I4B),
dimension(:),
contiguous,
intent(in) :: jainp
841 integer(I4B),
dimension(:),
contiguous,
intent(in) :: ihcinp
842 real(DP),
dimension(:),
contiguous,
intent(in) :: cl12inp
843 real(DP),
dimension(:),
contiguous,
intent(in) :: hwvainp
844 real(DP),
dimension(:),
contiguous,
intent(in) :: angldegxinp
845 integer(I4B),
intent(in) :: iangledegx
847 integer(I4B),
dimension(:),
allocatable :: ihctemp
848 real(DP),
dimension(:),
allocatable :: cl12temp
849 real(DP),
dimension(:),
allocatable :: hwvatemp
850 real(DP),
dimension(:),
allocatable :: angldegxtemp
851 integer(I4B) :: nr, nu, mr, mu, ipos, iposr, ierror
852 integer(I4B),
dimension(:),
allocatable :: rowmaxnnz
856 call this%allocate_scalars(name_model)
860 this%ianglex = iangledegx
864 if (nrsize == 0)
then
866 this%nja =
size(jainp)
867 this%njas = (this%nja - this%nodes) / 2
868 call this%allocate_arrays()
870 this%ia(nu) = iainp(nu)
872 do ipos = 1, this%nja
873 this%ja(ipos) = jainp(ipos)
878 call this%con_finalize(ihcinp, cl12inp, hwvainp, angldegxinp)
884 allocate (rowmaxnnz(this%nodes))
885 do nr = 1, this%nodes
887 rowmaxnnz(nr) = iainp(nu + 1) - iainp(nu)
889 call sparse%init(this%nodes, this%nodes, rowmaxnnz)
894 if (nr > 0)
call sparse%addconnection(nr, nr, 1)
895 do ipos = iainp(nu) + 1, iainp(nu + 1) - 1
900 call sparse%addconnection(nr, mr, 1)
903 this%nja = sparse%nnz
904 this%njas = (this%nja - this%nodes) / 2
907 call this%allocate_arrays()
911 call sparse%filliaja(this%ia, this%ja, ierror)
912 call sparse%destroy()
913 deallocate (rowmaxnnz)
916 allocate (ihctemp(this%nja))
917 allocate (cl12temp(this%nja))
918 allocate (hwvatemp(this%nja))
919 allocate (angldegxtemp(this%nja))
925 do ipos = iainp(nu), iainp(nu + 1) - 1
928 if (nr < 1 .or. mr < 1) cycle
929 ihctemp(iposr) = ihcinp(ipos)
930 cl12temp(iposr) = cl12inp(ipos)
931 hwvatemp(iposr) = hwvainp(ipos)
932 angldegxtemp(iposr) = angldegxinp(ipos)
938 call this%con_finalize(ihctemp, cl12temp, hwvatemp, angldegxtemp)
942 deallocate (cl12temp)
943 deallocate (hwvatemp)
944 deallocate (angldegxtemp)
949 call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
963 vertices, iavert, javert, &
964 cellxyz, cellfdc, nodereduced, nodeuser, &
973 character(len=*),
intent(in) :: name_model
974 integer(I4B),
intent(in) :: nodes
975 integer(I4B),
intent(in) :: nodesuser
976 integer(I4B),
intent(in) :: nrsize
977 integer(I4B),
intent(in) :: nvert
978 real(DP),
dimension(3, nvert),
intent(in) :: vertices
979 integer(I4B),
dimension(:),
intent(in) :: iavert
980 integer(I4B),
dimension(:),
intent(in) :: javert
981 real(DP),
dimension(3, nodesuser),
intent(in) :: cellxyz
983 real(DP),
dimension(nodesuser),
intent(in) :: cellfdc
984 integer(I4B),
dimension(:),
intent(in) :: nodereduced
985 integer(I4B),
dimension(:),
intent(in) :: nodeuser
986 real(DP),
dimension(:),
intent(in) :: reach_length
988 integer(I4B),
dimension(:),
allocatable :: itemp
989 integer(I4B),
dimension(:),
allocatable :: iavertcells
990 integer(I4B),
dimension(:),
allocatable :: javertcells
992 integer(I4B) :: n, m, i, j, ierror
995 call this%allocate_scalars(name_model)
1004 allocate (itemp(nvert))
1008 call vertcellspm%init(nvert, nodes, itemp)
1011 do i = iavert(j), iavert(j + 1) - 1
1012 call vertcellspm%addconnection(javert(i), j, 1)
1015 call vertcellspm%sort()
1016 allocate (iavertcells(nvert + 1))
1017 allocate (javertcells(vertcellspm%nnz))
1018 call vertcellspm%filliaja(iavertcells, javertcells, ierror)
1019 call vertcellspm%destroy()
1023 iavertcells, javertcells, nodereduced)
1026 this%nja = sparse%nnz
1027 this%njas = (this%nja - this%nodes) / 2
1030 deallocate (iavertcells)
1031 deallocate (javertcells)
1034 call this%allocate_arrays()
1038 call sparse%filliaja(this%ia, this%ja, ierror)
1039 call sparse%destroy()
1042 call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
1043 call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
1048 this%ihc, this%cl1, this%cl2)
1064 call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
1074 integer(I4B),
dimension(:),
intent(in) :: ia
1075 integer(I4B),
dimension(:),
intent(in) :: ja
1076 integer(I4B),
dimension(:),
intent(in) :: jas
1077 real(DP),
dimension(:),
intent(in) :: reach_length
1078 integer(I4B),
dimension(:),
intent(out) :: ihc
1079 real(DP),
dimension(:),
intent(out) :: cl1
1080 real(DP),
dimension(:),
intent(out) :: cl2
1084 integer(I4B) :: ipos
1085 integer(I4B) :: isympos
1088 do n = 1,
size(reach_length)
1089 do ipos = ia(n) + 1, ia(n + 1) - 1
1094 cl1(isympos) =
dhalf * reach_length(n)
1095 cl2(isympos) =
dhalf * reach_length(m)
1103 subroutine iajausr(this, nrsize, nodesuser, nodereduced, nodeuser)
1108 integer(I4B),
intent(in) :: nrsize
1109 integer(I4B),
intent(in) :: nodesuser
1110 integer(I4B),
dimension(:),
intent(in) :: nodereduced
1111 integer(I4B),
dimension(:),
intent(in) :: nodeuser
1113 integer(I4B) :: n, nr, ipos
1117 if (nrsize > 0)
then
1121 call mem_reallocate(this%iausr, nodesuser + 1,
'IAUSR', this%memoryPath)
1122 this%iausr(nodesuser + 1) = this%ia(this%nodes + 1)
1123 do n = nodesuser, 1, -1
1126 this%iausr(n) = this%iausr(n + 1)
1128 this%iausr(n) = this%ia(nr)
1134 call mem_reallocate(this%jausr, this%nja,
'JAUSR', this%memoryPath)
1135 do ipos = 1, this%nja
1138 this%jausr(ipos) = n
1144 call mem_setptr(this%iausr,
'IA', this%memoryPath)
1145 call mem_setptr(this%jausr,
'JA', this%memoryPath)
1167 integer(I4B),
intent(in) :: node1, node2
1172 if (node1 < 1 .or. node1 > this%nodes .or. node2 < 1 .or. &
1173 node2 > this%nodes)
then
1179 if (node1 == node2)
then
1185 do i = this%ia(node1) + 1, this%ia(node1 + 1) - 1
1186 if (this%ja(i) == node2)
then
1202 integer(I4B),
intent(in) :: neq
1203 integer(I4B),
intent(in) :: nja
1204 integer(I4B),
intent(inout),
dimension(nja) :: isym
1206 integer(I4B),
intent(in),
dimension(neq + 1) :: ia
1207 integer(I4B),
intent(in),
dimension(nja) :: ja
1208 integer(I4B) :: n, m, ii, jj
1211 do ii = ia(n), ia(n + 1) - 1
1215 search:
do jj = ia(m), ia(m + 1) - 1
1216 if (ja(jj) == n)
then
1235 integer(I4B),
intent(in) :: neq
1236 integer(I4B),
intent(in) :: nja
1237 integer(I4B),
intent(in),
dimension(neq + 1) :: ia
1238 integer(I4B),
intent(in),
dimension(nja) :: ja
1239 integer(I4B),
intent(in),
dimension(nja) :: isym
1240 integer(I4B),
intent(inout),
dimension(nja) :: jas
1242 integer(I4B) :: n, m, ii, ipos
1248 do ii = ia(n) + 1, ia(n + 1) - 1
1259 do ii = ia(n), ia(n + 1) - 1
1262 jas(ii) = jas(isym(ii))
1274 vertcellspm, cell1, cell2, nodereduced)
1279 integer(I4B),
intent(in) :: nodes
1280 integer(I4B),
intent(in) :: nrsize
1281 integer(I4B),
intent(in) :: maxnnz
1282 integer(I4B),
intent(in) :: nlay
1283 integer(I4B),
intent(in) :: ncpl
1286 integer(I4B),
dimension(:),
intent(in) :: nodereduced
1289 integer(I4B),
dimension(:),
allocatable :: rowmaxnnz
1290 integer(I4B) :: i, j, k, kk, nr, mr, j1, j2, icol1, icol2, nvert
1293 allocate (rowmaxnnz(nodes))
1295 rowmaxnnz(i) = maxnnz
1297 call sparse%init(nodes, nodes, rowmaxnnz)
1298 deallocate (rowmaxnnz)
1304 nr =
get_node(k, 1, j, nlay, 1, ncpl)
1305 if (nrsize > 0) nr = nodereduced(nr)
1309 call sparse%addconnection(nr, nr, 1)
1313 do kk = k - 1, 1, -1
1314 mr =
get_node(kk, 1, j, nlay, 1, ncpl)
1315 if (nrsize > 0) mr = nodereduced(mr)
1319 call sparse%addconnection(nr, mr, 1)
1326 mr =
get_node(kk, 1, j, nlay, 1, ncpl)
1327 if (nrsize > 0) mr = nodereduced(mr)
1331 call sparse%addconnection(nr, mr, 1)
1339 nvert = vertcellspm%nrow
1341 do icol1 = 1, vertcellspm%row(i)%nnz
1342 j1 = vertcellspm%row(i)%icolarray(icol1)
1344 nr =
get_node(k, 1, j1, nlay, 1, ncpl)
1345 if (nrsize > 0) nr = nodereduced(nr)
1347 call cell1%set_nodered(nr)
1348 do icol2 = 1, vertcellspm%row(i)%nnz
1349 j2 = vertcellspm%row(i)%icolarray(icol2)
1351 mr =
get_node(k, 1, j2, nlay, 1, ncpl)
1352 if (nrsize > 0) mr = nodereduced(mr)
1354 call cell2%set_nodered(mr)
1355 if (cell1%shares_edge(cell2))
then
1356 call sparse%addconnection(nr, mr, 1)
1371 iavertcells, javertcells, &
1377 integer(I4B),
intent(in) :: nodes
1378 integer(I4B),
intent(in) :: nrsize
1379 integer(I4B),
intent(in) :: maxnnz
1380 integer(I4B),
intent(in) :: nodeuser
1382 integer(I4B),
dimension(:),
intent(in) :: nodereduced
1383 integer(I4B),
dimension(:),
intent(in) :: iavertcells
1384 integer(I4B),
dimension(:),
intent(in) :: javertcells
1386 integer(I4B),
dimension(:),
allocatable :: rowmaxnnz
1387 integer(I4B) :: i, k, nr, mr, nvert
1391 allocate (rowmaxnnz(nodes))
1393 rowmaxnnz(i) = maxnnz
1395 call sparse%init(nodes, nodes, rowmaxnnz)
1396 deallocate (rowmaxnnz)
1401 if (nrsize > 0) mr = nodereduced(mr)
1403 call sparse%addconnection(mr, mr, 1)
1408 nvert =
size(iavertcells) - 1
1412 do k = iavertcells(i), iavertcells(i + 1) - 2
1414 do con = k + 1, iavertcells(i + 1) - 1
1416 if (nrsize > 0) nr = nodereduced(nr)
1418 mr = javertcells(con)
1419 if (nrsize > 0) mr = nodereduced(mr)
1421 call sparse%addconnection(nr, mr, 1)
1422 call sparse%addconnection(mr, nr, 1)
1439 integer(I4B),
intent(in) :: ipos
1440 integer(I4B),
intent(in) :: maskval
1445 if (
associated(this%mask, this%ja))
then
1446 call mem_allocate(this%mask, this%nja,
'MASK', this%memoryPath)
1454 this%mask(ipos) = maskval
1464 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: iac
1465 integer(I4B),
dimension(:),
contiguous,
intent(inout) :: ia
1467 integer(I4B) :: n, nodes
1473 if (n <
size(ia))
then
1474 ia(n) = iac(n) + ia(n - 1)
1476 ia(n) = ia(n) + ia(n - 1)
1479 do n = nodes + 1, 2, -1
1480 ia(n) = ia(n - 1) + 1
This module contains block parser methods.
subroutine disvconnections(this, name_model, nodes, ncpl, nlay, nrsize, nvert, vertex, iavert, javert, cellxy, top, bot, nodereduced, nodeuser)
Construct the connectivity arrays using cell disv information.
subroutine iajausr(this, nrsize, nodesuser, nodereduced, nodeuser)
Fill iausr and jausr if reduced grid, otherwise point them to ia and ja.
subroutine set_mask(this, ipos, maskval)
routine to set a value in the mask array (which has the same shape as thisja)
subroutine read_connectivity_from_block(this, name_model, nodes, nja, iout)
Read and process IAC and JA from an an input block called CONNECTIONDATA.
subroutine fill_disv1d_symarrays(ia, ja, jas, reach_length, ihc, cl1, cl2)
Fill symmetric connection arrays for disv1d.
subroutine, public iac_to_ia(iac, ia)
Convert an iac array into an ia array.
subroutine con_da(this)
Deallocate connection variables.
subroutine set_cl1_cl2_from_fleng(this, fleng)
Using a vector of cell lengths, calculate the cl1 and cl2 arrays.
subroutine vertexconnect(nodes, nrsize, maxnnz, nlay, ncpl, sparse, vertcellspm, cell1, cell2, nodereduced)
Routine to make cell connections from vertices.
subroutine, public fillisym(neq, nja, ia, ja, isym)
Function to fill the isym array.
subroutine allocate_arrays(this)
Allocate arrays for ConnectionsType.
integer(i4b) function getjaindex(this, node1, node2)
Get the index in the JA array corresponding to the connection between two nodes of interest.
subroutine disv1dconnections_verts(this, name_model, nodes, nodesuser, nrsize, nvert, vertices, iavert, javert, cellxyz, cellfdc, nodereduced, nodeuser, reach_length)
Fill the connections object for a disv1d package from vertices.
subroutine disconnections(this, name_model, nodes, ncol, nrow, nlay, nrsize, delr, delc, top, bot, nodereduced, nodeuser)
Construct the connectivity arrays for a structured three-dimensional grid.
subroutine, public filljas(neq, nja, ia, ja, isym, jas)
Function to fill the jas array.
subroutine vertexconnectl(nodes, nrsize, maxnnz, nodeuser, sparse, iavertcells, javertcells, nodereduced)
Routine to make cell connections from vertices for a linear network.
subroutine con_finalize(this, ihctemp, cl12temp, hwvatemp, angldegx)
Finalize connection data.
subroutine allocate_scalars(this, name_model)
Allocate scalars for ConnectionsType.
subroutine disuconnections(this, name_model, nodes, nodesuser, nrsize, nodereduced, nodeuser, iainp, jainp, ihcinp, cl12inp, hwvainp, angldegxinp, iangledegx)
Construct the connectivity arrays using disu information. Grid may be reduced.
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
real(dp), parameter dnodata
real no data constant
real(dp), parameter dhalf
real constant 1/2
real(dp), parameter dpi
real constant
real(dp), parameter dpio180
real constant
real(dp), parameter dzero
real constant zero
real(dp), parameter dtwo
real constant 2
integer(i4b), parameter lenmempath
maximum length of the memory path
real(dp), parameter dthree
real constant 3
real(dp), parameter done
real constant 1
integer(i4b) function, public get_node(ilay, irow, icol, nlay, nrow, ncol)
Get node number, given layer, row, and column indices for a structured grid. If any argument is inval...
This module defines variable data types.
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
Store and issue logging messages to output units.
subroutine, public write_message(text, iunit, fmt, skipbefore, skipafter, advance)
Write a message to an output unit.
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