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

Data Types

type  connectionstype
 

Functions/Subroutines

subroutine con_da (this)
 Deallocate connection variables. More...
 
subroutine allocate_scalars (this, name_model)
 Allocate scalars for ConnectionsType. More...
 
subroutine allocate_arrays (this)
 Allocate arrays for ConnectionsType. More...
 
subroutine con_finalize (this, ihctemp, cl12temp, hwvatemp, angldegx)
 Finalize connection data. More...
 
subroutine read_connectivity_from_block (this, name_model, nodes, nja, iout)
 Read and process IAC and JA from an an input block called CONNECTIONDATA. More...
 
subroutine set_cl1_cl2_from_fleng (this, fleng)
 Using a vector of cell lengths, calculate the cl1 and cl2 arrays. More...
 
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. More...
 
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. More...
 
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. More...
 
subroutine disv1dconnections_verts (this, name_model, nodes, nodesuser, nrsize, nvert, vertices, iavert, javert, cellxy, cellfdc, nodereduced, nodeuser, reach_length)
 Fill the connections object for a disv1d package from vertices. More...
 
subroutine fill_disv1d_symarrays (ia, ja, jas, cell_length, ihc, cl1, cl2, nrsize, nodereduced, nodeuser, fdc, iavert, javert, iavertcells, javertcells)
 Fill symmetric connection arrays for disv1d. More...
 
subroutine iajausr (this, nrsize, nodesuser, nodereduced, nodeuser)
 Fill iausr and jausr if reduced grid, otherwise point them to ia and ja. More...
 
integer(i4b) function getjaindex (this, node1, node2)
 Get the index in the JA array corresponding to the connection between two nodes of interest. More...
 
subroutine, public fillisym (neq, nja, ia, ja, isym)
 Function to fill the isym array. More...
 
subroutine, public filljas (neq, nja, ia, ja, isym, jas)
 Function to fill the jas array. More...
 
subroutine vertexconnect (nodes, nrsize, maxnnz, nlay, ncpl, sparse, vertcellspm, cell1, cell2, nodereduced)
 Routine to make cell connections from vertices. More...
 
subroutine vertexconnectl (nodes, nrsize, maxnnz, nodeuser, sparse, iavertcells, javertcells, nodereduced)
 Routine to make cell connections from vertices for a linear network. More...
 
subroutine set_mask (this, ipos, maskval)
 routine to set a value in the mask array (which has the same shape as thisja) More...
 
subroutine, public iac_to_ia (iac, ia)
 Convert an iac array into an ia array. More...
 
logical(lgp) function connected_down_n (nu, mu, iavert, javert, iavertcells, javertcells)
 Is cell m is connected to the downstream end of cell n. More...
 

Function/Subroutine Documentation

◆ allocate_arrays()

