14 integer(I4B) :: nodeusr
15 integer(I4B) :: nodered
22 real(dp),
pointer,
dimension(:) :: top_grid => null()
23 real(dp),
pointer,
dimension(:) :: bot_grid => null()
24 integer(I4B),
pointer,
dimension(:) :: iavert => null()
25 integer(I4B),
pointer,
dimension(:) :: javert => null()
26 real(dp),
pointer,
dimension(:, :) :: vertex_grid => null()
27 real(dp),
pointer,
dimension(:, :) :: cellxy_grid => null()
28 integer(I4B),
pointer,
dimension(:, :) :: nodereduced => null()
29 integer(I4B),
pointer,
dimension(:) :: nodeuser => null()
50 subroutine init(this, nlay, ncpl, nodes, top_grid, bot_grid, iavert, &
51 javert, vertex_grid, cellxy_grid, nodereduced, nodeuser)
54 integer(I4B),
intent(in) :: nlay
55 integer(I4B),
intent(in) :: ncpl
56 integer(I4B),
intent(in) :: nodes
57 real(DP),
dimension(nodes),
target :: top_grid
58 real(DP),
dimension(nodes),
target :: bot_grid
59 integer(I4B),
dimension(:),
target :: iavert
60 integer(I4B),
dimension(:),
target :: javert
61 real(DP),
dimension(:, :),
target :: vertex_grid
62 real(DP),
dimension(:, :),
target :: cellxy_grid
63 integer(I4B),
dimension(ncpl, nlay),
target :: nodereduced
64 integer(I4B),
dimension(nodes),
target :: nodeuser
66 integer(I4B) :: nodesuser
71 this%top_grid => top_grid
72 this%bot_grid => bot_grid
75 this%vertex_grid => vertex_grid
76 this%cellxy_grid => cellxy_grid
77 this%nodereduced => nodereduced
78 this%nodeuser => nodeuser
79 nodesuser = ncpl * nlay
81 if (nodes < nodesuser)
then
84 this%reduced = .false.
93 integer(I4B),
intent(in) :: k
94 integer(I4B),
intent(in) :: j
98 this%nodeusr =
get_node(k, 1, j, this%nlay, 1, this%ncpl)
99 if (this%reduced)
then
100 this%nodered = this%nodereduced(k, j)
102 this%nodered = this%nodeusr
104 call this%cell_setup()
115 integer(I4B),
intent(in) :: nodered
117 this%nodered = nodered
119 if (this%reduced)
then
120 this%nodeusr = this%nodeuser(nodered)
122 this%nodeusr = nodered
125 call get_jk(this%nodeusr, this%ncpl, this%nlay, this%j, this%k)
126 call this%cell_setup()
138 this%top = this%top_grid(this%nodered)
139 this%bot = this%bot_grid(this%nodered)
142 subroutine cprops(this, cell2, hwva, cl1, cl2, ax, ihc)
148 real(DP),
intent(out) :: hwva
149 real(DP),
intent(out) :: cl1
150 real(DP),
intent(out) :: cl2
151 real(DP),
intent(out) :: ax
152 integer(I4B),
intent(out) :: ihc
154 integer(I4B) :: ivert1, ivert2
155 integer(I4B) :: istart1, istart2, istop1, istop2
156 real(DP) :: x0, y0, x1, y1, x2, y2
158 if (this%j == cell2%j)
then
162 hwva = this%get_area()
163 cl1 =
dhalf * (this%top - this%bot)
164 cl2 =
dhalf * (cell2%top - cell2%bot)
170 istart1 = this%iavert(this%j)
171 istop1 = this%iavert(this%j + 1) - 1
172 istart2 = cell2%iavert(cell2%j)
173 istop2 = this%iavert(cell2%j + 1) - 1
175 this%javert(istart2:istop2), &
177 if (ivert1 == 0 .or. ivert2 == 0)
then
184 x1 = this%vertex_grid(1, ivert1)
185 y1 = this%vertex_grid(2, ivert1)
186 x2 = this%vertex_grid(1, ivert2)
187 y2 = this%vertex_grid(2, ivert2)
191 x0 = this%cellxy_grid(1, this%j)
192 y0 = this%cellxy_grid(2, this%j)
196 x0 = this%cellxy_grid(1, cell2%j)
197 y0 = this%cellxy_grid(2, cell2%j)
201 x1 = this%vertex_grid(1, ivert1)
202 y1 = this%vertex_grid(2, ivert1)
203 x2 = this%vertex_grid(1, ivert2)
204 y2 = this%vertex_grid(2, ivert2)
205 ax =
anglex(x1, y1, x2, y2)
221 real(DP),
intent(out) :: xcomp
222 real(DP),
intent(out) :: ycomp
224 integer(I4B) :: ivert1, ivert2
225 integer(I4B) :: istart1, istart2, istop1, istop2
226 real(DP) :: x1, y1, x2, y2
228 istart1 = this%iavert(this%j)
229 istop1 = this%iavert(this%j + 1) - 1
230 istart2 = cell2%iavert(cell2%j)
231 istop2 = this%iavert(cell2%j + 1) - 1
233 this%javert(istart2:istop2), &
235 x1 = this%vertex_grid(1, ivert1)
236 y1 = this%vertex_grid(2, ivert1)
237 x2 = this%vertex_grid(1, ivert2)
238 y2 = this%vertex_grid(2, ivert2)
251 ycomp, zcomp, conlen)
257 logical,
intent(in) :: nozee
258 real(DP),
intent(in) :: satn
259 real(DP),
intent(in) :: satm
260 real(DP),
intent(out) :: xcomp
261 real(DP),
intent(out) :: ycomp
262 real(DP),
intent(out) :: zcomp
263 real(DP),
intent(out) :: conlen
265 real(DP) :: x1, y1, z1, x2, y2, z2
267 x1 = this%cellxy_grid(1, this%j)
268 y1 = this%cellxy_grid(2, this%j)
269 x2 = this%cellxy_grid(1, cell2%j)
270 y2 = this%cellxy_grid(2, cell2%j)
275 z1 = this%bot +
dhalf * satn * (this%top - this%bot)
276 z2 = cell2%bot +
dhalf * satm * (cell2%top - cell2%bot)
295 integer(I4B) :: istart1, istop1, istart2, istop2
296 integer(I4B) :: ivert1, ivert2
298 istart1 = this%iavert(this%j)
299 istop1 = this%iavert(this%j + 1) - 1
300 istart2 = cell2%iavert(cell2%j)
301 istop2 = this%iavert(cell2%j + 1) - 1
303 this%javert(istart2:istop2), &
306 if (ivert1 == 0 .or. ivert2 == 0)
then
322 integer(I4B),
dimension(:) :: ivlist1
323 integer(I4B),
dimension(:) :: ivlist2
324 integer(I4B),
intent(out) :: ivert1
325 integer(I4B),
intent(out) :: ivert2
338 outerloop:
do il1 = 1, nv1 - 1
340 if (ivlist1(il1) == ivlist2(il2) .and. &
341 ivlist1(il1 + 1) == ivlist2(il2 - 1))
then
343 ivert1 = ivlist1(il1)
344 ivert2 = ivlist1(il1 + 1)
365 integer(I4B) :: ivert
366 integer(I4B) :: nvert
367 integer(I4B) :: icount
375 nvert = this%iavert(this%j + 1) - this%iavert(this%j)
377 iv1 = this%javert(this%iavert(this%j))
378 x1 = this%vertex_grid(1, iv1)
379 y1 = this%vertex_grid(2, iv1)
380 do ivert = this%iavert(this%j), this%iavert(this%j + 1) - 1
381 x = this%vertex_grid(1, this%javert(ivert))
382 if (icount < nvert)
then
383 y = this%vertex_grid(2, this%javert(ivert + 1))
385 y = this%vertex_grid(2, this%javert(this%iavert(this%j)))
387 area = area + (x - x1) * (y - y1)
392 do ivert = this%iavert(this%j), this%iavert(this%j + 1) - 1
393 y = this%vertex_grid(2, this%javert(ivert))
394 if (icount < nvert)
then
395 x = this%vertex_grid(1, this%javert(ivert + 1))
397 x = this%vertex_grid(1, this%javert(this%iavert(this%j)))
399 area = area - (x - x1) * (y - y1)
403 area = abs(area) *
dhalf
415 function anglex(x1, y1, x2, y2)
result(ax)
419 real(dp),
intent(in) :: x1
420 real(dp),
intent(in) :: x2
421 real(dp),
intent(in) :: y1
422 real(dp),
intent(in) :: y2
442 real(dp),
intent(in) :: x1
443 real(dp),
intent(in) :: x2
444 real(dp),
intent(in) :: y1
445 real(dp),
intent(in) :: y2
449 d = (x1 - x2)**2 + (y1 - y2)**2
461 real(dp),
intent(in) :: x0
462 real(dp),
intent(in) :: y0
463 real(dp),
intent(in) :: x1
464 real(dp),
intent(in) :: y1
465 real(dp),
intent(in) :: x2
466 real(dp),
intent(in) :: y2
470 d = abs((x2 - x1) * (y1 - y0) - (x1 - x0) * (y2 - y1))
482 real(DP),
intent(in) :: x0
483 real(DP),
intent(in) :: y0
484 real(DP),
intent(in) :: x1
485 real(DP),
intent(in) :: y1
486 real(DP),
intent(out) :: xcomp
487 real(DP),
intent(out) :: ycomp
489 real(DP) :: dx, dy, vmag
493 vmag = sqrt(dx**2 + dy**2)
507 xcomp, ycomp, zcomp, vmag)
509 real(dp),
intent(in) :: x0
510 real(dp),
intent(in) :: y0
511 real(dp),
intent(in) :: z0
512 real(dp),
intent(in) :: x1
513 real(dp),
intent(in) :: y1
514 real(dp),
intent(in) :: z1
515 real(dp),
intent(out) :: xcomp
516 real(dp),
intent(out) :: ycomp
517 real(dp),
intent(out) :: zcomp
520 real(dp) :: dx, dy, dz
525 vmag = sqrt(dx**2 + dy**2 + dz**2)
This module contains simulation constants.
real(dp), parameter dhalf
real constant 1/2
real(dp), parameter dpi
real constant
real(dp), parameter dzero
real constant zero
real(dp), parameter dtwo
real constant 2
real(dp), parameter done
real constant 1
subroutine cprops(this, cell2, hwva, cl1, cl2, ax, ihc)
real(dp) function anglex(x1, y1, x2, y2)
Calculate the angle that the x-axis makes with a line that is normal to the two points.
subroutine line_unit_normal(x0, y0, x1, y1, xcomp, ycomp)
Calculate the normal vector components (xcomp and ycomp) for a line defined by two points,...
subroutine cell_setup(this)
Set top and bottom elevations of grid cell.
subroutine shared_edge(ivlist1, ivlist2, ivert1, ivert2)
Find two common vertices shared by cell1 and cell2.
real(dp) function distance_normal(x0, y0, x1, y1, x2, y2)
Calculate normal distance from point (x0, y0) to line defined by two points, (x1, y1),...
subroutine set_nodered(this, nodered)
Set reduced node number.
logical function shares_edge(this, cell2)
Return true if this shares a horizontal edge with cell2.
subroutine edge_normal(this, cell2, xcomp, ycomp)
Return the x and y components of an outward normal facing vector.
real(dp) function get_area(this)
Calculate and return the area of the cell.
subroutine connection_vector(this, cell2, nozee, satn, satm, xcomp, ycomp, zcomp, conlen)
Return the x y and z components of a unit vector that points from from the center of this to the cent...
subroutine, public line_unit_vector(x0, y0, z0, x1, y1, z1, xcomp, ycomp, zcomp, vmag)
Calculate the vector components (xcomp, ycomp, and zcomp) for a line defined by two points,...
real(dp) function distance(x1, y1, x2, y2)
Calculate distance between two points.
subroutine set_kj(this, k, j)
Set node IDs.
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...
subroutine, public get_jk(nodenumber, ncpl, nlay, icpl, ilay)
Get layer index and within-layer node index from node number and grid dimensions. If nodenumber is in...
This module defines variable data types.