subroutine connectionsmodule::allocate_arrays ( class(connectionstype this)

Definition at line 135 of file Connections.f90.

137  ! -- dummy
138  class(ConnectionsType) :: this
139  !
140  ! -- allocate space for connection arrays
141  call mem_allocate(this%ia, this%nodes + 1, 'IA', this%memoryPath)
142  call mem_allocate(this%ja, this%nja, 'JA', this%memoryPath)
143  call mem_allocate(this%isym, this%nja, 'ISYM', this%memoryPath)
144  call mem_allocate(this%jas, this%nja, 'JAS', this%memoryPath)
145  call mem_allocate(this%hwva, this%njas, 'HWVA', this%memoryPath)
146  call mem_allocate(this%anglex, this%njas, 'ANGLEX', this%memoryPath)
147  call mem_allocate(this%ihc, this%njas, 'IHC', this%memoryPath)
148  call mem_allocate(this%cl1, this%njas, 'CL1', this%memoryPath)
149  call mem_allocate(this%cl2, this%njas, 'CL2', this%memoryPath)
150  call mem_allocate(this%iausr, 1, 'IAUSR', this%memoryPath)
151  call mem_allocate(this%jausr, 1, 'JAUSR', this%memoryPath)
152  !
153  ! -- let mask point to ja, which is always nonzero,
154  ! until someone decides to do a 'set_mask'
155  this%mask => this%ja

◆ allocate_scalars()

subroutine connectionsmodule::allocate_scalars ( class(connectionstype this,
character(len=*), intent(in)  name_model 
)

Definition at line 110 of file Connections.f90.

111  ! -- modules
114  ! -- dummy
115  class(ConnectionsType) :: this
116  character(len=*), intent(in) :: name_model
117  !
118  ! -- allocate
119  allocate (this%name_model)
120  !
121  this%memoryPath = create_mem_path(name_model, 'CON')
122  call mem_allocate(this%nodes, 'NODES', this%memoryPath)
123  call mem_allocate(this%nja, 'NJA', this%memoryPath)
124  call mem_allocate(this%njas, 'NJAS', this%memoryPath)
125  call mem_allocate(this%ianglex, 'IANGLEX', this%memoryPath)
126  this%name_model = name_model
127  this%nodes = 0
128  this%nja = 0
129  this%njas = 0
130  this%ianglex = 0
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
Here is the call graph for this function:

◆ con_da()

subroutine connectionsmodule::con_da ( class(connectionstype this)
private

Definition at line 61 of file Connections.f90.

62  ! -- modules
64  ! -- dummy
65  class(ConnectionsType) :: this
66  !
67  ! -- Strings
68  deallocate (this%name_model)
69  !
70  ! -- Scalars
71  call mem_deallocate(this%nodes)
72  call mem_deallocate(this%nja)
73  call mem_deallocate(this%njas)
74  call mem_deallocate(this%ianglex)
75  !
76  ! -- iausr and jausr
77  if (associated(this%iausr, this%ia)) then
78  nullify (this%iausr)
79  else
80  call mem_deallocate(this%iausr)
81  end if
82  if (associated(this%jausr, this%ja)) then
83  nullify (this%jausr)
84  else
85  call mem_deallocate(this%jausr)
86  end if
87  ! -- mask
88  if (associated(this%mask, this%ja)) then
89  nullify (this%mask)
90  else
91  call mem_deallocate(this%mask)
92  end if
93  !
94  ! -- Arrays
95  call mem_deallocate(this%ia)
96  if (size(this%ja) > 0) then
97  call mem_deallocate(this%ja)
98  call mem_deallocate(this%isym)
99  call mem_deallocate(this%jas)
100  call mem_deallocate(this%hwva)
101  call mem_deallocate(this%anglex)
102  call mem_deallocate(this%ihc)
103  call mem_deallocate(this%cl1)
104  call mem_deallocate(this%cl2)
105  end if

◆ con_finalize()

subroutine connectionsmodule::con_finalize ( class(connectionstype this,
integer(i4b), dimension(:), intent(in)  ihctemp,
real(dp), dimension(:), intent(in)  cl12temp,
real(dp), dimension(:), intent(in)  hwvatemp,
real(dp), dimension(:), intent(in)  angldegx 
)

Definition at line 160 of file Connections.f90.

161  ! -- modules
164  ! -- dummy
165  class(ConnectionsType) :: this
166  integer(I4B), dimension(:), intent(in) :: ihctemp
167  real(DP), dimension(:), intent(in) :: cl12temp
168  real(DP), dimension(:), intent(in) :: hwvatemp
169  real(DP), dimension(:), intent(in) :: angldegx
170  ! -- local
171  integer(I4B) :: ii, n, m
172  integer(I4B), parameter :: nname = 6
173  character(len=24), dimension(nname) :: aname(nname)
174  ! -- formats
175  character(len=*), parameter :: fmtsymerr = &
176  &"('Error in array: ',a,'.', &
177  &' Array is not symmetric in positions: ',i0,' and ',i0,'.', &
178  &' Values in these positions are: ',1pg15.6,' and ', 1pg15.6, &
179  &' For node ',i0,' connected to node ',i0)"
180  character(len=*), parameter :: fmtsymerrja = &
181  &"('Error in array: ',a,'.', &
182  &' Array does not have symmetric counterpart in position ',i0, &
183  &' for cell ',i0,' connected to cell ',i0)"
184  character(len=*), parameter :: fmtjanmerr = &
185  &"('Error in array: ',a,'.', &
186  &' First value for cell : ',i0,' must equal ',i0,'.', &
187  &' Found ',i0,' instead.')"
188  character(len=*), parameter :: fmtjasorterr = &
189  &"('Error in array: ',a,'.', &
190  &' Entries not sorted for row: ',i0,'.', &
191  &' Offending entries are: ',i0,' and ',i0)"
192  character(len=*), parameter :: fmtihcerr = &
193  "('IHC must be 0, 1, or 2. Found: ',i0)"
194  ! -- data
195  data aname(1)/' IAC'/
196  data aname(2)/' JA'/
197  data aname(3)/' IHC'/
198  data aname(4)/' CL12'/
199  data aname(5)/' HWVA'/
200  data aname(6)/' ANGLDEGX'/
201  !
202  ! -- Convert any negative ja numbers to positive
203  do ii = 1, this%nja
204  if (this%ja(ii) < 0) this%ja(ii) = -this%ja(ii)
205  end do
206  !
207  ! -- Ensure ja is sorted with the row column listed first
208  do n = 1, this%nodes
209  m = this%ja(this%ia(n))
210  if (n /= m) then
211  write (errmsg, fmtjanmerr) trim(adjustl(aname(2))), n, n, m
212  call store_error(errmsg)
213  end if
214  do ii = this%ia(n) + 1, this%ia(n + 1) - 2
215  m = this%ja(ii)
216  if (m > this%ja(ii + 1)) then
217  write (errmsg, fmtjasorterr) trim(adjustl(aname(2))), n, &
218  m, this%ja(ii + 1)
219  call store_error(errmsg)
220  end if
221  end do
222  end do
223  if (count_errors() > 0) then
224  call this%parser%StoreErrorUnit()
225  end if
226  !
227  ! -- fill the isym arrays
228  call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
229  !
230  ! -- check for symmetry in ja (isym value of zero indicates there is no
231  ! symmetric connection
232  do n = 1, this%nodes
233  do ii = this%ia(n), this%ia(n + 1) - 1
234  m = this%ja(ii)
235  if (this%isym(ii) == 0) then
236  write (errmsg, fmtsymerrja) trim(adjustl(aname(2))), ii, n, m
237  call store_error(errmsg)
238  end if
239  end do
240  end do
241  if (count_errors() > 0) then
242  call this%parser%StoreErrorUnit()
243  end if
244  !
245  ! -- Fill the jas array, which maps any connection to upper triangle
246  call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
247  !
248  ! -- Put into symmetric array
249  do n = 1, this%nodes
250  do ii = this%ia(n) + 1, this%ia(n + 1) - 1
251  m = this%ja(ii)
252  if (ihctemp(ii) /= ihctemp(this%isym(ii))) then
253  write (errmsg, fmtsymerr) trim(adjustl(aname(3))), ii, this%isym(ii), &
254  ihctemp(ii), ihctemp(this%isym(ii)), n, m
255  call store_error(errmsg)
256  else
257  this%ihc(this%jas(ii)) = ihctemp(ii)
258  end if
259  end do
260  end do
261  if (count_errors() > 0) then
262  call this%parser%StoreErrorUnit()
263  end if
264  !
265  ! -- Put cl12 into symmetric arrays cl1 and cl2
266  do n = 1, this%nodes
267  do ii = this%ia(n) + 1, this%ia(n + 1) - 1
268  m = this%ja(ii)
269  if (m > n) then
270  this%cl1(this%jas(ii)) = cl12temp(ii)
271  elseif (n > m) then
272  this%cl2(this%jas(ii)) = cl12temp(ii)
273  end if
274  end do
275  end do
276  !
277  ! -- Put HWVA into symmetric array based on the value of IHC
278  ! IHC = 0, vertical connection, HWVA is vertical flow area
279  ! IHC = 1, horizontal connection, HWVA is the width perpendicular to
280  ! flow
281  ! IHC = 2, horizontal connection for a vertically staggered grid.
282  ! HWVA is the width perpendicular to flow.
283  do n = 1, this%nodes
284  do ii = this%ia(n) + 1, this%ia(n + 1) - 1
285  m = this%ja(ii)
286  if (hwvatemp(ii) /= hwvatemp(this%isym(ii))) then
287  write (errmsg, fmtsymerr) trim(adjustl(aname(5))), ii, this%isym(ii), &
288  hwvatemp(ii), hwvatemp(this%isym(ii)), n, m
289  call store_error(errmsg)
290  end if
291  if (ihctemp(ii) < 0 .or. ihctemp(ii) > 2) then
292  write (errmsg, fmtihcerr) ihctemp(ii)
293  call store_error(errmsg)
294  end if
295  this%hwva(this%jas(ii)) = hwvatemp(ii)
296  end do
297  end do
298  if (count_errors() > 0) then
299  call this%parser%StoreErrorUnit()
300  end if
301  !
302  ! -- Put anglextemp into this%anglex; store only upper triangle
303  if (this%ianglex /= 0) then
304  do n = 1, this%nodes
305  do ii = this%ia(n) + 1, this%ia(n + 1) - 1
306  m = this%ja(ii)
307  if (n > m) cycle
308  this%anglex(this%jas(ii)) = angldegx(ii) * dpio180
309  end do
310  end do
311  else
312  do n = 1, size(this%anglex)
313  this%anglex(n) = dnodata
314  end do
315  end if
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
real(dp), parameter dnodata
real no data constant
Definition: Constants.f90:95
real(dp), parameter dhalf
real constant 1/2
Definition: Constants.f90:68
real(dp), parameter dpio180
real constant
Definition: Constants.f90:130
real(dp), parameter done
real constant 1
Definition: Constants.f90:76
This module contains simulation methods.
Definition: Sim.f90:10
subroutine, public store_error(msg, terminate)
Store an error message.
Definition: Sim.f90:92
integer(i4b) function, public count_errors()
Return number of errors.
Definition: Sim.f90:59
subroutine, public store_error_unit(iunit, terminate)
Store the file unit number.
Definition: Sim.f90:168
Here is the call graph for this function:

◆ connected_down_n()

logical(lgp) function connectionsmodule::connected_down_n ( integer(i4b), intent(in)  nu,
integer(i4b), intent(in)  mu,
integer(i4b), dimension(:), intent(in)  iavert,
integer(i4b), dimension(:), intent(in)  javert,
integer(i4b), dimension(:), intent(in)  iavertcells,
integer(i4b), dimension(:), intent(in)  javertcells 
)
private

Definition at line 1477 of file Connections.f90.

1479  ! dummy
1480  integer(I4B), intent(in) :: nu ! user nodenumber for cell n
1481  integer(I4B), intent(in) :: mu ! user nodenumber for connected cell m
1482  integer(I4B), dimension(:), intent(in) :: iavert ! csr index array of size (nodeuser + 1) for javert
1483  integer(I4B), dimension(:), intent(in) :: javert ! csr array containing vertex numbers that define each cell
1484  integer(I4B), dimension(:), intent(in) :: iavertcells ! csr index array of size (nvert + 1) for javert
1485  integer(I4B), dimension(:), intent(in) :: javertcells ! csr array containing cells numbers that referenced for a vertex
1486  ! return
1487  logical(LGP) :: connected_down
1488  ! - local
1489  integer(I4B) :: ipos
1490  integer(I4B) :: ivert_down
1491 
1492  ! ivert_down is the last vertex for node nu; the last vertex is considered
1493  ! to be the dowstream end of node nu
1494  ivert_down = javert(iavert(nu + 1) - 1)
1495 
1496  ! look through vertex ivert_down, and see if mu is in it; if so, then
1497  ! that means mu is connected to the downstream end of nu
1498  connected_down = .false.
1499  do ipos = iavertcells(ivert_down), iavertcells(ivert_down + 1) - 1
1500  if (javertcells(ipos) == mu) then
1501  connected_down = .true.
1502  exit
1503  end if
1504  end do
1505 
Here is the caller graph for this function:

◆ disconnections()

subroutine connectionsmodule::disconnections ( class(connectionstype this,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  ncol,
integer(i4b), intent(in)  nrow,
integer(i4b), intent(in)  nlay,
integer(i4b), intent(in)  nrsize,
real(dp), dimension(ncol), intent(in)  delr,
real(dp), dimension(nrow), intent(in)  delc,
real(dp), dimension(nodes), intent(in)  top,
real(dp), dimension(nodes), intent(in)  bot,
integer(i4b), dimension(:), intent(in), target  nodereduced,
integer(i4b), dimension(:), intent(in)  nodeuser 
)

Definition at line 463 of file Connections.f90.

466  ! -- modules
467  use constantsmodule, only: dhalf, dzero, dthree, dtwo, dpi
468  use sparsemodule, only: sparsematrix
469  ! -- dummy
470  class(ConnectionsType) :: this
471  character(len=*), intent(in) :: name_model
472  integer(I4B), intent(in) :: nodes
473  integer(I4B), intent(in) :: ncol
474  integer(I4B), intent(in) :: nrow
475  integer(I4B), intent(in) :: nlay
476  integer(I4B), intent(in) :: nrsize
477  real(DP), dimension(ncol), intent(in) :: delr
478  real(DP), dimension(nrow), intent(in) :: delc
479  real(DP), dimension(nodes), intent(in) :: top
480  real(DP), dimension(nodes), intent(in) :: bot
481  integer(I4B), dimension(:), target, intent(in) :: nodereduced
482  integer(I4B), dimension(:), intent(in) :: nodeuser
483  ! -- local
484  integer(I4B), dimension(:, :, :), pointer :: nrdcd_ptr => null() !non-contiguous because is a slice
485  integer(I4B), dimension(:), allocatable :: rowmaxnnz
486  type(sparsematrix) :: sparse
487  integer(I4B) :: i, j, k, kk, ierror, isympos, nodesuser
488  integer(I4B) :: nr, mr
489  !
490  ! -- Allocate scalars
491  call this%allocate_scalars(name_model)
492  !
493  ! -- Set scalars
494  this%nodes = nodes
495  this%ianglex = 1
496  !
497  ! -- Setup the sparse matrix object
498  allocate (rowmaxnnz(this%nodes))
499  do i = 1, this%nodes
500  rowmaxnnz(i) = 6
501  end do
502  call sparse%init(this%nodes, this%nodes, rowmaxnnz)
503  !
504  ! -- Create a 3d pointer to nodereduced for easier processing
505  if (nrsize /= 0) then
506  nrdcd_ptr(1:ncol, 1:nrow, 1:nlay) => nodereduced
507  end if
508  !
509  ! -- Add connections to sparse
510  do k = 1, nlay
511  do i = 1, nrow
512  do j = 1, ncol
513  !
514  ! -- Find the reduced node number and then cycle if the
515  ! node is always inactive
516  if (nrsize == 0) then
517  nr = get_node(k, i, j, nlay, nrow, ncol)
518  else
519  nr = nrdcd_ptr(j, i, k)
520  end if
521  if (nr <= 0) cycle
522  !
523  ! -- Process diagonal
524  call sparse%addconnection(nr, nr, 1)
525  !
526  ! -- Up direction
527  if (k > 1) then
528  do kk = k - 1, 1, -1
529  if (nrsize == 0) then
530  mr = get_node(kk, i, j, nlay, nrow, ncol)
531  else
532  mr = nrdcd_ptr(j, i, kk)
533  end if
534  if (mr >= 0) exit
535  end do
536  if (mr > 0) then
537  call sparse%addconnection(nr, mr, 1)
538  end if
539  end if
540  !
541  ! -- Back direction
542  if (i > 1) then
543  if (nrsize == 0) then
544  mr = get_node(k, i - 1, j, nlay, nrow, ncol)
545  else
546  mr = nrdcd_ptr(j, i - 1, k)
547  end if
548  if (mr > 0) then
549  call sparse%addconnection(nr, mr, 1)
550  end if
551  end if
552  !
553  ! -- Left direction
554  if (j > 1) then
555  if (nrsize == 0) then
556  mr = get_node(k, i, j - 1, nlay, nrow, ncol)
557  else
558  mr = nrdcd_ptr(j - 1, i, k)
559  end if
560  if (mr > 0) then
561  call sparse%addconnection(nr, mr, 1)
562  end if
563  end if
564  !
565  ! -- Right direction
566  if (j < ncol) then
567  if (nrsize == 0) then
568  mr = get_node(k, i, j + 1, nlay, nrow, ncol)
569  else
570  mr = nrdcd_ptr(j + 1, i, k)
571  end if
572  if (mr > 0) then
573  call sparse%addconnection(nr, mr, 1)
574  end if
575  end if
576  !
577  ! -- Front direction
578  if (i < nrow) then !front
579  if (nrsize == 0) then
580  mr = get_node(k, i + 1, j, nlay, nrow, ncol)
581  else
582  mr = nrdcd_ptr(j, i + 1, k)
583  end if
584  if (mr > 0) then
585  call sparse%addconnection(nr, mr, 1)
586  end if
587  end if
588  !
589  ! -- Down direction
590  if (k < nlay) then
591  do kk = k + 1, nlay
592  if (nrsize == 0) then
593  mr = get_node(kk, i, j, nlay, nrow, ncol)
594  else
595  mr = nrdcd_ptr(j, i, kk)
596  end if
597  if (mr >= 0) exit
598  end do
599  if (mr > 0) then
600  call sparse%addconnection(nr, mr, 1)
601  end if
602  end if
603  end do
604  end do
605  end do
606  this%nja = sparse%nnz
607  this%njas = (this%nja - this%nodes) / 2
608  !
609  ! -- Allocate index arrays of size nja and symmetric arrays
610  call this%allocate_arrays()
611  !
612  ! -- Fill the IA and JA arrays from sparse, then destroy sparse
613  call sparse%filliaja(this%ia, this%ja, ierror)
614  call sparse%destroy()
615  !
616  ! -- fill the isym and jas arrays
617  call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
618  call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
619  !
620  ! -- Fill symmetric discretization arrays (ihc,cl1,cl2,hwva,anglex)
621  isympos = 1
622  do k = 1, nlay
623  do i = 1, nrow
624  do j = 1, ncol
625  !
626  ! -- cycle if node is always inactive
627  if (nrsize == 0) then
628  nr = get_node(k, i, j, nlay, nrow, ncol)
629  else
630  nr = nrdcd_ptr(j, i, k)
631  end if
632  if (nr <= 0) cycle
633  !
634  ! -- right connection
635  if (j < ncol) then
636  if (nrsize == 0) then
637  mr = get_node(k, i, j + 1, nlay, nrow, ncol)
638  else
639  mr = nrdcd_ptr(j + 1, i, k)
640  end if
641  if (mr > 0) then
642  this%ihc(isympos) = 1
643  this%cl1(isympos) = dhalf * delr(j)
644  this%cl2(isympos) = dhalf * delr(j + 1)
645  this%hwva(isympos) = delc(i)
646  this%anglex(isympos) = dzero
647  isympos = isympos + 1
648  end if
649  end if
650  !
651  ! -- front connection
652  if (i < nrow) then
653  if (nrsize == 0) then
654  mr = get_node(k, i + 1, j, nlay, nrow, ncol)
655  else
656  mr = nrdcd_ptr(j, i + 1, k)
657  end if
658  if (mr > 0) then
659  this%ihc(isympos) = 1
660  this%cl1(isympos) = dhalf * delc(i)
661  this%cl2(isympos) = dhalf * delc(i + 1)
662  this%hwva(isympos) = delr(j)
663  this%anglex(isympos) = dthree / dtwo * dpi
664  isympos = isympos + 1
665  end if
666  end if
667  !
668  ! -- down connection
669  if (k < nlay) then
670  do kk = k + 1, nlay
671  if (nrsize == 0) then
672  mr = get_node(kk, i, j, nlay, nrow, ncol)
673  else
674  mr = nrdcd_ptr(j, i, kk)
675  end if
676  if (mr >= 0) exit
677  end do
678  if (mr > 0) then
679  this%ihc(isympos) = 0
680  this%cl1(isympos) = dhalf * (top(nr) - bot(nr))
681  this%cl2(isympos) = dhalf * (top(mr) - bot(mr))
682  this%hwva(isympos) = delr(j) * delc(i)
683  this%anglex(isympos) = dzero
684  isympos = isympos + 1
685  end if
686  end if
687  end do
688  end do
689  end do
690  !
691  ! -- Deallocate temporary arrays
692  deallocate (rowmaxnnz)
693  !
694  ! -- If reduced system, then need to build iausr and jausr, otherwise point
695  ! them to ia and ja.
696  nodesuser = nlay * nrow * ncol
697  call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
real(dp), parameter dpi
real constant
Definition: Constants.f90:128
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65
real(dp), parameter dtwo
real constant 2
Definition: Constants.f90:79
real(dp), parameter dthree
real constant 3
Definition: Constants.f90:80
Here is the call graph for this function:

◆ disuconnections()

subroutine connectionsmodule::disuconnections ( class(connectionstype this,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nodesuser,
integer(i4b), intent(in)  nrsize,
integer(i4b), dimension(:), intent(in), contiguous  nodereduced,
integer(i4b), dimension(:), intent(in), contiguous  nodeuser,
integer(i4b), dimension(:), intent(in), contiguous  iainp,
integer(i4b), dimension(:), intent(in), contiguous  jainp,
integer(i4b), dimension(:), intent(in), contiguous  ihcinp,
real(dp), dimension(:), intent(in), contiguous  cl12inp,
real(dp), dimension(:), intent(in), contiguous  hwvainp,
real(dp), dimension(:), intent(in), contiguous  angldegxinp,
integer(i4b), intent(in)  iangledegx 
)

Definition at line 801 of file Connections.f90.

805  ! -- modules
806  use constantsmodule, only: dhalf, dzero, dthree, dtwo, dpi
807  use sparsemodule, only: sparsematrix
809  ! -- dummy
810  class(ConnectionsType) :: this
811  character(len=*), intent(in) :: name_model
812  integer(I4B), intent(in) :: nodes
813  integer(I4B), intent(in) :: nodesuser
814  integer(I4B), intent(in) :: nrsize
815  integer(I4B), dimension(:), contiguous, intent(in) :: nodereduced
816  integer(I4B), dimension(:), contiguous, intent(in) :: nodeuser
817  integer(I4B), dimension(:), contiguous, intent(in) :: iainp
818  integer(I4B), dimension(:), contiguous, intent(in) :: jainp
819  integer(I4B), dimension(:), contiguous, intent(in) :: ihcinp
820  real(DP), dimension(:), contiguous, intent(in) :: cl12inp
821  real(DP), dimension(:), contiguous, intent(in) :: hwvainp
822  real(DP), dimension(:), contiguous, intent(in) :: angldegxinp
823  integer(I4B), intent(in) :: iangledegx
824  ! -- local
825  integer(I4B), dimension(:), allocatable :: ihctemp
826  real(DP), dimension(:), allocatable :: cl12temp
827  real(DP), dimension(:), allocatable :: hwvatemp
828  real(DP), dimension(:), allocatable :: angldegxtemp
829  integer(I4B) :: nr, nu, mr, mu, ipos, iposr, ierror
830  integer(I4B), dimension(:), allocatable :: rowmaxnnz
831  type(sparsematrix) :: sparse
832  !
833  ! -- Allocate scalars
834  call this%allocate_scalars(name_model)
835  !
836  ! -- Set scalars
837  this%nodes = nodes
838  this%ianglex = iangledegx
839  !
840  ! -- If not a reduced grid, then copy and finalize, otherwise more
841  ! processing is required
842  if (nrsize == 0) then
843  this%nodes = nodes
844  this%nja = size(jainp)
845  this%njas = (this%nja - this%nodes) / 2
846  call this%allocate_arrays()
847  do nu = 1, nodes + 1
848  this%ia(nu) = iainp(nu)
849  end do
850  do ipos = 1, this%nja
851  this%ja(ipos) = jainp(ipos)
852  end do
853  !
854  ! -- Call con_finalize to check inp arrays and push larger arrays
855  ! into compressed symmetric arrays
856  call this%con_finalize(ihcinp, cl12inp, hwvainp, angldegxinp)
857  !
858  else
859  ! -- reduced system requires more work
860  !
861  ! -- Setup the sparse matrix object
862  allocate (rowmaxnnz(this%nodes))
863  do nr = 1, this%nodes
864  nu = nodeuser(nr)
865  rowmaxnnz(nr) = iainp(nu + 1) - iainp(nu)
866  end do
867  call sparse%init(this%nodes, this%nodes, rowmaxnnz)
868  !
869  ! -- go through user connectivity and create sparse
870  do nu = 1, nodesuser
871  nr = nodereduced(nu)
872  if (nr > 0) call sparse%addconnection(nr, nr, 1)
873  do ipos = iainp(nu) + 1, iainp(nu + 1) - 1
874  mu = jainp(ipos)
875  mr = nodereduced(mu)
876  if (nr < 1) cycle
877  if (mr < 1) cycle
878  call sparse%addconnection(nr, mr, 1)
879  end do
880  end do
881  this%nja = sparse%nnz
882  this%njas = (this%nja - this%nodes) / 2
883  !
884  ! -- Allocate index arrays of size nja and symmetric arrays
885  call this%allocate_arrays()
886  !
887  ! -- Fill the IA and JA arrays from sparse, then destroy sparse
888  call sparse%sort()
889  call sparse%filliaja(this%ia, this%ja, ierror)
890  call sparse%destroy()
891  deallocate (rowmaxnnz)
892  !
893  ! -- At this point, need to reduce ihc, cl12, hwva, and angldegx
894  allocate (ihctemp(this%nja))
895  allocate (cl12temp(this%nja))
896  allocate (hwvatemp(this%nja))
897  allocate (angldegxtemp(this%nja))
898  !
899  ! -- Compress user arrays into reduced arrays
900  iposr = 1
901  do nu = 1, nodesuser
902  nr = nodereduced(nu)
903  do ipos = iainp(nu), iainp(nu + 1) - 1
904  mu = jainp(ipos)
905  mr = nodereduced(mu)
906  if (nr < 1 .or. mr < 1) cycle
907  ihctemp(iposr) = ihcinp(ipos)
908  cl12temp(iposr) = cl12inp(ipos)
909  hwvatemp(iposr) = hwvainp(ipos)
910  angldegxtemp(iposr) = angldegxinp(ipos)
911  iposr = iposr + 1
912  end do
913  end do
914  !
915  ! -- call finalize
916  call this%con_finalize(ihctemp, cl12temp, hwvatemp, angldegxtemp)
917  !
918  ! -- deallocate temporary arrays
919  deallocate (ihctemp)
920  deallocate (cl12temp)
921  deallocate (hwvatemp)
922  deallocate (angldegxtemp)
923  end if
924  !
925  ! -- If reduced system, then need to build iausr and jausr, otherwise point
926  ! them to ia and ja.
927  call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)

◆ disv1dconnections_verts()

subroutine connectionsmodule::disv1dconnections_verts ( class(connectionstype this,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nodesuser,
integer(i4b), intent(in)  nrsize,
integer(i4b), intent(in)  nvert,
real(dp), dimension(3, nvert), intent(in)  vertices,
integer(i4b), dimension(:), intent(in)  iavert,
integer(i4b), dimension(:), intent(in)  javert,
real(dp), dimension(2, nodesuser), intent(in)  cellxy,
real(dp), dimension(nodesuser), intent(in)  cellfdc,
integer(i4b), dimension(:), intent(in)  nodereduced,
integer(i4b), dimension(:), intent(in)  nodeuser,
real(dp), dimension(:), intent(in)  reach_length 
)

Note that nothing is done for hwva

Parameters
[in]reach_lengthlength of each reach

Definition at line 935 of file Connections.f90.

940  ! modules
941  use constantsmodule, only: dhalf, dzero, dthree, dtwo, dpi
942  use sparsemodule, only: sparsematrix
943  use geomutilmodule, only: get_node
945  ! dummy
946  class(ConnectionsType) :: this
947  character(len=*), intent(in) :: name_model
948  integer(I4B), intent(in) :: nodes
949  integer(I4B), intent(in) :: nodesuser
950  integer(I4B), intent(in) :: nrsize
951  integer(I4B), intent(in) :: nvert
952  real(DP), dimension(3, nvert), intent(in) :: vertices
953  integer(I4B), dimension(:), intent(in) :: iavert
954  integer(I4B), dimension(:), intent(in) :: javert
955  real(DP), dimension(2, nodesuser), intent(in) :: cellxy
956  real(DP), dimension(nodesuser), intent(in) :: cellfdc
957  integer(I4B), dimension(:), intent(in) :: nodereduced
958  integer(I4B), dimension(:), intent(in) :: nodeuser
959  real(DP), dimension(:), intent(in) :: reach_length !< length of each reach
960  ! local
961  integer(I4B), dimension(:), allocatable :: itemp
962  integer(I4B), dimension(:), allocatable :: iavertcells
963  integer(I4B), dimension(:), allocatable :: javertcells
964  type(sparsematrix) :: sparse, vertcellspm
965  integer(I4B) :: i, j, ierror
966 
967  ! Allocate scalars
968  call this%allocate_scalars(name_model)
969 
970  ! Set scalars
971  this%nodes = nodes
972  this%ianglex = 1
973 
974  ! Create a sparse matrix array with a row for each vertex. The columns
975  ! in the sparse matrix contains the cells that include that vertex.
976  ! This array will be used to determine cell connectivity.
977  allocate (itemp(nvert))
978  do i = 1, nvert
979  itemp(i) = 4
980  end do
981  call vertcellspm%init(nvert, nodesuser, itemp)
982  deallocate (itemp)
983  do j = 1, nodesuser
984  do i = iavert(j), iavert(j + 1) - 1
985  call vertcellspm%addconnection(javert(i), j, 1)
986  end do
987  end do
988  call vertcellspm%sort()
989  allocate (iavertcells(nvert + 1))
990  allocate (javertcells(vertcellspm%nnz))
991  call vertcellspm%filliaja(iavertcells, javertcells, ierror)
992  call vertcellspm%destroy()
993 
994  ! Call routine to build a sparse matrix of the connections
995  call vertexconnectl(this%nodes, nrsize, 6, nodesuser, sparse, &
996  iavertcells, javertcells, nodereduced)
997  this%nja = sparse%nnz
998  this%njas = (this%nja - this%nodes) / 2
999 
1000  ! Allocate index arrays of size nja and symmetric arrays
1001  call this%allocate_arrays()
1002 
1003  ! Fill the IA and JA arrays from sparse, then destroy sparse
1004  call sparse%sort()
1005  call sparse%filliaja(this%ia, this%ja, ierror)
1006  call sparse%destroy()
1007 
1008  ! fill the isym and jas arrays
1009  call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
1010  call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
1011 
1012  ! fill the ihc, cl1, and cl2 arrays
1013  call fill_disv1d_symarrays(this%ia, this%ja, this%jas, reach_length, &
1014  this%ihc, this%cl1, this%cl2, &
1015  nrsize, nodereduced, nodeuser, cellfdc, &
1016  iavert, javert, iavertcells, javertcells)
1017 
1018  ! cleanup memory
1019  deallocate (iavertcells)
1020  deallocate (javertcells)
1021 
1022  ! If reduced system, then need to build iausr and jausr, otherwise point
1023  ! them to ia and ja.
1024  call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
1025 
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...
Definition: GeomUtil.f90:83
Here is the call graph for this function:

◆ disvconnections()

subroutine connectionsmodule::disvconnections ( class(connectionstype this,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  ncpl,
integer(i4b), intent(in)  nlay,
integer(i4b), intent(in)  nrsize,
integer(i4b), intent(in)  nvert,
real(dp), dimension(2, nvert), intent(in)  vertex,
integer(i4b), dimension(:), intent(in)  iavert,
integer(i4b), dimension(:), intent(in)  javert,
real(dp), dimension(2, ncpl), intent(in)  cellxy,
real(dp), dimension(nodes), intent(in)  top,
real(dp), dimension(nodes), intent(in)  bot,
integer(i4b), dimension(:), intent(in)  nodereduced,
integer(i4b), dimension(:), intent(in)  nodeuser 
)

Definition at line 702 of file Connections.f90.

705  ! -- modules
706  use constantsmodule, only: dhalf, dzero, dthree, dtwo, dpi
707  use sparsemodule, only: sparsematrix
708  use disvgeom, only: disvgeomtype
710  ! -- dummy
711  class(ConnectionsType) :: this
712  character(len=*), intent(in) :: name_model
713  integer(I4B), intent(in) :: nodes
714  integer(I4B), intent(in) :: ncpl
715  integer(I4B), intent(in) :: nlay
716  integer(I4B), intent(in) :: nrsize
717  integer(I4B), intent(in) :: nvert
718  real(DP), dimension(2, nvert), intent(in) :: vertex
719  integer(I4B), dimension(:), intent(in) :: iavert
720  integer(I4B), dimension(:), intent(in) :: javert
721  real(DP), dimension(2, ncpl), intent(in) :: cellxy
722  real(DP), dimension(nodes), intent(in) :: top
723  real(DP), dimension(nodes), intent(in) :: bot
724  integer(I4B), dimension(:), intent(in) :: nodereduced
725  integer(I4B), dimension(:), intent(in) :: nodeuser
726  ! -- local
727  integer(I4B), dimension(:), allocatable :: itemp
728  type(sparsematrix) :: sparse, vertcellspm
729  integer(I4B) :: n, m, ipos, i, j, ierror, nodesuser
730  type(DisvGeomType) :: cell1, cell2
731  !
732  ! -- Allocate scalars
733  call this%allocate_scalars(name_model)
734  !
735  ! -- Set scalars
736  this%nodes = nodes
737  this%ianglex = 1
738  !
739  ! -- Initialize DisvGeomType objects
740  call cell1%init(nlay, ncpl, nodes, top, bot, iavert, javert, vertex, &
741  cellxy, nodereduced, nodeuser)
742  call cell2%init(nlay, ncpl, nodes, top, bot, iavert, javert, vertex, &
743  cellxy, nodereduced, nodeuser)
744  !
745  ! -- Create a sparse matrix array with a row for each vertex. The columns
746  ! in the sparse matrix contains the cells that include that vertex.
747  ! This array will be used to determine horizontal cell connectivity.
748  allocate (itemp(nvert))
749  do i = 1, nvert
750  itemp(i) = 4
751  end do
752  call vertcellspm%init(nvert, ncpl, itemp)
753  deallocate (itemp)
754  do j = 1, ncpl
755  do i = iavert(j), iavert(j + 1) - 1
756  call vertcellspm%addconnection(javert(i), j, 1)
757  end do
758  end do
759  !
760  ! -- Call routine to build a sparse matrix of the connections
761  call vertexconnect(this%nodes, nrsize, 6, nlay, ncpl, sparse, &
762  vertcellspm, cell1, cell2, nodereduced)
763  this%nja = sparse%nnz
764  this%njas = (this%nja - this%nodes) / 2
765  !
766  ! -- Allocate index arrays of size nja and symmetric arrays
767  call this%allocate_arrays()
768  !
769  ! -- Fill the IA and JA arrays from sparse, then destroy sparse
770  call sparse%sort()
771  call sparse%filliaja(this%ia, this%ja, ierror)
772  call sparse%destroy()
773  !
774  ! -- fill the isym and jas arrays
775  call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
776  call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, this%jas)
777  !
778  ! -- Fill symmetric discretization arrays (ihc,cl1,cl2,hwva,anglex)
779  do n = 1, this%nodes
780  call cell1%set_nodered(n)
781  do ipos = this%ia(n) + 1, this%ia(n + 1) - 1
782  m = this%ja(ipos)
783  if (m < n) cycle
784  call cell2%set_nodered(m)
785  call cell1%cprops(cell2, this%hwva(this%jas(ipos)), &
786  this%cl1(this%jas(ipos)), this%cl2(this%jas(ipos)), &
787  this%anglex(this%jas(ipos)), &
788  this%ihc(this%jas(ipos)))
789  end do
790  end do
791  !
792  ! -- If reduced system, then need to build iausr and jausr, otherwise point
793  ! them to ia and ja.
794  nodesuser = nlay * ncpl
795  call this%iajausr(nrsize, nodesuser, nodereduced, nodeuser)
Here is the call graph for this function:

◆ fill_disv1d_symarrays()

subroutine connectionsmodule::fill_disv1d_symarrays ( integer(i4b), dimension(:), intent(in)  ia,
integer(i4b), dimension(:), intent(in)  ja,
integer(i4b), dimension(:), intent(in)  jas,
real(dp), dimension(:), intent(in)  cell_length,
integer(i4b), dimension(:), intent(out)  ihc,
real(dp), dimension(:), intent(out)  cl1,
real(dp), dimension(:), intent(out)  cl2,
integer(i4b), intent(in)  nrsize,
integer(i4b), dimension(:), intent(in)  nodereduced,
integer(i4b), dimension(:), intent(in)  nodeuser,
real(dp), dimension(:), intent(in)  fdc,
integer(i4b), dimension(:), intent(in)  iavert,
integer(i4b), dimension(:), intent(in)  javert,
integer(i4b), dimension(:), intent(in)  iavertcells,
integer(i4b), dimension(:), intent(in)  javertcells 
)
Parameters
[in]iacsr pointer array
[in]jacsr array
[in]jascsr symmetric array
[in]cell_lengthlength of each cell (all cells)
[out]ihchorizontal connection flag
[out]cl1distance from n to shared face with m
[out]cl2distance from m to shared face with n
[in]nrsizegreat than zero indicated reduced nodes present
[in]nodereducedmap user node to reduced node number
[in]nodeusermap user reduced node to user node number
[in]fdcfractional distance along cell to reach cell center

Definition at line 1030 of file Connections.f90.

1033  ! dummy
1034  integer(I4B), dimension(:), intent(in) :: ia !< csr pointer array
1035  integer(I4B), dimension(:), intent(in) :: ja !< csr array
1036  integer(I4B), dimension(:), intent(in) :: jas !< csr symmetric array
1037  real(DP), dimension(:), intent(in) :: cell_length !< length of each cell (all cells)
1038  integer(I4B), dimension(:), intent(out) :: ihc !< horizontal connection flag
1039  real(DP), dimension(:), intent(out) :: cl1 !< distance from n to shared face with m
1040  real(DP), dimension(:), intent(out) :: cl2 !< distance from m to shared face with n
1041  integer(I4B), intent(in) :: nrsize !< great than zero indicated reduced nodes present
1042  integer(I4B), dimension(:), intent(in) :: nodereduced !< map user node to reduced node number
1043  integer(I4B), dimension(:), intent(in) :: nodeuser !< map user reduced node to user node number
1044  real(DP), dimension(:), intent(in) :: fdc !< fractional distance along cell to reach cell center
1045  integer(I4B), dimension(:), intent(in) :: iavert ! csr index array of size (nodeuser + 1) for javert
1046  integer(I4B), dimension(:), intent(in) :: javert ! csr array containing vertex numbers that define each cell
1047  integer(I4B), dimension(:), intent(in) :: iavertcells ! csr index array of size (nvert + 1) for javert
1048  integer(I4B), dimension(:), intent(in) :: javertcells ! csr array containing cells numbers that referenced for a vertex
1049 
1050  ! local
1051  integer(I4B) :: nr, nu
1052  integer(I4B) :: mr, mu
1053  integer(I4B) :: ipos
1054  integer(I4B) :: isympos
1055  real(DP) :: f
1056 
1057  ! loop through and set array values
1058  do nu = 1, size(cell_length)
1059 
1060  ! find reduced node number and cycle if cell does not exist
1061  nr = nu
1062  if (nrsize > 0) nr = nodereduced(nu)
1063  if (nr <= 0) cycle
1064 
1065  ! visit each cell connected to reduced cell nr
1066  do ipos = ia(nr) + 1, ia(nr + 1) - 1
1067 
1068  ! process upper triangle by skipping mr cells less than nr
1069  mr = ja(ipos)
1070  if (mr < nr) cycle
1071 
1072  mu = mr
1073  if (nrsize > 0) mu = nodeuser(mr)
1074 
1075  isympos = jas(ipos)
1076  ihc(isympos) = 1
1077 
1078  ! if cell m is connected to the downstream end of cell n, then use
1079  ! 1 - fdc times the cell length, otherwise use fdc * length
1080  if (fdc(nu) == dhalf) then
1081  f = dhalf
1082  else
1083  if (connected_down_n(nu, mu, iavert, javert, iavertcells, &
1084  javertcells)) then
1085  f = (done - fdc(nu))
1086  else
1087  f = fdc(nu)
1088  end if
1089  end if
1090  cl1(isympos) = f * cell_length(nu)
1091 
1092  ! do the opposite for the cl2 distance as it is relative to m
1093  if (fdc(mu) == dhalf) then
1094  f = dhalf
1095  else
1096  if (connected_down_n(mu, nu, iavert, javert, iavertcells, &
1097  javertcells)) then
1098  f = (done - fdc(mu))
1099  else
1100  f = fdc(mu)
1101  end if
1102  end if
1103  cl2(isympos) = f * cell_length(mu)
1104 
1105  end do
1106  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ fillisym()

subroutine, public connectionsmodule::fillisym ( integer(i4b), intent(in)  neq,
integer(i4b), intent(in)  nja,
integer(i4b), dimension(neq + 1), intent(in)  ia,
integer(i4b), dimension(nja), intent(in)  ja,
integer(i4b), dimension(nja), intent(inout)  isym 
)

Definition at line 1205 of file Connections.f90.

1206  ! -- dummy
1207  integer(I4B), intent(in) :: neq
1208  integer(I4B), intent(in) :: nja
1209  integer(I4B), intent(inout), dimension(nja) :: isym
1210  ! -- local
1211  integer(I4B), intent(in), dimension(neq + 1) :: ia
1212  integer(I4B), intent(in), dimension(nja) :: ja
1213  integer(I4B) :: n, m, ii, jj
1214  !
1215  do n = 1, neq
1216  do ii = ia(n), ia(n + 1) - 1
1217  m = ja(ii)
1218  if (m /= n) then
1219  isym(ii) = 0
1220  search: do jj = ia(m), ia(m + 1) - 1
1221  if (ja(jj) == n) then
1222  isym(ii) = jj
1223  exit search
1224  end if
1225  end do search
1226  else
1227  isym(ii) = ii
1228  end if
1229  end do
1230  end do
Here is the caller graph for this function:

◆ filljas()

subroutine, public connectionsmodule::filljas ( integer(i4b), intent(in)  neq,
integer(i4b), intent(in)  nja,
integer(i4b), dimension(neq + 1), intent(in)  ia,
integer(i4b), dimension(nja), intent(in)  ja,
integer(i4b), dimension(nja), intent(in)  isym,
integer(i4b), dimension(nja), intent(inout)  jas 
)

Definition at line 1235 of file Connections.f90.

1236  ! -- dummy
1237  integer(I4B), intent(in) :: neq
1238  integer(I4B), intent(in) :: nja
1239  integer(I4B), intent(in), dimension(neq + 1) :: ia
1240  integer(I4B), intent(in), dimension(nja) :: ja
1241  integer(I4B), intent(in), dimension(nja) :: isym
1242  integer(I4B), intent(inout), dimension(nja) :: jas
1243  ! -- local
1244  integer(I4B) :: n, m, ii, ipos
1245  !
1246  ! -- set diagonal to zero and fill upper
1247  ipos = 1
1248  do n = 1, neq
1249  jas(ia(n)) = 0
1250  do ii = ia(n) + 1, ia(n + 1) - 1
1251  m = ja(ii)
1252  if (m > n) then
1253  jas(ii) = ipos
1254  ipos = ipos + 1
1255  end if
1256  end do
1257  end do
1258  !
1259  ! -- fill lower
1260  do n = 1, neq
1261  do ii = ia(n), ia(n + 1) - 1
1262  m = ja(ii)
1263  if (m < n) then
1264  jas(ii) = jas(isym(ii))
1265  end if
1266  end do
1267  end do
Here is the caller graph for this function:

◆ getjaindex()

integer(i4b) function connectionsmodule::getjaindex ( class(connectionstype this,
integer(i4b), intent(in)  node1,
integer(i4b), intent(in)  node2 
)

Node1 is used as the index in the IA array, and IA(Node1) is the row index in the (nodes x nodes) matrix represented by the compressed sparse row format.

-1 is returned if either node number is invalid. 0 is returned if the two nodes are not connected.

Definition at line 1168 of file Connections.f90.

1169  ! -- return
1170  integer(I4B) :: getjaindex
1171  ! -- dummy
1172  class(ConnectionsType) :: this
1173  integer(I4B), intent(in) :: node1, node2 ! nodes of interest
1174  ! -- local
1175  integer(I4B) :: i
1176  !
1177  ! -- error checking
1178  if (node1 < 1 .or. node1 > this%nodes .or. node2 < 1 .or. &
1179  node2 > this%nodes) then
1180  getjaindex = -1 ! indicates error (an invalid node number)
1181  return
1182  end if
1183  !
1184  ! -- If node1==node2, just return the position for the diagonal.
1185  if (node1 == node2) then
1186  getjaindex = this%ia(node1)
1187  return
1188  end if
1189  !
1190  ! -- Look for connection among nonzero elements defined for row node1.
1191  do i = this%ia(node1) + 1, this%ia(node1 + 1) - 1
1192  if (this%ja(i) == node2) then
1193  getjaindex = i
1194  return
1195  end if
1196  end do
1197  !
1198  ! -- If execution reaches here, no connection exists
1199  ! between nodes of interest.
1200  getjaindex = 0 ! indicates no connection exists

◆ iac_to_ia()

subroutine, public connectionsmodule::iac_to_ia ( integer(i4b), dimension(:), intent(in), pointer, contiguous  iac,
integer(i4b), dimension(:), intent(inout), contiguous  ia 
)

Definition at line 1452 of file Connections.f90.

1453  ! -- dummy
1454  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: iac
1455  integer(I4B), dimension(:), contiguous, intent(inout) :: ia
1456  ! -- local
1457  integer(I4B) :: n, nodes
1458  !
1459  ! -- Convert iac to ia
1460  nodes = size(iac)
1461  ia(1) = iac(1)
1462  do n = 2, size(ia) ! size(ia) == size(iac) + 1
1463  if (n < size(ia)) then
1464  ia(n) = iac(n) + ia(n - 1)
1465  else
1466  ia(n) = ia(n) + ia(n - 1)
1467  end if
1468  end do
1469  do n = nodes + 1, 2, -1
1470  ia(n) = ia(n - 1) + 1
1471  end do
1472  ia(1) = 1
Here is the caller graph for this function:

◆ iajausr()

subroutine connectionsmodule::iajausr ( class(connectionstype this,
integer(i4b), intent(in)  nrsize,
integer(i4b), intent(in)  nodesuser,
integer(i4b), dimension(:), intent(in)  nodereduced,
integer(i4b), dimension(:), intent(in)  nodeuser 
)
private

Definition at line 1112 of file Connections.f90.

1113  ! -- modules
1115  ! -- dummy
1116  class(ConnectionsType) :: this
1117  integer(I4B), intent(in) :: nrsize
1118  integer(I4B), intent(in) :: nodesuser
1119  integer(I4B), dimension(:), intent(in) :: nodereduced
1120  integer(I4B), dimension(:), intent(in) :: nodeuser
1121  ! -- local
1122  integer(I4B) :: n, nr, ipos
1123  !
1124  ! -- If reduced system, then need to build iausr and jausr, otherwise point
1125  ! them to ia and ja.
1126  if (nrsize > 0) then
1127  !
1128  ! -- Create the iausr array of size nodesuser + 1. For excluded cells,
1129  ! iausr(n) and iausr(n + 1) should be equal to indicate no connections.
1130  call mem_reallocate(this%iausr, nodesuser + 1, 'IAUSR', this%memoryPath)
1131  this%iausr(nodesuser + 1) = this%ia(this%nodes + 1)
1132  do n = nodesuser, 1, -1
1133  nr = nodereduced(n)
1134  if (nr < 1) then
1135  this%iausr(n) = this%iausr(n + 1)
1136  else
1137  this%iausr(n) = this%ia(nr)
1138  end if
1139  end do
1140  !
1141  ! -- Create the jausr array, which is the same size as ja, but it
1142  ! contains user node numbers instead of reduced node numbers
1143  call mem_reallocate(this%jausr, this%nja, 'JAUSR', this%memoryPath)
1144  do ipos = 1, this%nja
1145  nr = this%ja(ipos)
1146  n = nodeuser(nr)
1147  this%jausr(ipos) = n
1148  end do
1149  else
1150  ! -- iausr and jausr will be pointers
1151  call mem_deallocate(this%iausr)
1152  call mem_deallocate(this%jausr)
1153  call mem_setptr(this%iausr, 'IA', this%memoryPath)
1154  call mem_setptr(this%jausr, 'JA', this%memoryPath)
1155  end if

◆ read_connectivity_from_block()

subroutine connectionsmodule::read_connectivity_from_block ( class(connectionstype this,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nja,
integer(i4b), intent(in)  iout 
)

Definition at line 321 of file Connections.f90.

322  ! -- modules
323  use constantsmodule, only: linelength
325  ! -- dummy
326  class(ConnectionsType) :: this
327  character(len=*), intent(in) :: name_model
328  integer(I4B), intent(in) :: nodes
329  integer(I4B), intent(in) :: nja
330  integer(I4B), intent(in) :: iout
331  ! -- local
332  character(len=LINELENGTH) :: line
333  character(len=LINELENGTH) :: keyword
334  integer(I4B) :: ii, n, m
335  integer(I4B) :: ierr
336  logical :: isfound, endOfBlock
337  integer(I4B), parameter :: nname = 2
338  logical, dimension(nname) :: lname
339  character(len=24), dimension(nname) :: aname(nname)
340  ! -- formats
341  character(len=*), parameter :: fmtsymerr = &
342  &"(/,'Error in array: ',(a),/, &
343  &'Array is not symmetric in positions: ',2i9,/, &
344  &'Values in these positions are: ', 2(1pg15.6))"
345  character(len=*), parameter :: fmtihcerr = &
346  &"(/,'IHC must be 0, 1, or 2. Found: ',i0)"
347  ! -- data
348  data aname(1)/' IAC'/
349  data aname(2)/' JA'/
350  !
351  ! -- Allocate and initialize dimensions
352  call this%allocate_scalars(name_model)
353  this%nodes = nodes
354  this%nja = nja
355  this%njas = (this%nja - this%nodes) / 2
356  !
357  ! -- Allocate space for connection arrays
358  call this%allocate_arrays()
359  !
360  ! -- get connectiondata block
361  call this%parser%GetBlock('CONNECTIONDATA', isfound, ierr)
362  lname(:) = .false.
363  if (isfound) then
364  write (iout, '(1x,a)') 'PROCESSING CONNECTIONDATA'
365  do
366  call this%parser%GetNextLine(endofblock)
367  if (endofblock) exit
368  call this%parser%GetStringCaps(keyword)
369  select case (keyword)
370  case ('IAC')
371  call readarray(this%parser%iuactive, this%ia, aname(1), 1, &
372  this%nodes, iout, 0)
373  lname(1) = .true.
374  case ('JA')
375  call readarray(this%parser%iuactive, this%ja, aname(2), 1, &
376  this%nja, iout, 0)
377  lname(2) = .true.
378  case default
379  write (errmsg, '(a,a)') &
380  'Unknown CONNECTIONDATA tag: ', trim(keyword)
381  call store_error(errmsg)
382  call this%parser%StoreErrorUnit()
383  end select
384  end do
385  write (iout, '(1x,a)') 'END PROCESSING CONNECTIONDATA'
386  else
387  call store_error('Required CONNECTIONDATA block not found.')
388  call this%parser%StoreErrorUnit()
389  end if
390  !
391  ! -- verify all items were read
392  do n = 1, nname
393  if (.not. lname(n)) then
394  write (errmsg, '(a,a)') &
395  'Required input was not specified: ', aname(n)
396  call store_error(errmsg)
397  end if
398  end do
399  if (count_errors() > 0) then
400  call this%parser%StoreErrorUnit()
401  end if
402  !
403  ! -- Convert iac to ia
404  do n = 2, this%nodes + 1
405  this%ia(n) = this%ia(n) + this%ia(n - 1)
406  end do
407  do n = this%nodes + 1, 2, -1
408  this%ia(n) = this%ia(n - 1) + 1
409  end do
410  this%ia(1) = 1
411  !
412  ! -- Convert any negative ja numbers to positive
413  do ii = 1, this%nja
414  if (this%ja(ii) < 0) this%ja(ii) = -this%ja(ii)
415  end do
416  !
417  ! -- fill the isym and jas arrays
418  call fillisym(this%nodes, this%nja, this%ia, this%ja, this%isym)
419  call filljas(this%nodes, this%nja, this%ia, this%ja, this%isym, &
420  this%jas)
421  !
422  ! -- check for symmetry in ja
423  do n = 1, this%nodes
424  do ii = this%ia(n), this%ia(n + 1) - 1
425  m = this%ja(ii)
426  if (n /= this%ja(this%isym(ii))) then
427  write (line, fmtsymerr) aname(2), ii, this%isym(ii)
428  call write_message(line)
429  call this%parser%StoreErrorUnit()
430  end if
431  end do
432  end do
433  !
434  if (count_errors() > 0) then
435  call this%parser%StoreErrorUnit()
436  end if
Here is the call graph for this function:

◆ set_cl1_cl2_from_fleng()

subroutine connectionsmodule::set_cl1_cl2_from_fleng ( class(connectionstype this,
real(dp), dimension(:), intent(in)  fleng 
)

Definition at line 441 of file Connections.f90.

442  ! -- modules
443  use constantsmodule, only: dhalf
444  ! -- dummy
445  class(ConnectionsType) :: this
446  real(DP), dimension(:), intent(in) :: fleng
447  ! -- local
448  integer(I4B) :: n, m, ii
449  !
450  ! -- Fill symmetric arrays cl1 and cl2 from fleng of the node
451  do n = 1, this%nodes
452  do ii = this%ia(n) + 1, this%ia(n + 1) - 1
453  m = this%ja(ii)
454  this%cl1(this%jas(ii)) = fleng(n) * dhalf
455  this%cl2(this%jas(ii)) = fleng(m) * dhalf
456  end do
457  end do

◆ set_mask()

subroutine connectionsmodule::set_mask ( class(connectionstype this,
integer(i4b), intent(in)  ipos,
integer(i4b), intent(in)  maskval 
)

Definition at line 1427 of file Connections.f90.

1428  ! -- modules
1430  ! -- dummy
1431  class(ConnectionsType) :: this
1432  integer(I4B), intent(in) :: ipos
1433  integer(I4B), intent(in) :: maskval
1434  ! -- local
1435  integer(I4B) :: i
1436  !
1437  ! if we still point to this%ja, we first need to allocate space
1438  if (associated(this%mask, this%ja)) then
1439  call mem_allocate(this%mask, this%nja, 'MASK', this%memoryPath)
1440  ! and initialize with unmasked
1441  do i = 1, this%nja
1442  this%mask(i) = 1
1443  end do
1444  end if
1445  !
1446  ! -- set the mask value
1447  this%mask(ipos) = maskval

◆ vertexconnect()

subroutine connectionsmodule::vertexconnect ( integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nrsize,
integer(i4b), intent(in)  maxnnz,
integer(i4b), intent(in)  nlay,
integer(i4b), intent(in)  ncpl,
type(sparsematrix), intent(inout)  sparse,
type(sparsematrix), intent(inout)  vertcellspm,
type(disvgeomtype), intent(inout)  cell1,
type(disvgeomtype), intent(inout)  cell2,
integer(i4b), dimension(:), intent(in)  nodereduced 
)
private

Definition at line 1272 of file Connections.f90.

1274  ! -- modules
1275  use sparsemodule, only: sparsematrix
1276  use disvgeom, only: disvgeomtype
1277  ! -- dummy
1278  integer(I4B), intent(in) :: nodes
1279  integer(I4B), intent(in) :: nrsize
1280  integer(I4B), intent(in) :: maxnnz
1281  integer(I4B), intent(in) :: nlay
1282  integer(I4B), intent(in) :: ncpl
1283  type(SparseMatrix), intent(inout) :: sparse
1284  type(SparseMatrix), intent(inout) :: vertcellspm
1285  integer(I4B), dimension(:), intent(in) :: nodereduced
1286  type(DisvGeomType), intent(inout) :: cell1, cell2
1287  ! -- local
1288  integer(I4B), dimension(:), allocatable :: rowmaxnnz
1289  integer(I4B) :: i, j, k, kk, nr, mr, j1, j2, icol1, icol2, nvert
1290  !
1291  ! -- Allocate and fill the ia and ja arrays
1292  allocate (rowmaxnnz(nodes))
1293  do i = 1, nodes
1294  rowmaxnnz(i) = maxnnz
1295  end do
1296  call sparse%init(nodes, nodes, rowmaxnnz)
1297  deallocate (rowmaxnnz)
1298  do k = 1, nlay
1299  do j = 1, ncpl
1300  !
1301  ! -- Find the reduced node number and then cycle if the
1302  ! node is always inactive
1303  nr = get_node(k, 1, j, nlay, 1, ncpl)
1304  if (nrsize > 0) nr = nodereduced(nr)
1305  if (nr <= 0) cycle
1306  !
1307  ! -- Process diagonal
1308  call sparse%addconnection(nr, nr, 1)
1309  !
1310  ! -- Up direction
1311  if (k > 1) then
1312  do kk = k - 1, 1, -1
1313  mr = get_node(kk, 1, j, nlay, 1, ncpl)
1314  if (nrsize > 0) mr = nodereduced(mr)
1315  if (mr >= 0) exit
1316  end do
1317  if (mr > 0) then
1318  call sparse%addconnection(nr, mr, 1)
1319  end if
1320  end if
1321  !
1322  ! -- Down direction
1323  if (k < nlay) then
1324  do kk = k + 1, nlay
1325  mr = get_node(kk, 1, j, nlay, 1, ncpl)
1326  if (nrsize > 0) mr = nodereduced(mr)
1327  if (mr >= 0) exit
1328  end do
1329  if (mr > 0) then
1330  call sparse%addconnection(nr, mr, 1)
1331  end if
1332  end if
1333  end do
1334  end do
1335  !
1336  ! -- Go through each vertex and connect up all the cells that use
1337  ! this vertex in their definition and share an edge.
1338  nvert = vertcellspm%nrow
1339  do i = 1, nvert
1340  do icol1 = 1, vertcellspm%row(i)%nnz
1341  j1 = vertcellspm%row(i)%icolarray(icol1)
1342  do k = 1, nlay
1343  nr = get_node(k, 1, j1, nlay, 1, ncpl)
1344  if (nrsize > 0) nr = nodereduced(nr)
1345  if (nr <= 0) cycle
1346  call cell1%set_nodered(nr)
1347  do icol2 = 1, vertcellspm%row(i)%nnz
1348  j2 = vertcellspm%row(i)%icolarray(icol2)
1349  if (j1 == j2) cycle
1350  mr = get_node(k, 1, j2, nlay, 1, ncpl)
1351  if (nrsize > 0) mr = nodereduced(mr)
1352  if (mr <= 0) cycle
1353  call cell2%set_nodered(mr)
1354  if (cell1%shares_edge(cell2)) then
1355  call sparse%addconnection(nr, mr, 1)
1356  end if
1357  end do
1358  end do
1359  end do
1360  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vertexconnectl()

subroutine connectionsmodule::vertexconnectl ( integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nrsize,
integer(i4b), intent(in)  maxnnz,
integer(i4b), intent(in)  nodeuser,
type(sparsematrix), intent(inout)  sparse,
integer(i4b), dimension(:), intent(in)  iavertcells,
integer(i4b), dimension(:), intent(in)  javertcells,
integer(i4b), dimension(:), intent(in)  nodereduced 
)
Parameters
[in]nodesnumber of active nodes
[in]nrsizeif > 0 then reduced grid
[in]maxnnzmax number of non zeros
[in]nodeusernumber of user nodes
[in,out]sparsesparse matrix object
[in]nodereducedmap from user to reduced node
[in]iavertcellscsr ia index array for vertices
[in]javertcellscsr list of cells that use each vertex

Definition at line 1366 of file Connections.f90.

1369  ! modules
1370  use sparsemodule, only: sparsematrix
1371  use geomutilmodule, only: get_node
1372  ! dummy
1373  integer(I4B), intent(in) :: nodes !< number of active nodes
1374  integer(I4B), intent(in) :: nrsize !< if > 0 then reduced grid
1375  integer(I4B), intent(in) :: maxnnz !< max number of non zeros
1376  integer(I4B), intent(in) :: nodeuser !< number of user nodes
1377  type(SparseMatrix), intent(inout) :: sparse !< sparse matrix object
1378  integer(I4B), dimension(:), intent(in) :: nodereduced !< map from user to reduced node
1379  integer(I4B), dimension(:), intent(in) :: iavertcells !< csr ia index array for vertices
1380  integer(I4B), dimension(:), intent(in) :: javertcells !< csr list of cells that use each vertex
1381  ! local
1382  integer(I4B), dimension(:), allocatable :: rowmaxnnz
1383  integer(I4B) :: i, k, nr, mr, nvert
1384  integer(I4B) :: con
1385 
1386  ! Setup a sparse object
1387  allocate (rowmaxnnz(nodes))
1388  do i = 1, nodes
1389  rowmaxnnz(i) = maxnnz
1390  end do
1391  call sparse%init(nodes, nodes, rowmaxnnz)
1392  deallocate (rowmaxnnz)
1393 
1394  ! Fill the sparse diagonal
1395  do nr = 1, nodeuser
1396  mr = nr
1397  if (nrsize > 0) mr = nodereduced(mr)
1398  if (mr <= 0) cycle
1399  call sparse%addconnection(mr, mr, 1)
1400  end do
1401 
1402  ! Go through each vertex and connect up all the cells that use
1403  ! this vertex in their definition.
1404  nvert = size(iavertcells) - 1
1405  do i = 1, nvert
1406  ! loop through cells that share the vertex
1407  do k = iavertcells(i), iavertcells(i + 1) - 2
1408  ! loop again through connected cells that share vertex
1409  do con = k + 1, iavertcells(i + 1) - 1
1410  nr = javertcells(k)
1411  if (nrsize > 0) nr = nodereduced(nr)
1412  if (nr <= 0) cycle
1413  mr = javertcells(con)
1414  if (nrsize > 0) mr = nodereduced(mr)
1415  if (mr <= 0) cycle
1416  call sparse%addconnection(nr, mr, 1)
1417  call sparse%addconnection(mr, nr, 1)
1418  end do
1419  end do
1420  end do
1421 
Here is the call graph for this function:
Here is the caller graph for this function: