MODFLOW 6  version 6.5.0.dev2
MODFLOW 6 Code Documentation
Data Types | Functions/Subroutines
uzfcellgroupmodule Module Reference

Data Types

type  uzfcellgrouptype
 

Functions/Subroutines

subroutine init (this, ncells, nwav, memory_path)
 Allocate and set uzf object variables. More...
 
subroutine dealloc (this)
 Deallocate uzf object variables. More...
 
subroutine setdata (this, icell, area, top, bot, surfdep, vks, thtr, thts, thti, eps, ntrail, landflag, ivertcon)
 Set uzf object material properties. More...
 
subroutine sethead (this, icell, hgwf)
 Set initial head for uzf object. More...
 
subroutine setdatafinf (this, icell, finf)
 Set infiltration. More...
 
subroutine setdatauzfarea (this, icell, areamult)
 Set uzfarea using cellarea and areamult. More...
 
subroutine setdataet (this, icell, jbelow, pet, extdp)
 Set unsaturated ET-related variables. More...
 
subroutine setgwpet (this, icell)
 Subtract aet from pet to calculate residual et for gw. More...
 
subroutine setbelowpet (this, icell, jbelow)
 Subtract aet from pet to calculate residual et for deeper cells. More...
 
subroutine setdataetwc (this, icell, jbelow, extwc)
 Set extinction water content. More...
 
subroutine setdataetha (this, icell, jbelow, ha, hroot, rootact)
 Set variables for head-based unsaturated flow. More...
 
subroutine advance (this, icell)
 Set variables to advance to new time step. nothing yet. More...
 
subroutine solve (this, thiswork, jbelow, icell, totfluxtot, ietflag, issflag, iseepflag, hgwf, qfrommvr, ierr, reset_state, trhs, thcof, deriv, watercontent)
 Formulate the unsaturated flow object, calculate terms for gwf equation. More...
 
subroutine addrech (this, icell, jbelow, hgwf, trhs, thcof, deriv, delt)
 Add recharge or infiltration to cells. More...
 
subroutine rejfinf (this, icell, deriv, hgwf, trhs, thcof, finfact)
 Reject applied infiltration due to low vks. More...
 
subroutine gwseep (this, icell, deriv, scale, hgwf, trhs, thcof, seep)
 Calculate groudwater discharge to land surface. More...
 
subroutine simgwet (this, igwetflag, icell, hgwf, trhs, thcof, det)
 Calculate gwf et using residual uzf pet. More...
 
real(dp) function etfunc_lin (s, x, c, det, trhs, thcof, hgwf, celtop, celbot)
 Calculate gwf et using linear ET function from mf-2005. More...
 
real(dp) function etfunc_nlin (s, x, c, det, trhs, thcof, hgwf)
 Square-wave ET function with smoothing at extinction depth. More...
 
subroutine uz_rise (this, icell, totfluxtot)
 Calculate recharge due to a rise in the gwf head. More...
 
subroutine setwaves (this, icell)
 Reset waves to default values at start of simulation. More...
 
subroutine routewaves (this, totfluxtot, delt, ietflag, icell, ierr)
 Prepare and route waves over time step. More...
 
subroutine wave_shift (this, this2, icell, icell2, shft, strt, stp, cntr)
 Copy waves or shift waves in arrays. More...
 
subroutine uzflow (this, thick, thickold, delt, ietflag, icell, ierr)
 Method of Characteristics solution for kinematic wave equation. More...
 
subroutine factors (feps1, feps2)
 Calculate unit specific tolerances. More...
 
subroutine trailwav (this, icell, ierr)
 Create and set trail waves. More...
 
subroutine leadwav (this, time, itester, itrailflg, thetab, fluxb, ffcheck, feps2, delt, icell)
 Create a lead wave and route over time step. More...
 
real(dp) function leadspeed (theta1, theta2, flux1, flux2, thts, thtr, eps, vks)
 Calculates waves speed from dflux/dtheta. More...
 
real(dp) function unsat_stor (this, icell, d1)
 Sums up mobile water over depth interval. More...
 
subroutine update_wav (this, icell, delt, iss, itest)
 Update to new state of uz at end of time step. More...
 
subroutine uzet (this, icell, delt, ietflag, ierr)
 Remove water from uz due to et. More...
 
real(dp) function caph (this, icell, tho)
 Calculate capillary pressure head from B-C equation. More...
 
real(dp) function rate_et_z (this, icell, factor, fktho, h)
 Calculate capillary pressure-based uz et. More...
 
real(dp) function get_water_content_at_depth (this, icell, depth)
 Determine the water content at a specific depth. More...
 
real(dp) function get_wcnew (this, icell)
 Calculate and return the cell-based water content value. More...
 

Function/Subroutine Documentation

◆ addrech()

subroutine uzfcellgroupmodule::addrech ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  jbelow,
real(dp), intent(in)  hgwf,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(inout)  deriv,
real(dp), intent(in)  delt 
)

Definition at line 735 of file UzfCellGroup.f90.

736  ! -- dummy
737  class(UzfCellGroupType) :: this
738  integer(I4B), intent(in) :: icell
739  integer(I4B), intent(in) :: jbelow
740  real(DP), intent(inout) :: trhs
741  real(DP), intent(inout) :: thcof
742  real(DP), intent(inout) :: deriv
743  real(DP), intent(in) :: delt
744  real(DP), intent(in) :: hgwf
745  ! -- local
746  real(DP) :: fcheck
747  real(DP) :: x, scale, range
748  !
749  ! -- initialize
750  range = dem5
751  deriv = dzero
752  thcof = dzero
753  trhs = this%uzfarea(icell) * this%totflux(icell) / delt
754  if (this%totflux(icell) < dem14) return
755  scale = done
756  !
757  ! -- smoothly reduce flow between cells when head close to cell top
758  x = hgwf - (this%celbot(icell) - range)
759  call sscurve(x, range, deriv, scale)
760  deriv = this%uzfarea(icell) * deriv * this%totflux(icell) / delt
761  this%finf(jbelow) = (done - scale) * this%totflux(icell) / delt
762  fcheck = this%finf(jbelow) - this%vks(jbelow)
763  !
764  ! -- reduce flow between cells when vks is too small
765  if (fcheck < dem14) fcheck = dzero
766  this%finf(jbelow) = this%finf(jbelow) - fcheck
767  this%surfluxbelow(icell) = this%finf(jbelow)
768  this%totflux(icell) = scale * this%totflux(icell) + fcheck * delt
769  trhs = this%uzfarea(icell) * this%totflux(icell) / delt
770  !
771  ! -- Return
772  return
Here is the call graph for this function:

◆ advance()

subroutine uzfcellgroupmodule::advance ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell 
)
private

Definition at line 580 of file UzfCellGroup.f90.

581  ! -- dummy
582  class(UzfCellGroupType) :: this
583  integer(I4B), intent(in) :: icell
584  !
585  ! -- set variables
586  this%surfseep(icell) = dzero
587  !
588  ! -- Return
589  return

◆ caph()

real(dp) function uzfcellgroupmodule::caph ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  tho 
)
private

Definition at line 2042 of file UzfCellGroup.f90.

2043  ! -- dummy
2044  class(UzfCellGroupType) :: this
2045  integer(I4B), intent(in) :: icell
2046  real(DP), intent(in) :: tho
2047  ! -- local
2048  real(DP) :: caph, lambda, star
2049  !
2050  caph = -dem6
2051  star = (tho - this%thtr(icell)) / (this%thts(icell) - this%thtr(icell))
2052  if (star < dem15) star = dem15
2053  lambda = dtwo / (this%eps(icell) - dthree)
2054  if (star > dem15) then
2055  if (tho - this%thts(icell) < dem15) then
2056  caph = this%ha(icell) * star**(-done / lambda)
2057  else
2058  caph = dzero
2059  end if
2060  end if
2061  !
2062  ! -- Return
2063  return

◆ dealloc()

subroutine uzfcellgroupmodule::dealloc ( class(uzfcellgrouptype)  this)

Definition at line 243 of file UzfCellGroup.f90.

244  ! -- modules
245  use memorymanagermodule, only: mem_deallocate
246  ! -- dummy
247  class(UzfCellGroupType) :: this
248  !
249  ! -- deallocate based on whether or not memory manager was used
250  if (this%imem_manager == 0) then
251  deallocate (this%uzdpst)
252  deallocate (this%uzthst)
253  deallocate (this%uzflst)
254  deallocate (this%uzspst)
255  deallocate (this%nwavst)
256  deallocate (this%thtr)
257  deallocate (this%thts)
258  deallocate (this%thti)
259  deallocate (this%eps)
260  deallocate (this%ha)
261  deallocate (this%hroot)
262  deallocate (this%rootact)
263  deallocate (this%extwc)
264  deallocate (this%etact)
265  deallocate (this%nwav)
266  deallocate (this%ntrail)
267  deallocate (this%totflux)
268  deallocate (this%sinf)
269  deallocate (this%finf)
270  deallocate (this%finf_rej)
271  deallocate (this%gwet)
272  deallocate (this%uzfarea)
273  deallocate (this%cellarea)
274  deallocate (this%celtop)
275  deallocate (this%celbot)
276  deallocate (this%landtop)
277  deallocate (this%watab)
278  deallocate (this%watabold)
279  deallocate (this%surfdep)
280  deallocate (this%vks)
281  deallocate (this%surflux)
282  deallocate (this%surfluxbelow)
283  deallocate (this%surfseep)
284  deallocate (this%gwpet)
285  deallocate (this%pet)
286  deallocate (this%petmax)
287  deallocate (this%extdp)
288  deallocate (this%extdpuz)
289  deallocate (this%landflag)
290  deallocate (this%ivertcon)
291  else
292  call mem_deallocate(this%uzdpst)
293  call mem_deallocate(this%uzthst)
294  call mem_deallocate(this%uzflst)
295  call mem_deallocate(this%uzspst)
296  call mem_deallocate(this%nwavst)
297  call mem_deallocate(this%thtr)
298  call mem_deallocate(this%thts)
299  call mem_deallocate(this%thti)
300  call mem_deallocate(this%eps)
301  call mem_deallocate(this%ha)
302  call mem_deallocate(this%hroot)
303  call mem_deallocate(this%rootact)
304  call mem_deallocate(this%extwc)
305  call mem_deallocate(this%etact)
306  call mem_deallocate(this%nwav)
307  call mem_deallocate(this%ntrail)
308  call mem_deallocate(this%totflux)
309  call mem_deallocate(this%sinf)
310  call mem_deallocate(this%finf)
311  call mem_deallocate(this%finf_rej)
312  call mem_deallocate(this%gwet)
313  call mem_deallocate(this%uzfarea)
314  call mem_deallocate(this%cellarea)
315  call mem_deallocate(this%celtop)
316  call mem_deallocate(this%celbot)
317  call mem_deallocate(this%landtop)
318  call mem_deallocate(this%watab)
319  call mem_deallocate(this%watabold)
320  call mem_deallocate(this%surfdep)
321  call mem_deallocate(this%vks)
322  call mem_deallocate(this%surflux)
323  call mem_deallocate(this%surfluxbelow)
324  call mem_deallocate(this%surfseep)
325  call mem_deallocate(this%gwpet)
326  call mem_deallocate(this%pet)
327  call mem_deallocate(this%petmax)
328  call mem_deallocate(this%extdp)
329  call mem_deallocate(this%extdpuz)
330  call mem_deallocate(this%landflag)
331  call mem_deallocate(this%ivertcon)
332  end if
333  !
334  ! -- Return
335  return
memorymanagermodule::mem_deallocate
Definition: MemoryManager.f90:121

◆ etfunc_lin()

real(dp) function uzfcellgroupmodule::etfunc_lin ( real(dp), intent(in)  s,
real(dp), intent(in)  x,
real(dp), intent(in)  c,
real(dp), intent(inout)  det,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(in)  hgwf,
real(dp), intent(in)  celtop,
real(dp), intent(in)  celbot 
)
private

Definition at line 891 of file UzfCellGroup.f90.

892  ! -- Return
893  real(DP) :: etfunc_lin
894  ! -- dummy
895  real(DP), intent(in) :: s
896  real(DP), intent(in) :: x
897  real(DP), intent(in) :: c
898  real(DP), intent(inout) :: det
899  real(DP), intent(inout) :: trhs
900  real(DP), intent(inout) :: thcof
901  real(DP), intent(in) :: hgwf
902  real(DP), intent(in) :: celtop
903  real(DP), intent(in) :: celbot
904  ! -- local
905  real(DP) :: etgw
906  real(DP) :: range
907  real(DP) :: depth, scale, thick
908  !
909  ! -- Between ET surface and extinction depth
910  if (hgwf > (s - x) .and. hgwf < s) THEN
911  etgw = (c * (hgwf - (s - x)) / x)
912  if (etgw > c) then
913  etgw = c
914  else
915  trhs = c - c * s / x
916  thcof = -c / x
917  etgw = trhs - (thcof * hgwf)
918  end if
919  !
920  ! -- Above land surface
921  else if (hgwf >= s) then
922  trhs = c
923  etgw = c
924  !
925  ! Below extinction depth
926  else
927  etgw = dzero
928  end if
929  !
930  ! -- calculate rate
931  depth = hgwf - (s - x)
932  thick = celtop - celbot
933  if (depth > thick) depth = thick
934  if (depth < dzero) depth = dzero
935  range = dem4 * x
936  call scubic(depth, range, det, scale)
937  trhs = scale * trhs
938  thcof = scale * thcof
939  etgw = trhs - (thcof * hgwf)
940  det = -det * etgw
941  etfunc_lin = etgw
942  !
943  ! -- Return
944  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ etfunc_nlin()

real(dp) function uzfcellgroupmodule::etfunc_nlin ( real(dp), intent(in)  s,
real(dp), intent(in)  x,
real(dp), intent(in)  c,
real(dp), intent(inout)  det,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(in)  hgwf 
)
private

Definition at line 949 of file UzfCellGroup.f90.

950  ! -- Return
951  real(DP) :: etfunc_nlin
952  ! -- dummy
953  real(DP), intent(in) :: s
954  real(DP), intent(in) :: x
955  real(DP), intent(in) :: c
956  real(DP), intent(inout) :: det
957  real(DP), intent(inout) :: trhs
958  real(DP), intent(inout) :: thcof
959  real(DP), intent(in) :: hgwf
960  ! -- local
961  real(DP) :: etgw
962  real(DP) :: range
963  real(DP) :: depth, scale
964  !
965  depth = hgwf - (s - x)
966  if (depth < dzero) depth = dzero
967  etgw = c
968  range = dem3 * x
969  call scubic(depth, range, det, scale)
970  etgw = etgw * scale
971  trhs = -etgw
972  det = -det * etgw
973  etfunc_nlin = etgw
974  !
975  ! -- Return
976  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ factors()

subroutine uzfcellgroupmodule::factors ( real(dp), intent(out)  feps1,
real(dp), intent(out)  feps2 
)
private

Definition at line 1217 of file UzfCellGroup.f90.

1218  ! -- dummy
1219  real(DP), intent(out) :: feps1
1220  real(DP), intent(out) :: feps2
1221  real(DP) :: factor1
1222  real(DP) :: factor2
1223  !
1224  ! calculate constants for uzflow
1225  factor1 = done
1226  factor2 = done
1227  feps1 = dem9
1228  feps2 = dem9
1229  if (itmuni == 1) then
1230  factor1 = done / 86400.d0
1231  else if (itmuni == 2) then
1232  factor1 = done / 1440.d0
1233  else if (itmuni == 3) then
1234  factor1 = done / 24.0d0
1235  else if (itmuni == 5) then
1236  factor1 = 365.0d0
1237  end if
1238  factor2 = done / 0.3048
1239  feps1 = feps1 * factor1 * factor2
1240  feps2 = feps2 * factor1 * factor2
1241  !
1242  ! -- Return
1243  return
Here is the caller graph for this function:

◆ get_water_content_at_depth()

real(dp) function uzfcellgroupmodule::get_water_content_at_depth ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  depth 
)
private

Because UZF-calculated waves are internal to UZF objects, different water contents exists at different depths.

Parameters
[in]icelluzf cell containing depth
[in]depthdepth within the cell

Definition at line 2087 of file UzfCellGroup.f90.

2088  ! -- dummy
2089  class(UzfCellGroupType) :: this
2090  integer(I4B), intent(in) :: icell !< uzf cell containing depth
2091  real(DP), intent(in) :: depth !< depth within the cell
2092  ! -- return
2093  real(DP) :: theta_at_depth
2094  ! -- local
2095  real(DP) :: d1
2096  real(DP) :: d2
2097  real(DP) :: f1
2098  real(DP) :: f2
2099  !
2100  if (this%watab(icell) < this%celtop(icell)) then
2101  if (this%celtop(icell) - depth > this%watab(icell)) then
2102  d1 = depth - dem3
2103  d2 = depth + dem3
2104  f1 = this%unsat_stor(icell, d1)
2105  f2 = this%unsat_stor(icell, d2)
2106  theta_at_depth = this%thtr(icell) + (f2 - f1) / (d2 - d1)
2107  else
2108  theta_at_depth = this%thts(icell)
2109  end if
2110  else
2111  theta_at_depth = this%thts(icell)
2112  end if
2113  !
2114  ! -- Return
2115  return

◆ get_wcnew()

real(dp) function uzfcellgroupmodule::get_wcnew ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell 
)
private
Parameters
[in]icelluzf cell containing depth

Definition at line 2120 of file UzfCellGroup.f90.

2121  ! -- dummy
2122  class(UzfCellGroupType) :: this
2123  integer(I4B), intent(in) :: icell !< uzf cell containing depth
2124  ! -- return
2125  real(DP) :: watercontent
2126  ! -- local
2127  real(DP) :: top
2128  real(DP) :: bot
2129  real(DP) :: theta_r
2130  real(DP) :: thk
2131  real(DP) :: hgwf
2132  real(DP) :: fm
2133  real(DP) :: d
2134  !
2135  hgwf = this%watab(icell)
2136  top = this%celtop(icell)
2137  bot = this%celbot(icell)
2138  thk = top - max(bot, hgwf)
2139  if (thk > dzero) then
2140  theta_r = this%thtr(icell)
2141  d = thk
2142  fm = this%unsat_stor(icell, d)
2143  watercontent = fm / thk
2144  watercontent = watercontent + theta_r
2145  else
2146  watercontent = dzero
2147  end if
2148  !
2149  ! -- Return
2150  return

◆ gwseep()

subroutine uzfcellgroupmodule::gwseep ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(inout)  deriv,
real(dp), intent(out)  scale,
real(dp), intent(in)  hgwf,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(inout)  seep 
)
private

Definition at line 808 of file UzfCellGroup.f90.

809  ! -- dummy
810  class(UzfCellGroupType) :: this
811  integer(I4B), intent(in) :: icell
812  real(DP), intent(inout) :: deriv
813  real(DP), intent(inout) :: trhs
814  real(DP), intent(inout) :: thcof
815  real(DP), intent(inout) :: seep
816  real(DP), intent(out) :: scale
817  real(DP), intent(in) :: hgwf
818  ! -- local
819  real(DP) :: x, range, y, deriv1, d1, d2, Q
820  !
821  seep = dzero
822  deriv = dzero
823  deriv1 = dzero
824  d1 = dzero
825  d2 = dzero
826  scale = dzero
827  q = this%uzfarea(icell) * this%vks(icell)
828  range = this%surfdep(icell)
829  x = hgwf - this%celtop(icell)
830  call scubiclinear(x, range, deriv1, y)
831  scale = y
832  seep = scale * q * (hgwf - this%celtop(icell)) / range
833  trhs = scale * q * this%celtop(icell) / range
834  thcof = -scale * q / range
835  d1 = -deriv1 * q * x / range
836  d2 = -scale * q / range
837  deriv = d1 + d2
838  if (seep < dzero) then
839  seep = dzero
840  deriv = dzero
841  trhs = dzero
842  thcof = dzero
843  end if
844  !
845  ! -- Return
846  return
Here is the call graph for this function:

◆ init()

subroutine uzfcellgroupmodule::init ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  ncells,
integer(i4b), intent(in)  nwav,
character(len=*), intent(in), optional  memory_path 
)
private

Definition at line 97 of file UzfCellGroup.f90.

98  ! -- modules
99  use memorymanagermodule, only: mem_allocate
100  ! -- dummy
101  class(UzfCellGroupType) :: this
102  integer(I4B), intent(in) :: nwav
103  integer(I4B), intent(in) :: ncells
104  character(len=*), intent(in), optional :: memory_path
105  ! -- local
106  integer(I4B) :: icell
107  !
108  ! -- Use mem_allocate if memory path is passed in, otherwise it's a temp object
109  if (present(memory_path)) then
110  this%imem_manager = 1
111  call mem_allocate(this%uzdpst, nwav, ncells, 'UZDPST', memory_path)
112  call mem_allocate(this%uzthst, nwav, ncells, 'UZTHST', memory_path)
113  call mem_allocate(this%uzflst, nwav, ncells, 'UZFLST', memory_path)
114  call mem_allocate(this%uzspst, nwav, ncells, 'UZSPST', memory_path)
115  call mem_allocate(this%nwavst, ncells, 'NWAVST', memory_path)
116  call mem_allocate(this%thtr, ncells, 'THTR', memory_path)
117  call mem_allocate(this%thts, ncells, 'THTS', memory_path)
118  call mem_allocate(this%thti, ncells, 'THTI', memory_path)
119  call mem_allocate(this%eps, ncells, 'EPS', memory_path)
120  call mem_allocate(this%ha, ncells, 'HA', memory_path)
121  call mem_allocate(this%hroot, ncells, 'HROOT', memory_path)
122  call mem_allocate(this%rootact, ncells, 'ROOTACT', memory_path)
123  call mem_allocate(this%extwc, ncells, 'EXTWC', memory_path)
124  call mem_allocate(this%etact, ncells, 'ETACT', memory_path)
125  call mem_allocate(this%nwav, ncells, 'NWAV', memory_path)
126  call mem_allocate(this%ntrail, ncells, 'NTRAIL', memory_path)
127  call mem_allocate(this%totflux, ncells, 'TOTFLUX', memory_path)
128  call mem_allocate(this%sinf, ncells, 'SINF', memory_path)
129  call mem_allocate(this%finf, ncells, 'FINF', memory_path)
130  call mem_allocate(this%finf_rej, ncells, 'FINF_REJ', memory_path)
131  call mem_allocate(this%gwet, ncells, 'GWET', memory_path)
132  call mem_allocate(this%uzfarea, ncells, 'UZFAREA', memory_path)
133  call mem_allocate(this%cellarea, ncells, 'CELLAREA', memory_path)
134  call mem_allocate(this%celtop, ncells, 'CELTOP', memory_path)
135  call mem_allocate(this%celbot, ncells, 'CELBOT', memory_path)
136  call mem_allocate(this%landtop, ncells, 'LANDTOP', memory_path)
137  call mem_allocate(this%watab, ncells, 'WATAB', memory_path)
138  call mem_allocate(this%watabold, ncells, 'WATABOLD', memory_path)
139  call mem_allocate(this%surfdep, ncells, 'SURFDEP', memory_path)
140  call mem_allocate(this%vks, ncells, 'VKS', memory_path)
141  call mem_allocate(this%surflux, ncells, 'SURFLUX', memory_path)
142  call mem_allocate(this%surfluxbelow, ncells, 'SURFLUXBELOW', memory_path)
143  call mem_allocate(this%surfseep, ncells, 'SURFSEEP', memory_path)
144  call mem_allocate(this%gwpet, ncells, 'GWPET', memory_path)
145  call mem_allocate(this%pet, ncells, 'PET', memory_path)
146  call mem_allocate(this%petmax, ncells, 'PETMAX', memory_path)
147  call mem_allocate(this%extdp, ncells, 'EXTDP', memory_path)
148  call mem_allocate(this%extdpuz, ncells, 'EXTDPUZ', memory_path)
149  call mem_allocate(this%landflag, ncells, 'LANDFLAG', memory_path)
150  call mem_allocate(this%ivertcon, ncells, 'IVERTCON', memory_path)
151  else
152  this%imem_manager = 0
153  allocate (this%uzdpst(nwav, ncells))
154  allocate (this%uzthst(nwav, ncells))
155  allocate (this%uzflst(nwav, ncells))
156  allocate (this%uzspst(nwav, ncells))
157  allocate (this%nwavst(ncells))
158  allocate (this%thtr(ncells))
159  allocate (this%thts(ncells))
160  allocate (this%thti(ncells))
161  allocate (this%eps(ncells))
162  allocate (this%ha(ncells))
163  allocate (this%hroot(ncells))
164  allocate (this%rootact(ncells))
165  allocate (this%extwc(ncells))
166  allocate (this%etact(ncells))
167  allocate (this%nwav(ncells))
168  allocate (this%ntrail(ncells))
169  allocate (this%totflux(ncells))
170  allocate (this%sinf(ncells))
171  allocate (this%finf(ncells))
172  allocate (this%finf_rej(ncells))
173  allocate (this%gwet(ncells))
174  allocate (this%uzfarea(ncells))
175  allocate (this%cellarea(ncells))
176  allocate (this%celtop(ncells))
177  allocate (this%celbot(ncells))
178  allocate (this%landtop(ncells))
179  allocate (this%watab(ncells))
180  allocate (this%watabold(ncells))
181  allocate (this%surfdep(ncells))
182  allocate (this%vks(ncells))
183  allocate (this%surflux(ncells))
184  allocate (this%surfluxbelow(ncells))
185  allocate (this%surfseep(ncells))
186  allocate (this%gwpet(ncells))
187  allocate (this%pet(ncells))
188  allocate (this%petmax(ncells))
189  allocate (this%extdp(ncells))
190  allocate (this%extdpuz(ncells))
191  allocate (this%landflag(ncells))
192  allocate (this%ivertcon(ncells))
193  end if
194  do icell = 1, ncells
195  this%uzdpst(:, icell) = dzero
196  this%uzthst(:, icell) = dzero
197  this%uzflst(:, icell) = dzero
198  this%uzspst(:, icell) = dzero
199  this%nwavst(icell) = 1
200  this%thtr(icell) = dzero
201  this%thts(icell) = dzero
202  this%thti(icell) = dzero
203  this%eps(icell) = dzero
204  this%ha(icell) = dzero
205  this%hroot(icell) = dzero
206  this%rootact(icell) = dzero
207  this%extwc(icell) = dzero
208  this%etact(icell) = dzero
209  this%nwav(icell) = nwav
210  this%ntrail(icell) = 0
211  this%totflux(icell) = dzero
212  this%sinf(icell) = dzero
213  this%finf(icell) = dzero
214  this%finf_rej(icell) = dzero
215  this%gwet(icell) = dzero
216  this%uzfarea(icell) = dzero
217  this%cellarea(icell) = dzero
218  this%celtop(icell) = dzero
219  this%celbot(icell) = dzero
220  this%landtop(icell) = dzero
221  this%watab(icell) = dzero
222  this%watabold(icell) = dzero
223  this%surfdep(icell) = dzero
224  this%vks(icell) = dzero
225  this%surflux(icell) = dzero
226  this%surfluxbelow(icell) = dzero
227  this%surfseep(icell) = dzero
228  this%gwpet(icell) = dzero
229  this%pet(icell) = dzero
230  this%petmax(icell) = dzero
231  this%extdp(icell) = dzero
232  this%extdpuz(icell) = dzero
233  this%landflag(icell) = 0
234  this%ivertcon(icell) = 0
235  end do
236  !
237  ! -- Return
238  return
memorymanagermodule::mem_allocate
Definition: MemoryManager.f90:53

◆ leadspeed()

real(dp) function uzfcellgroupmodule::leadspeed ( real(dp), intent(in)  theta1,
real(dp), intent(in)  theta2,
real(dp), intent(in)  flux1,
real(dp), intent(inout)  flux2,
real(dp), intent(in)  thts,
real(dp), intent(in)  thtr,
real(dp), intent(in)  eps,
real(dp), intent(in)  vks 
)
private

Definition at line 1541 of file UzfCellGroup.f90.

1542  ! -- Return
1543  real(DP) :: leadspeed
1544  ! -- dummy
1545  real(DP), intent(in) :: theta1
1546  real(DP), intent(in) :: theta2
1547  real(DP), intent(in) :: flux1
1548  real(DP), intent(inout) :: flux2
1549  real(DP), intent(in) :: thts
1550  real(DP), intent(in) :: thtr
1551  real(DP), intent(in) :: eps
1552  real(DP), intent(in) :: vks
1553  ! -- local
1554  real(DP) :: comp1, comp2, thsrinv, epsfksths
1555  real(DP) :: eps_m1, fhold, comp3
1556  !
1557  eps_m1 = eps - done
1558  thsrinv = done / (thts - thtr)
1559  epsfksths = eps * vks * thsrinv
1560  comp1 = theta2 - theta1
1561  comp2 = abs(flux2 - flux1)
1562  comp3 = theta1 - thtr
1563  if (comp2 < dem15) flux2 = flux1 + dem15
1564  if (abs(comp1) < dem30) then
1565  if (comp3 > dem30) fhold = (comp3 * thsrinv)**eps
1566  if (fhold < dem30) fhold = dem30
1567  leadspeed = epsfksths * (fhold**eps_m1)
1568  else
1569  leadspeed = (flux2 - flux1) / (theta2 - theta1)
1570  end if
1571  if (leadspeed < dem30) leadspeed = dem30
1572  !
1573  ! -- Return
1574  return
Here is the caller graph for this function:

◆ leadwav()

subroutine uzfcellgroupmodule::leadwav ( class(uzfcellgrouptype)  this,
real(dp), intent(inout)  time,
integer(i4b), intent(inout)  itester,
integer(i4b), intent(inout)  itrailflg,
real(dp), intent(inout)  thetab,
real(dp), intent(inout)  fluxb,
real(dp), intent(inout)  ffcheck,
real(dp), intent(in)  feps2,
real(dp), intent(in)  delt,
integer(i4b), intent(in)  icell 
)
private

Definition at line 1337 of file UzfCellGroup.f90.

1339  ! -- dummy
1340  class(UzfCellGroupType) :: this
1341  real(DP), intent(inout) :: thetab
1342  real(DP), intent(inout) :: fluxb
1343  real(DP), intent(in) :: feps2
1344  real(DP), intent(inout) :: time
1345  integer(I4B), intent(inout) :: itester
1346  integer(I4B), intent(inout) :: itrailflg
1347  real(DP), intent(inout) :: ffcheck
1348  real(DP), intent(in) :: delt
1349  integer(I4B), intent(in) :: icell
1350  ! -- local
1351  real(DP) :: bottomtime, shortest, fcheck
1352  real(DP) :: eps_m1, timenew, bottom, timedt
1353  real(DP) :: thtsrinv, diff, fluxhld2
1354  real(DP) :: flux1, flux2, theta1, theta2, ftest
1355  real(DP), allocatable, dimension(:) :: checktime
1356  integer(I4B) :: iflx, iremove, j, l
1357  integer(I4B) :: nwavp1, jshort
1358  integer(I4B), allocatable, dimension(:) :: more
1359  !
1360  allocate (checktime(this%nwavst(icell)))
1361  allocate (more(this%nwavst(icell)))
1362  ftest = dzero
1363  eps_m1 = dble(this%eps(icell)) - done
1364  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1365  !
1366  ! -- initialize new wave
1367  if (itrailflg == 0) then
1368  if (ffcheck > feps2) then
1369  this%uzflst(this%nwavst(icell), icell) = this%surflux(icell)
1370  if (this%uzflst(this%nwavst(icell), icell) < dem30) &
1371  this%uzflst(this%nwavst(icell), icell) = dzero
1372  this%uzthst(this%nwavst(icell), icell) = &
1373  (((this%uzflst(this%nwavst(icell), icell) / this%vks(icell))** &
1374  (done / this%eps(icell))) * (this%thts(icell) - this%thtr(icell))) &
1375  + this%thtr(icell)
1376  theta2 = this%uzthst(this%nwavst(icell), icell)
1377  flux2 = this%uzflst(this%nwavst(icell), icell)
1378  flux1 = this%uzflst(this%nwavst(icell) - 1, icell)
1379  theta1 = this%uzthst(this%nwavst(icell) - 1, icell)
1380  this%uzspst(this%nwavst(icell), icell) = &
1381  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1382  this%thtr(icell), this%eps(icell), this%vks(icell))
1383  this%uzdpst(this%nwavst(icell), icell) = dzero
1384  end if
1385  end if
1386  !
1387  ! -- route all waves and interception of waves over times step
1388  diff = done
1389  timedt = dzero
1390  iflx = 0
1391  fluxhld2 = this%uzflst(1, icell)
1392  if (this%nwavst(icell) == 0) itester = 1
1393  if (itester /= 1) then
1394  do while (diff > dem6)
1395  nwavp1 = this%nwavst(icell) + 1
1396  timedt = delt - time
1397  do j = 1, this%nwavst(icell)
1398  checktime(j) = dep20
1399  more(j) = 0
1400  end do
1401  shortest = timedt
1402  if (this%nwavst(icell) > 2) then
1403  j = 2
1404  !
1405  ! -- calculate time until wave overtakes wave ahead
1406  nwavp1 = this%nwavst(icell) + 1
1407  do while (j < nwavp1)
1408  ftest = this%uzspst(j - 1, icell) - this%uzspst(j, icell)
1409  if (abs(ftest) > dem30) then
1410  checktime(j) = (this%uzdpst(j, icell) - &
1411  this%uzdpst(j - 1, icell)) / ftest
1412  if (checktime(j) < dem30) checktime(j) = dep20
1413  end if
1414  j = j + 1
1415  end do
1416  end if
1417  !
1418  ! - calc time until wave reaches bottom of cell
1419  bottomtime = dep20
1420  if (this%nwavst(icell) > 1) then
1421  if (this%uzspst(2, icell) > dzero) then
1422  bottom = this%uzspst(2, icell)
1423  if (bottom < dem15) bottom = dem15
1424  bottomtime = (this%uzdpst(1, icell) - this%uzdpst(2, icell)) / bottom
1425  if (bottomtime < dzero) bottomtime = dem12
1426  end if
1427  end if
1428  !
1429  ! -- calc time for wave interception
1430  jshort = 0
1431  do j = this%nwavst(icell), 3, -1
1432  if (shortest - checktime(j) > -dem9) then
1433  more(j) = 1
1434  jshort = j
1435  shortest = checktime(j)
1436  end if
1437  end do
1438  do j = 3, this%nwavst(icell)
1439  if (shortest - checktime(j) < dem9) then
1440  if (j /= jshort) more(j) = 0
1441  end if
1442  end do
1443  !
1444  ! -- what happens first, waves hits bottom or interception
1445  iremove = 0
1446  timenew = time
1447  fcheck = (time + shortest) - delt
1448  if (shortest < dem7) fcheck = -done
1449  if (bottomtime < shortest .AND. time + bottomtime < delt) then
1450  j = 2
1451  do while (j < nwavp1)
1452  !
1453  ! -- route waves
1454  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1455  this%uzspst(j, icell) * bottomtime
1456  j = j + 1
1457  end do
1458  fluxb = this%uzflst(2, icell)
1459  thetab = this%uzthst(2, icell)
1460  iflx = 1
1461  call this%wave_shift(this, icell, icell, 1, 1, &
1462  this%nwavst(icell) - 1, 1)
1463  iremove = 1
1464  timenew = time + bottomtime
1465  this%uzspst(1, icell) = dzero
1466  !
1467  ! -- do waves intercept before end of time step
1468  else if (fcheck < dzero .AND. this%nwavst(icell) > 2) then
1469  j = 2
1470  do while (j < nwavp1)
1471  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1472  this%uzspst(j, icell) * shortest
1473  j = j + 1
1474  end do
1475  !
1476  ! -- combine waves that intercept, remove a wave
1477  j = 3
1478  l = j
1479  do while (j < this%nwavst(icell) + 1)
1480  if (more(j) == 1) then
1481  l = j
1482  theta2 = this%uzthst(j, icell)
1483  flux2 = this%uzflst(j, icell)
1484  if (j == 3) then
1485  flux1 = fluxb
1486  theta1 = thetab
1487  else
1488  flux1 = this%uzflst(j - 2, icell)
1489  theta1 = this%uzthst(j - 2, icell)
1490  end if
1491  this%uzspst(j, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1492  this%thts(icell), &
1493  this%thtr(icell), &
1494  this%eps(icell), this%vks(icell))
1495  !
1496  ! -- update waves
1497  call this%wave_shift(this, icell, icell, 1, l - 1, &
1498  this%nwavst(icell) - 1, 1)
1499  l = this%nwavst(icell) + 1
1500  iremove = iremove + 1
1501  end if
1502  j = j + 1
1503  end do
1504  timenew = timenew + shortest
1505  !
1506  ! -- calc. total flux to bottom during remaining time in step
1507  else
1508  j = 2
1509  do while (j < nwavp1)
1510  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1511  this%uzspst(j, icell) * timedt
1512  j = j + 1
1513  end do
1514  timenew = delt
1515  end if
1516  this%totflux(icell) = this%totflux(icell) + fluxhld2 * (timenew - time)
1517  if (iflx == 1) then
1518  fluxhld2 = this%uzflst(1, icell)
1519  iflx = 0
1520  end if
1521  !
1522  ! -- remove dead waves
1523  this%nwavst(icell) = this%nwavst(icell) - iremove
1524  time = timenew
1525  diff = delt - time
1526  if (this%nwavst(icell) == 1) then
1527  itester = 1
1528  exit
1529  end if
1530  end do
1531  end if
1532  deallocate (checktime)
1533  deallocate (more)
1534  !
1535  ! -- Return
1536  return
Here is the call graph for this function:

◆ rate_et_z()

real(dp) function uzfcellgroupmodule::rate_et_z ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  factor,
real(dp), intent(in)  fktho,
real(dp), intent(in)  h 
)
private

Definition at line 2067 of file UzfCellGroup.f90.

2068  ! -- Return
2069  real(DP) :: rate_et_z
2070  ! -- dummy
2071  class(UzfCellGroupType) :: this
2072  integer(I4B), intent(in) :: icell
2073  real(DP), intent(in) :: factor, fktho, h
2074  !
2075  rate_et_z = factor * fktho * (h - this%hroot(icell))
2076  if (rate_et_z < dzero) rate_et_z = dzero
2077  !
2078  ! -- Return
2079  return

◆ rejfinf()

subroutine uzfcellgroupmodule::rejfinf ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(inout)  deriv,
real(dp), intent(in)  hgwf,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(inout)  finfact 
)
private

Definition at line 777 of file UzfCellGroup.f90.

778  ! -- dummy
779  class(UzfCellGroupType) :: this
780  integer(I4B), intent(in) :: icell
781  real(DP), intent(inout) :: deriv
782  real(DP), intent(inout) :: finfact
783  real(DP), intent(inout) :: thcof
784  real(DP), intent(inout) :: trhs
785  real(DP), intent(in) :: hgwf
786  ! -- local
787  real(DP) :: x, range, scale, q
788  !
789  range = this%surfdep(icell)
790  q = this%surflux(icell)
791  finfact = q
792  trhs = finfact * this%uzfarea(icell)
793  x = this%celtop(icell) - hgwf
794  call slinear(x, range, deriv, scale)
795  deriv = -q * deriv * this%uzfarea(icell) * scale
796  if (scale < done) then
797  finfact = q * scale
798  trhs = finfact * this%uzfarea(icell) * this%celtop(icell) / range
799  thcof = finfact * this%uzfarea(icell) / range
800  end if
801  !
802  ! -- Return
803  return
Here is the call graph for this function:

◆ routewaves()

subroutine uzfcellgroupmodule::routewaves ( class(uzfcellgrouptype)  this,
real(dp), intent(inout)  totfluxtot,
real(dp), intent(in)  delt,
integer(i4b), intent(in)  ietflag,
integer(i4b), intent(in)  icell,
integer(i4b), intent(inout)  ierr 
)
private

Definition at line 1056 of file UzfCellGroup.f90.

1057  ! -- dummy
1058  class(UzfCellGroupType) :: this
1059  real(DP), intent(inout) :: totfluxtot
1060  real(DP), intent(in) :: delt
1061  integer(I4B), intent(in) :: ietflag
1062  integer(I4B), intent(in) :: icell
1063  integer(I4B), intent(inout) :: ierr
1064  ! -- local
1065  real(DP) :: thick, thickold
1066  integer(I4B) :: idelt, iwav, ik
1067  !
1068  ! -- initialize
1069  this%totflux(icell) = dzero
1070  this%etact(icell) = dzero
1071  thick = this%celtop(icell) - this%watab(icell)
1072  thickold = this%celtop(icell) - this%watabold(icell)
1073  !
1074  ! -- no uz, clear waves
1075  if (thickold < dzero) then
1076  do iwav = 1, 6
1077  this%uzthst(iwav, icell) = this%thtr(icell)
1078  this%uzdpst(iwav, icell) = dzero
1079  this%uzspst(iwav, icell) = dzero
1080  this%uzflst(iwav, icell) = dzero
1081  this%nwavst(icell) = 1
1082  end do
1083  end if
1084  idelt = 1
1085  do ik = 1, idelt
1086  call this%uzflow(thick, thickold, delt, ietflag, icell, ierr)
1087  if (ierr > 0) return
1088  totfluxtot = totfluxtot + this%totflux(icell)
1089  end do
1090  !
1091  ! -- Return
1092  return

◆ setbelowpet()

subroutine uzfcellgroupmodule::setbelowpet ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  jbelow 
)

Definition at line 511 of file UzfCellGroup.f90.

512  ! -- modules
513  use tdismodule, only: delt
514  ! -- dummy
515  class(UzfCellGroupType) :: this
516  integer(I4B), intent(in) :: icell
517  integer(I4B), intent(in) :: jbelow
518  ! -- dummy
519  real(DP) :: pet
520  !
521  pet = dzero
522  !
523  ! -- transfer unmet pet to lower cell
524  !
525  if (this%extdpuz(jbelow) > dem3) then
526  pet = this%pet(icell) - this%etact(icell) / delt - &
527  this%gwet(icell) / this%uzfarea(icell)
528  if (pet < dzero) pet = dzero
529  end if
530  this%pet(jbelow) = pet
531  !
532  ! -- Return
533  return
tdismodule::delt
real(dp), pointer, public delt
length of the current time step
Definition: tdis.f90:29

◆ setdata()

subroutine uzfcellgroupmodule::setdata ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  area,
real(dp), intent(in)  top,
real(dp), intent(in)  bot,
real(dp), intent(in)  surfdep,
real(dp), intent(in)  vks,
real(dp), intent(in)  thtr,
real(dp), intent(in)  thts,
real(dp), intent(in)  thti,
real(dp), intent(in)  eps,
integer(i4b), intent(in)  ntrail,
integer(i4b), intent(in)  landflag,
integer(i4b), intent(in)  ivertcon 
)

Definition at line 340 of file UzfCellGroup.f90.

342  ! -- dummy
343  class(UzfCellGroupType) :: this
344  integer(I4B), intent(in) :: icell
345  real(DP), intent(in) :: area
346  real(DP), intent(in) :: top
347  real(DP), intent(in) :: bot
348  real(DP), intent(in) :: surfdep
349  real(DP), intent(in) :: vks
350  real(DP), intent(in) :: thtr
351  real(DP), intent(in) :: thts
352  real(DP), intent(in) :: thti
353  real(DP), intent(in) :: eps
354  integer(I4B), intent(in) :: ntrail
355  integer(I4B), intent(in) :: landflag
356  integer(I4B), intent(in) :: ivertcon
357  !
358  ! -- set the values for uzf cell icell
359  this%landflag(icell) = landflag
360  this%ivertcon(icell) = ivertcon
361  this%surfdep(icell) = surfdep
362  this%uzfarea(icell) = area
363  this%cellarea(icell) = area
364  if (this%landflag(icell) == 1) then
365  this%celtop(icell) = top - dhalf * this%surfdep(icell)
366  else
367  this%celtop(icell) = top
368  end if
369  this%celbot(icell) = bot
370  this%vks(icell) = vks
371  this%thtr(icell) = thtr
372  this%thts(icell) = thts
373  this%thti(icell) = thti
374  this%eps(icell) = eps
375  this%ntrail(icell) = ntrail
376  this%pet(icell) = dzero
377  this%extdp(icell) = dzero
378  this%extwc(icell) = dzero
379  this%ha(icell) = dzero
380  this%hroot(icell) = dzero

◆ setdataet()

subroutine uzfcellgroupmodule::setdataet ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  jbelow,
real(dp), intent(in)  pet,
real(dp), intent(in)  extdp 
)
private

Definition at line 442 of file UzfCellGroup.f90.

443  ! -- dummy
444  class(UzfCellGroupType) :: this
445  integer(I4B), intent(in) :: icell
446  integer(I4B), intent(in) :: jbelow
447  real(DP), intent(in) :: pet
448  real(DP), intent(in) :: extdp
449  ! -- local
450  real(DP) :: thick
451  !
452  if (this%landflag(icell) == 1) then
453  this%pet(icell) = pet
454  this%gwpet(icell) = pet
455  else
456  this%pet(icell) = dzero
457  this%gwpet(icell) = dzero
458  end if
459  thick = this%celtop(icell) - this%celbot(icell)
460  this%extdp(icell) = extdp
461  if (this%landflag(icell) > 0) then
462  this%landtop(icell) = this%celtop(icell)
463  this%petmax(icell) = this%pet(icell)
464  end if
465  !
466  ! -- set uz extinction depth
467  if (this%landtop(icell) - this%extdp(icell) < this%celbot(icell)) then
468  this%extdpuz(icell) = thick
469  else
470  this%extdpuz(icell) = this%celtop(icell) - &
471  (this%landtop(icell) - this%extdp(icell))
472  end if
473  if (this%extdpuz(icell) < dzero) this%extdpuz(icell) = dzero
474  if (this%extdpuz(icell) > dem7 .and. this%extdp(icell) < dem7) &
475  this%extdp(icell) = this%extdpuz(icell)
476  !
477  ! -- set pet for underlying cell
478  if (jbelow > 0) then
479  this%landtop(jbelow) = this%landtop(icell)
480  this%petmax(jbelow) = this%petmax(icell)
481  end if
482  !
483  ! -- Return
484  return

◆ setdataetha()

subroutine uzfcellgroupmodule::setdataetha ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  jbelow,
real(dp), intent(in)  ha,
real(dp), intent(in)  hroot,
real(dp), intent(in)  rootact 
)
private

Definition at line 555 of file UzfCellGroup.f90.

556  ! -- dummy
557  class(UzfCellGroupType) :: this
558  integer(I4B), intent(in) :: icell
559  integer(I4B), intent(in) :: jbelow
560  real(DP), intent(in) :: ha
561  real(DP), intent(in) :: hroot
562  real(DP), intent(in) :: rootact
563  !
564  ! -- set variables
565  this%ha(icell) = ha
566  this%hroot(icell) = hroot
567  this%rootact(icell) = rootact
568  if (jbelow > 0) then
569  this%ha(jbelow) = ha
570  this%hroot(jbelow) = hroot
571  this%rootact(jbelow) = rootact
572  end if
573  !
574  ! -- Return
575  return

◆ setdataetwc()

subroutine uzfcellgroupmodule::setdataetwc ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  jbelow,
real(dp), intent(in)  extwc 
)

Definition at line 538 of file UzfCellGroup.f90.

539  ! -- dummy
540  class(UzfCellGroupType) :: this
541  integer(I4B), intent(in) :: icell
542  integer(I4B), intent(in) :: jbelow
543  real(DP), intent(in) :: extwc
544  !
545  ! -- set extinction water content
546  this%extwc(icell) = extwc
547  if (jbelow > 0) this%extwc(jbelow) = extwc
548  !
549  ! -- Return
550  return

◆ setdatafinf()

subroutine uzfcellgroupmodule::setdatafinf ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  finf 
)
private

Definition at line 404 of file UzfCellGroup.f90.

405  ! -- dummy
406  class(UzfCellGroupType) :: this
407  integer(I4B), intent(in) :: icell
408  real(DP), intent(in) :: finf
409  !
410  if (this%landflag(icell) == 1) then
411  this%sinf(icell) = finf
412  this%finf(icell) = finf
413  else
414  this%sinf(icell) = dzero
415  this%finf(icell) = dzero
416  end if
417  this%finf_rej(icell) = dzero
418  this%surflux(icell) = dzero
419  this%surfluxbelow(icell) = dzero
420  !
421  ! -- Return
422  return

◆ setdatauzfarea()

subroutine uzfcellgroupmodule::setdatauzfarea ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  areamult 
)
private

Definition at line 427 of file UzfCellGroup.f90.

428  ! -- dummy
429  class(UzfCellGroupType) :: this
430  integer(I4B), intent(in) :: icell
431  real(DP), intent(in) :: areamult
432  !
433  ! -- set uzf area
434  this%uzfarea(icell) = this%cellarea(icell) * areamult
435  !
436  ! -- Return
437  return

◆ setgwpet()

subroutine uzfcellgroupmodule::setgwpet ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell 
)
private

Definition at line 489 of file UzfCellGroup.f90.

490  ! -- modules
491  use tdismodule, only: delt
492  ! -- dummy
493  class(UzfCellGroupType) :: this
494  integer(I4B), intent(in) :: icell
495  ! -- dummy
496  real(DP) :: pet
497  !
498  pet = dzero
499  !
500  ! -- reduce pet for gw by uzet
501  pet = this%pet(icell) - this%etact(icell) / delt
502  if (pet < dzero) pet = dzero
503  this%gwpet(icell) = pet
504  !
505  ! -- Return
506  return

◆ sethead()

subroutine uzfcellgroupmodule::sethead ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  hgwf 
)
private

Definition at line 385 of file UzfCellGroup.f90.

386  ! -- dummy
387  class(UzfCellGroupType) :: this
388  integer(I4B), intent(in) :: icell
389  real(DP), intent(in) :: hgwf
390  !
391  ! -- set initial head
392  this%watab(icell) = this%celbot(icell)
393  if (hgwf > this%celbot(icell)) this%watab(icell) = hgwf
394  if (this%watab(icell) > this%celtop(icell)) &
395  this%watab(icell) = this%celtop(icell)
396  this%watabold(icell) = this%watab(icell)
397  !
398  ! -- Return
399  return

◆ setwaves()

subroutine uzfcellgroupmodule::setwaves ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell 
)
private

Definition at line 1004 of file UzfCellGroup.f90.

1005  ! -- dummy
1006  class(UzfCellGroupType) :: this
1007  ! -- local
1008  integer(I4B), intent(in) :: icell
1009  real(DP) :: bottom, top
1010  integer(I4B) :: jk
1011  real(DP) :: thick
1012  !
1013  ! -- initialize
1014  this%totflux(icell) = dzero
1015  this%nwavst(icell) = 1
1016  this%uzdpst(:, icell) = dzero
1017  thick = this%celtop(icell) - this%watab(icell)
1018  do jk = 1, this%nwav(icell)
1019  this%uzthst(jk, icell) = this%thtr(icell)
1020  end do
1021  !
1022  ! -- initialize waves for first stress period
1023  if (thick > dzero) then
1024  this%uzdpst(1, icell) = thick
1025  this%uzthst(1, icell) = this%thti(icell)
1026  top = this%uzthst(1, icell) - this%thtr(icell)
1027  if (top < dzero) top = dzero
1028  bottom = this%thts(icell) - this%thtr(icell)
1029  if (bottom < dzero) bottom = dzero
1030  this%uzflst(1, icell) = this%vks(icell) * (top / bottom)**this%eps(icell)
1031  if (this%uzthst(1, icell) < this%thtr(icell)) &
1032  this%uzthst(1, icell) = this%thtr(icell)
1033  !
1034  ! -- calculate water stored in the unsaturated zone
1035  if (top > dzero) then
1036  this%uzspst(1, icell) = dzero
1037  else
1038  this%uzflst(1, icell) = dzero
1039  this%uzspst(1, icell) = dzero
1040  end if
1041  !
1042  ! no unsaturated zone
1043  else
1044  this%uzflst(1, icell) = dzero
1045  this%uzdpst(1, icell) = dzero
1046  this%uzspst(1, icell) = dzero
1047  this%uzthst(1, icell) = this%thtr(icell)
1048  end if
1049  !
1050  ! -- Return
1051  return

◆ simgwet()

subroutine uzfcellgroupmodule::simgwet ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  igwetflag,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  hgwf,
real(dp), intent(inout)  trhs,
real(dp), intent(inout)  thcof,
real(dp), intent(inout)  det 
)
private

Definition at line 851 of file UzfCellGroup.f90.

852  ! -- dummy
853  class(UzfCellGroupType) :: this
854  integer(I4B), intent(in) :: igwetflag
855  integer(I4B), intent(in) :: icell
856  real(DP), intent(in) :: hgwf
857  real(DP), intent(inout) :: trhs
858  real(DP), intent(inout) :: thcof
859  real(DP), intent(inout) :: det
860  ! -- local
861  real(DP) :: s, x, c, b, et
862  !
863  this%gwet(icell) = dzero
864  trhs = dzero
865  thcof = dzero
866  det = dzero
867  s = this%landtop(icell)
868  x = this%extdp(icell)
869  c = this%gwpet(icell)
870  b = this%celbot(icell)
871  if (b > hgwf) return
872  if (x < dem6) return
873  if (igwetflag == 1) then
874  et = etfunc_lin(s, x, c, det, trhs, thcof, hgwf, &
875  this%celtop(icell), this%celbot(icell))
876  else if (igwetflag == 2) then
877  et = etfunc_nlin(s, x, c, det, trhs, thcof, hgwf)
878  end if
879  ! this%gwet(icell) = et * this%uzfarea(icell)
880  trhs = trhs * this%uzfarea(icell)
881  thcof = thcof * this%uzfarea(icell)
882  this%gwet(icell) = trhs - (thcof * hgwf)
883  ! write(99,*)'in group', icell, this%gwet(icell)
884  !
885  ! -- Return
886  return
Here is the call graph for this function:

◆ solve()

subroutine uzfcellgroupmodule::solve ( class(uzfcellgrouptype)  this,
type(uzfcellgrouptype)  thiswork,
integer(i4b), intent(in)  jbelow,
integer(i4b), intent(in)  icell,
real(dp), intent(inout)  totfluxtot,
integer(i4b), intent(in)  ietflag,
integer(i4b), intent(in)  issflag,
integer(i4b), intent(in)  iseepflag,
real(dp), intent(in)  hgwf,
real(dp), intent(in)  qfrommvr,
integer(i4b), intent(inout)  ierr,
logical, intent(in)  reset_state,
real(dp), intent(inout), optional  trhs,
real(dp), intent(inout), optional  thcof,
real(dp), intent(inout), optional  deriv,
real(dp), intent(inout), optional  watercontent 
)
private
Parameters
thisworkwork object for resetting wave state
[in]jbelownumber of underlying uzf object or 0 if none
[in]icellnumber of this uzf object
[in]ietflaget is off (0) or based one water content (1) or pressure (2)
[in]issflagsteady state flag
[in]iseepflagdischarge to land is active (1) or not (0)
[in]hgwfhead for cell icell
[in]qfrommvrwater inflow from mover
[in,out]ierrflag indicating not enough waves
[in]reset_stateflag indicating that waves should be reset after solution
[in,out]trhstotal uzf rhs contribution to GWF model
[in,out]thcoftotal uzf hcof contribution to GWF model
[in,out]derivderivate term for contribution to GWF model
[in,out]watercontentcalculated water content

Definition at line 595 of file UzfCellGroup.f90.

598  ! -- modules
599  use tdismodule, only: delt
600  ! -- dummy
601  class(UzfCellGroupType) :: this
602  type(UzfCellGroupType) :: thiswork !< work object for resetting wave state
603  integer(I4B), intent(in) :: jbelow !< number of underlying uzf object or 0 if none
604  integer(I4B), intent(in) :: icell !< number of this uzf object
605  real(DP), intent(inout) :: totfluxtot !<
606  integer(I4B), intent(in) :: ietflag !< et is off (0) or based one water content (1) or pressure (2)
607  integer(I4B), intent(in) :: issflag !< steady state flag
608  integer(I4B), intent(in) :: iseepflag !< discharge to land is active (1) or not (0)
609  real(DP), intent(in) :: hgwf !< head for cell icell
610  real(DP), intent(in) :: qfrommvr !< water inflow from mover
611  integer(I4B), intent(inout) :: ierr !< flag indicating not enough waves
612  logical, intent(in) :: reset_state !< flag indicating that waves should be reset after solution
613  real(DP), intent(inout), optional :: trhs !< total uzf rhs contribution to GWF model
614  real(DP), intent(inout), optional :: thcof !< total uzf hcof contribution to GWF model
615  real(DP), intent(inout), optional :: deriv !< derivate term for contribution to GWF model
616  real(DP), intent(inout), optional :: watercontent !< calculated water content
617  ! -- local
618  real(DP) :: test
619  real(DP) :: scale
620  real(DP) :: seep
621  real(DP) :: finfact
622  real(DP) :: derivfinf
623  real(DP) :: trhsfinf
624  real(DP) :: thcoffinf
625  real(DP) :: trhsseep
626  real(DP) :: thcofseep
627  real(DP) :: deriv1
628  real(DP) :: deriv2
629  !
630  ! -- initialize variables
631  totfluxtot = dzero
632  trhsfinf = dzero
633  thcoffinf = dzero
634  trhsseep = dzero
635  thcofseep = dzero
636  this%finf_rej(icell) = dzero
637  this%surflux(icell) = this%finf(icell) + qfrommvr / this%uzfarea(icell)
638  this%watab(icell) = hgwf
639  this%surfseep(icell) = dzero
640  seep = dzero
641  finfact = dzero
642  this%etact(icell) = dzero
643  this%surfluxbelow(icell) = dzero
644  if (this%ivertcon(icell) > 0) then
645  this%finf(jbelow) = dzero
646  end if
647  if (this%watab(icell) < this%celbot(icell)) &
648  this%watab(icell) = this%celbot(icell)
649  !
650  ! -- initialize derivative variables
651  deriv1 = dzero
652  deriv2 = dzero
653  derivfinf = dzero
654  !
655  ! -- save wave states for resetting after iteration.
656  if (reset_state) then
657  call thiswork%wave_shift(this, 1, icell, 0, 1, this%nwavst(icell), 1)
658  end if
659  !
660  if (this%watab(icell) > this%celtop(icell)) &
661  this%watab(icell) = this%celtop(icell)
662  !
663  ! -- add water from mover to applied infiltration.
664  if (this%surflux(icell) > this%vks(icell)) then
665  this%surflux(icell) = this%vks(icell)
666  end if
667  !
668  ! -- saturation excess rejected infiltration
669  if (this%landflag(icell) == 1) then
670  call this%rejfinf(icell, deriv1, hgwf, trhsfinf, thcoffinf, finfact)
671  this%surflux(icell) = finfact
672  end if
673  !
674  ! -- calculate rejected infiltration
675  this%finf_rej(icell) = this%finf(icell) + &
676  (qfrommvr / this%uzfarea(icell)) - this%surflux(icell)
677  !
678  ! -- calculate groundwater discharge
679  if (iseepflag > 0 .and. this%landflag(icell) == 1) then
680  call this%gwseep(icell, deriv2, scale, hgwf, trhsseep, thcofseep, seep)
681  this%surfseep(icell) = seep
682  end if
683  !
684  ! -- route water through unsat zone, calc. storage change and recharge
685  test = this%watab(icell)
686  if (this%watabold(icell) - test < -dem15) test = this%watabold(icell)
687  if (this%celtop(icell) - test > dem15) then
688  if (issflag == 0) then
689  call this%routewaves(totfluxtot, delt, ietflag, icell, ierr)
690  if (ierr > 0) return
691  call this%uz_rise(icell, totfluxtot)
692  this%totflux(icell) = totfluxtot
693  if (this%ivertcon(icell) > 0) then
694  call this%addrech(icell, jbelow, hgwf, trhsfinf, thcoffinf, &
695  derivfinf, delt)
696  end if
697  else
698  this%totflux(icell) = this%surflux(icell) * delt
699  totfluxtot = this%surflux(icell) * delt
700  end if
701  thcoffinf = dzero
702  trhsfinf = this%totflux(icell) * this%uzfarea(icell) / delt
703  if (.not. reset_state) then
704  call this%update_wav(icell, delt, issflag, 0)
705  end if
706  else
707  this%totflux(icell) = this%surflux(icell) * delt
708  totfluxtot = this%surflux(icell) * delt
709  if (.not. reset_state) then
710  call this%update_wav(icell, delt, issflag, 1)
711  end if
712  end if
713  !
714  ! -- If formulating, then these variables will be present
715  if (present(deriv)) deriv = deriv1 + deriv2 + derivfinf
716  if (present(trhs)) trhs = trhsfinf + trhsseep
717  if (present(thcof)) thcof = thcoffinf + thcofseep
718  !
719  ! -- Assign water content prior to resetting waves
720  if (present(watercontent)) then
721  watercontent = this%get_wcnew(icell)
722  end if
723  !
724  ! -- reset waves to previous state for next iteration
725  if (reset_state) then
726  call this%wave_shift(thiswork, icell, 1, 0, 1, thiswork%nwavst(1), 1)
727  end if
728  !
729  ! -- Return
730  return

◆ trailwav()

subroutine uzfcellgroupmodule::trailwav ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
integer(i4b), intent(inout)  ierr 
)
private

Definition at line 1248 of file UzfCellGroup.f90.

1249  ! -- dummy
1250  class(UzfCellGroupType) :: this
1251  integer(I4B), intent(in) :: icell
1252  integer(I4B), intent(inout) :: ierr
1253  ! -- local
1254  real(DP) :: smoist, smoistinc, ftrail, eps_m1
1255  real(DP) :: thtsrinv
1256  real(DP) :: flux1, flux2, theta1, theta2
1257  real(DP) :: fnuminc
1258  integer(I4B) :: j, jj, jk, nwavstm1
1259  !
1260  ! -- initialize
1261  eps_m1 = dble(this%eps(icell)) - done
1262  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1263  nwavstm1 = this%nwavst(icell) - 1
1264  !
1265  ! -- initialize trailwaves
1266  smoist = (((this%surflux(icell) / this%vks(icell))** &
1267  (done / this%eps(icell))) * &
1268  (this%thts(icell) - this%thtr(icell))) + this%thtr(icell)
1269  if (this%uzthst(nwavstm1, icell) - smoist > dem9) then
1270  fnuminc = dzero
1271  do jk = 1, this%ntrail(icell)
1272  fnuminc = fnuminc + float(jk)
1273  end do
1274  smoistinc = (this%uzthst(nwavstm1, icell) - smoist) / (fnuminc - done)
1275  jj = this%ntrail(icell)
1276  ftrail = dble(this%ntrail(icell)) + done
1277  do j = this%nwavst(icell), this%nwavst(icell) + this%ntrail(icell) - 1
1278  if (j > this%nwav(icell)) then
1279  ! -- too many waves error
1280  ierr = 1
1281  return
1282  end if
1283  if (j > this%nwavst(icell)) then
1284  this%uzthst(j, icell) = this%uzthst(j - 1, icell) &
1285  - ((ftrail - float(jj)) * smoistinc)
1286  else
1287  this%uzthst(j, icell) = this%uzthst(j - 1, icell) - dem9
1288  end if
1289  jj = jj - 1
1290  if (this%uzthst(j, icell) <= this%thtr(icell) + dem9) &
1291  this%uzthst(j, icell) = this%thtr(icell) + dem9
1292  this%uzflst(j, icell) = &
1293  this%vks(icell) * (((this%uzthst(j, icell) - this%thtr(icell)) * &
1294  thtsrinv)**this%eps(icell))
1295  theta2 = this%uzthst(j - 1, icell)
1296  flux2 = this%uzflst(j - 1, icell)
1297  flux1 = this%uzflst(j, icell)
1298  theta1 = this%uzthst(j, icell)
1299  this%uzspst(j, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1300  this%thts(icell), this%thtr(icell), &
1301  this%eps(icell), this%vks(icell))
1302  this%uzdpst(j, icell) = dzero
1303  if (j == this%nwavst(icell)) then
1304  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1305  (this%ntrail(icell) + 1) * dem9
1306  else
1307  this%uzdpst(j, icell) = this%uzdpst(j - 1, icell) - dem9
1308  end if
1309  end do
1310  this%nwavst(icell) = this%nwavst(icell) + this%ntrail(icell) - 1
1311  if (this%nwavst(icell) >= this%nwav(icell)) then
1312  ! -- too many waves error
1313  ierr = 1
1314  return
1315  end if
1316  else
1317  this%uzdpst(this%nwavst, icell) = dzero
1318  this%uzflst(this%nwavst, icell) = &
1319  this%vks(icell) * (((this%uzthst(this%nwavst, icell) - &
1320  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1321  this%uzthst(this%nwavst, icell) = smoist
1322  theta2 = this%uzthst(this%nwavst(icell) - 1, icell)
1323  flux2 = this%uzflst(this%nwavst(icell) - 1, icell)
1324  flux1 = this%uzflst(this%nwavst(icell), icell)
1325  theta1 = this%uzthst(this%nwavst(icell), icell)
1326  this%uzspst(this%nwavst(icell), icell) = &
1327  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1328  this%thtr(icell), this%eps(icell), this%vks(icell))
1329  end if
1330  !
1331  ! -- Return
1332  return
Here is the call graph for this function:

◆ unsat_stor()

real(dp) function uzfcellgroupmodule::unsat_stor ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(inout)  d1 
)
private

Definition at line 1579 of file UzfCellGroup.f90.

1580  ! -- Return
1581  real(DP) :: unsat_stor
1582  ! -- dummy
1583  class(UzfCellGroupType) :: this
1584  integer(I4B), intent(in) :: icell
1585  real(DP), intent(inout) :: d1
1586  ! -- local
1587  real(DP) :: fm
1588  integer(I4B) :: j, k, nwavm1, jj
1589  !
1590  fm = dzero
1591  j = this%nwavst(icell) + 1
1592  k = this%nwavst(icell)
1593  nwavm1 = k - 1
1594  if (d1 > this%uzdpst(1, icell)) d1 = this%uzdpst(1, icell)
1595  !
1596  ! -- find deepest wave above depth d1, counter held as j
1597  do while (k > 0)
1598  if (this%uzdpst(k, icell) - d1 < -dem30) j = k
1599  k = k - 1
1600  end do
1601  if (j > this%nwavst(icell)) then
1602  fm = fm + (this%uzthst(this%nwavst(icell), icell) - this%thtr(icell)) * d1
1603  elseif (this%nwavst(icell) > 1) then
1604  if (j > 1) then
1605  fm = fm + (this%uzthst(j - 1, icell) - this%thtr(icell)) &
1606  * (d1 - this%uzdpst(j, icell))
1607  end if
1608  do jj = j, nwavm1
1609  fm = fm + (this%uzthst(jj, icell) - this%thtr(icell)) &
1610  * (this%uzdpst(jj, icell) &
1611  - this%uzdpst(jj + 1, icell))
1612  end do
1613  fm = fm + (this%uzthst(this%nwavst(icell), icell) - this%thtr(icell)) &
1614  * (this%uzdpst(this%nwavst(icell), icell))
1615  else
1616  fm = fm + (this%uzthst(1, icell) - this%thtr(icell)) * d1
1617  end if
1618  unsat_stor = fm
1619  !
1620  ! -- Return
1621  return

◆ update_wav()

subroutine uzfcellgroupmodule::update_wav ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  delt,
integer(i4b), intent(in)  iss,
integer(i4b), intent(in)  itest 
)
private

Definition at line 1626 of file UzfCellGroup.f90.

1627  ! -- dummy
1628  class(UzfCellGroupType) :: this
1629  integer(I4B), intent(in) :: icell
1630  integer(I4B), intent(in) :: itest
1631  integer(I4B), intent(in) :: iss
1632  real(DP), intent(in) :: delt
1633  ! -- local
1634  real(DP) :: bot, depthsave, top
1635  real(DP) :: thick, thtsrinv
1636  integer(I4B) :: nwavhld, k, j
1637  !
1638  bot = this%watab(icell)
1639  top = this%celtop(icell)
1640  thick = top - bot
1641  nwavhld = this%nwavst(icell)
1642  if (itest == 1) then
1643  this%uzflst(1, icell) = dzero
1644  this%uzthst(1, icell) = this%thtr(icell)
1645  return
1646  end if
1647  if (iss == 1) then
1648  if (this%thts(icell) - this%thtr(icell) < dem7) then
1649  thtsrinv = done / dem7
1650  else
1651  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1652  end if
1653  this%totflux(icell) = this%surflux(icell) * delt
1654  this%watabold(icell) = this%watab(icell)
1655  this%uzthst(1, icell) = this%thti(icell)
1656  this%uzflst(1, icell) = &
1657  this%vks(icell) * (((this%uzthst(1, icell) - this%thtr(icell)) &
1658  * thtsrinv)**this%eps(icell))
1659  this%uzdpst(1, icell) = thick
1660  this%uzspst(1, icell) = thick
1661  this%nwavst(icell) = 1
1662  else
1663  !
1664  ! -- water table rises through waves
1665  if (this%watab(icell) - this%watabold(icell) > dem30) then
1666  depthsave = this%uzdpst(1, icell)
1667  j = 0
1668  k = this%nwavst(icell)
1669  do while (k > 0)
1670  if (this%uzdpst(k, icell) - thick < -dem30) j = k
1671  k = k - 1
1672  end do
1673  this%uzdpst(1, icell) = thick
1674  if (j > 1) then
1675  this%uzspst(1, icell) = dzero
1676  this%nwavst(icell) = this%nwavst(icell) - j + 2
1677  this%uzthst(1, icell) = this%uzthst(j - 1, icell)
1678  this%uzflst(1, icell) = this%uzflst(j - 1, icell)
1679  if (j > 2) call this%wave_shift(this, icell, icell, j - 2, 2, &
1680  nwavhld - (j - 2), 1)
1681  elseif (j == 0) then
1682  this%uzspst(1, icell) = dzero
1683  this%uzthst(1, icell) = this%uzthst(this%nwavst(icell), icell)
1684  this%uzflst(1, icell) = this%uzflst(this%nwavst(icell), icell)
1685  this%nwavst(icell) = 1
1686  end if
1687  end if
1688  !
1689  ! -- calculate new unsat. storage
1690  if (thick <= dzero) then
1691  this%uzspst(1, icell) = dzero
1692  this%nwavst(icell) = 1
1693  this%uzthst(1, icell) = this%thtr(icell)
1694  this%uzflst(1, icell) = dzero
1695  end if
1696  this%watabold(icell) = this%watab(icell)
1697  end if
1698  !
1699  ! -- Return
1700  return

◆ uz_rise()

subroutine uzfcellgroupmodule::uz_rise ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(inout)  totfluxtot 
)
private

Definition at line 981 of file UzfCellGroup.f90.

982  ! -- dummy
983  class(UzfCellGroupType) :: this
984  integer(I4B), intent(in) :: icell
985  real(DP), intent(inout) :: totfluxtot
986  ! -- local
987  real(DP) :: fm1, fm2, d1
988  !
989  ! -- additional recharge from a rising water table
990  if (this%watab(icell) - this%watabold(icell) > dem30) then
991  d1 = this%celtop(icell) - this%watabold(icell)
992  fm1 = this%unsat_stor(icell, d1)
993  d1 = this%celtop(icell) - this%watab(icell)
994  fm2 = this%unsat_stor(icell, d1)
995  totfluxtot = totfluxtot + (fm1 - fm2)
996  end if
997  !
998  ! -- Return
999  return

◆ uzet()

subroutine uzfcellgroupmodule::uzet ( class(uzfcellgrouptype)  this,
integer(i4b), intent(in)  icell,
real(dp), intent(in)  delt,
integer(i4b), intent(in)  ietflag,
integer(i4b), intent(inout)  ierr 
)
private

Definition at line 1705 of file UzfCellGroup.f90.

1706  ! -- dummy
1707  class(UzfCellGroupType) :: this
1708  integer(I4B), intent(in) :: icell
1709  real(DP), intent(in) :: delt
1710  integer(I4B), intent(in) :: ietflag
1711  integer(I4B), intent(inout) :: ierr
1712  ! -- local
1713  type(UzfCellGroupType) :: uzfktemp
1714  real(DP) :: diff
1715  real(DP) :: thetaout
1716  real(DP) :: fm
1717  real(DP) :: st
1718  real(DP) :: thtsrinv
1719  real(DP) :: epsfksthts
1720  real(DP) :: fmp
1721  real(DP) :: fktho
1722  real(DP) :: theta1
1723  real(DP) :: theta2
1724  real(DP) :: flux1
1725  real(DP) :: flux2
1726  real(DP) :: hcap
1727  real(DP) :: factor
1728  real(DP) :: tho
1729  real(DP) :: depth
1730  real(DP) :: extwc1
1731  real(DP) :: petsub
1732  integer(I4B) :: i
1733  integer(I4B) :: j
1734  integer(I4B) :: jhold
1735  integer(I4B) :: jk
1736  integer(I4B) :: kj
1737  integer(I4B) :: kk
1738  integer(I4B) :: numadd
1739  integer(I4B) :: k
1740  integer(I4B) :: nwv
1741  integer(I4B) :: itest
1742  !
1743  ! -- initialize
1744  this%etact(icell) = dzero
1745  if (this%extdpuz(icell) < dem7) return
1746  petsub = this%rootact(icell) * this%pet(icell) * &
1747  this%extdpuz(icell) / this%extdp(icell)
1748  thetaout = delt * petsub / this%extdp(icell)
1749  if (ietflag == 1) thetaout = delt * this%pet(icell) / this%extdp(icell)
1750  if (thetaout < dem10) return
1751  depth = this%uzdpst(1, icell)
1752  st = this%unsat_stor(icell, depth)
1753  if (st < dem4) return
1754  !
1755  ! -- allocate temporary wave storage.
1756  nwv = this%nwavst(icell)
1757  itest = 0
1758  call uzfktemp%init(1, nwv)
1759  !
1760  ! store original wave characteristics
1761  call uzfktemp%wave_shift(this, 1, icell, 0, 1, nwv, 1)
1762  factor = done
1763  this%etact(icell) = dzero
1764  if (this%thts(icell) - this%thtr(icell) < dem7) then
1765  thtsrinv = 1.0 / dem7
1766  else
1767  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1768  end if
1769  epsfksthts = this%eps(icell) * this%vks(icell) * thtsrinv
1770  this%etact(icell) = dzero
1771  fmp = dzero
1772  extwc1 = this%extwc(icell) - this%thtr(icell)
1773  if (extwc1 < dem6) extwc1 = dem7
1774  numadd = 0
1775  fm = st
1776  k = 0
1777  !
1778  ! -- loop for reducing aet to pet when et is head dependent
1779  do while (itest == 0)
1780  k = k + 1
1781  if (k > 1 .AND. abs(fmp - petsub) > dem5 * petsub) then
1782  factor = factor / (fm / petsub)
1783  end if
1784  !
1785  ! -- one wave shallower than extdp
1786  if (this%nwavst(icell) == 1 .AND. &
1787  this%uzdpst(1, icell) <= this%extdpuz(icell)) then
1788  if (ietflag == 2) then
1789  tho = this%uzthst(1, icell)
1790  fktho = this%uzflst(1, icell)
1791  hcap = this%caph(icell, tho)
1792  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1793  end if
1794  if ((this%uzthst(1, icell) - thetaout) > this%thtr(icell) + extwc1) then
1795  this%uzthst(1, icell) = this%uzthst(1, icell) - thetaout
1796  this%uzflst(1, icell) = &
1797  this%vks(icell) * (((this%uzthst(1, icell) - &
1798  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1799  else if (this%uzthst(1, icell) > this%thtr(icell) + extwc1) then
1800  this%uzthst(1, icell) = this%thtr(icell) + extwc1
1801  this%uzflst(1, icell) = &
1802  this%vks(icell) * (((this%uzthst(1, icell) - &
1803  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1804  end if
1805  !
1806  ! -- all waves shallower than extinction depth
1807  else if (this%nwavst(icell) > 1 .AND. &
1808  this%uzdpst(this%nwavst(icell), icell) > this%extdpuz(icell)) then
1809  if (ietflag == 2) then
1810  tho = this%uzthst(this%nwavst(icell), icell)
1811  fktho = this%uzflst(this%nwavst(icell), icell)
1812  hcap = this%caph(icell, tho)
1813  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1814  end if
1815  if (this%uzthst(this%nwavst(icell), icell) - thetaout > &
1816  this%thtr(icell) + extwc1) then
1817  this%uzthst(this%nwavst(icell) + 1, icell) = &
1818  this%uzthst(this%nwavst(icell), icell) - thetaout
1819  numadd = 1
1820  else if (this%uzthst(this%nwavst(icell), icell) > &
1821  this%thtr(icell) + extwc1) then
1822  this%uzthst(this%nwavst(icell) + 1, icell) = this%thtr(icell) + extwc1
1823  numadd = 1
1824  end if
1825  if (numadd == 1) then
1826  this%uzflst(this%nwavst(icell) + 1, icell) = &
1827  this%vks(icell) * &
1828  (((this%uzthst(this%nwavst(icell) + 1, icell) - &
1829  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1830  theta2 = this%uzthst(this%nwavst(icell) + 1, icell)
1831  flux2 = this%uzflst(this%nwavst(icell) + 1, icell)
1832  flux1 = this%uzflst(this%nwavst(icell), icell)
1833  theta1 = this%uzthst(this%nwavst(icell), icell)
1834  this%uzspst(this%nwavst(icell) + 1, icell) = &
1835  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1836  this%thtr(icell), this%eps(icell), this%vks(icell))
1837  this%uzdpst(this%nwavst(icell) + 1, icell) = this%extdpuz(icell)
1838  this%nwavst(icell) = this%nwavst(icell) + 1
1839  if (this%nwavst(icell) > this%nwav(icell)) then
1840  !
1841  ! -- too many waves error, deallocate temp arrays and return
1842  ierr = 1
1843  goto 500
1844  end if
1845  else
1846  numadd = 0
1847  end if
1848  !
1849  ! -- one wave below extinction depth
1850  else if (this%nwavst(icell) == 1) then
1851  if (ietflag == 2) then
1852  tho = this%uzthst(1, icell)
1853  fktho = this%uzflst(1, icell)
1854  hcap = this%caph(icell, tho)
1855  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1856  end if
1857  if ((this%uzthst(1, icell) - thetaout) > this%thtr(icell) + extwc1) then
1858  if (thetaout > dem30) then
1859  this%uzthst(2, icell) = this%uzthst(1, icell) - thetaout
1860  this%uzflst(2, icell) = &
1861  this%vks(icell) * (((this%uzthst(2, icell) - this%thtr(icell)) * &
1862  thtsrinv)**this%eps(icell))
1863  this%uzdpst(2, icell) = this%extdpuz(icell)
1864  theta2 = this%uzthst(2, icell)
1865  flux2 = this%uzflst(2, icell)
1866  flux1 = this%uzflst(1, icell)
1867  theta1 = this%uzthst(1, icell)
1868  this%uzspst(2, icell) = &
1869  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1870  this%thtr(icell), this%eps(icell), this%vks(icell))
1871  this%nwavst(icell) = this%nwavst(icell) + 1
1872  if (this%nwavst(icell) > this%nwav(icell)) then
1873  !
1874  ! -- too many waves error
1875  ierr = 1
1876  goto 500
1877  end if
1878  end if
1879  else if (this%uzthst(1, icell) > this%thtr(icell) + extwc1) then
1880  if (thetaout > dem30) then
1881  this%uzthst(2, icell) = this%thtr(icell) + extwc1
1882  this%uzflst(2, icell) = &
1883  this%vks(icell) * (((this%uzthst(2, icell) - &
1884  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1885  this%uzdpst(2, icell) = this%extdpuz(icell)
1886  theta2 = this%uzthst(2, icell)
1887  flux2 = this%uzflst(2, icell)
1888  flux1 = this%uzflst(1, icell)
1889  theta1 = this%uzthst(1, icell)
1890  this%uzspst(2, icell) = &
1891  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1892  this%thtr(icell), this%eps(icell), this%vks(icell))
1893  this%nwavst(icell) = this%nwavst(icell) + 1
1894  if (this%nwavst(icell) > this%nwav(icell)) then
1895  !
1896  ! -- too many waves error
1897  ierr = 1
1898  goto 500
1899  end if
1900  end if
1901  end if
1902  else
1903  !
1904  ! -- extinction depth splits waves
1905  if (this%uzdpst(1, icell) - this%extdpuz(icell) > dem7) then
1906  j = 2
1907  jk = 0
1908  !
1909  ! -- locate extinction depth between waves
1910  do while (jk == 0)
1911  diff = this%uzdpst(j, icell) - this%extdpuz(icell)
1912  if (diff > dzero) then
1913  j = j + 1
1914  else
1915  jk = 1
1916  end if
1917  end do
1918  kk = j
1919  if (this%uzthst(j, icell) > this%thtr(icell) + extwc1) then
1920  !
1921  ! -- create a wave at extinction depth
1922  if (abs(diff) > dem5) then
1923  call this%wave_shift(this, icell, icell, -1, &
1924  this%nwavst(icell) + 1, j, -1)
1925  this%uzdpst(j, icell) = this%extdpuz(icell)
1926  this%nwavst(icell) = this%nwavst(icell) + 1
1927  if (this%nwavst(icell) > this%nwav(icell)) then
1928  !
1929  ! -- too many waves error
1930  ierr = 1
1931  goto 500
1932  end if
1933  end if
1934  kk = j
1935  else
1936  jhold = this%nwavst(icell)
1937  i = j + 1
1938  do while (i < this%nwavst(icell))
1939  if (this%uzthst(i, icell) > this%thtr(icell) + extwc1) then
1940  jhold = i
1941  i = this%nwavst(icell) + 1
1942  end if
1943  i = i + 1
1944  end do
1945  j = jhold
1946  kk = jhold
1947  end if
1948  else
1949  kk = 1
1950  end if
1951  !
1952  ! -- all waves above extinction depth
1953  do while (kk <= this%nwavst(icell))
1954  if (ietflag == 2) then
1955  tho = this%uzthst(kk, icell)
1956  fktho = this%uzflst(kk, icell)
1957  hcap = this%caph(icell, tho)
1958  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1959  end if
1960  if (this%uzthst(kk, icell) > this%thtr(icell) + extwc1) then
1961  if (this%uzthst(kk, icell) - thetaout > &
1962  this%thtr(icell) + extwc1) then
1963  this%uzthst(kk, icell) = this%uzthst(kk, icell) - thetaout
1964  else if (this%uzthst(kk, icell) > this%thtr(icell) + extwc1) then
1965  this%uzthst(kk, icell) = this%thtr(icell) + extwc1
1966  end if
1967  if (kk == 1) then
1968  this%uzflst(kk, icell) = &
1969  this%vks(icell) * &
1970  (((this%uzthst(kk, icell) - &
1971  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1972  end if
1973  if (kk > 1) then
1974  flux1 = &
1975  this%vks(icell) * ((this%uzthst(kk - 1, icell) - &
1976  this%thtr(icell)) * thtsrinv)**this%eps(icell)
1977  flux2 = &
1978  this%vks(icell) * ((this%uzthst(kk, icell) - &
1979  this%thtr(icell)) * thtsrinv)**this%eps(icell)
1980  this%uzflst(kk, icell) = flux2
1981  theta2 = this%uzthst(kk, icell)
1982  theta1 = this%uzthst(kk - 1, icell)
1983  this%uzspst(kk, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1984  this%thts(icell), &
1985  this%thtr(icell), &
1986  this%eps(icell), this%vks(icell))
1987  end if
1988  end if
1989  kk = kk + 1
1990  end do
1991  end if
1992  !
1993  ! -- calculate aet
1994  kj = 1
1995  do while (kj <= this%nwavst(icell) - 1)
1996  if (abs(this%uzthst(kj, icell) - this%uzthst(kj + 1, icell)) < dem6) then
1997  call this%wave_shift(this, icell, icell, 1, kj + 1, &
1998  this%nwavst(icell) - 1, 1)
1999  kj = kj - 1
2000  this%nwavst(icell) = this%nwavst(icell) - 1
2001  end if
2002  kj = kj + 1
2003  end do
2004  depth = this%uzdpst(1, icell)
2005  fm = this%unsat_stor(icell, depth)
2006  this%etact(icell) = st - fm
2007  fm = this%etact(icell) / delt
2008  if (this%etact(icell) < dzero) then
2009  call this%wave_shift(uzfktemp, icell, 1, 0, 1, nwv, 1)
2010  this%nwavst(icell) = nwv
2011  this%etact(icell) = dzero
2012  elseif (petsub - fm < -dem15 .AND. ietflag == 2) then
2013  !
2014  ! -- aet greater than pet, reset and try again
2015  call this%wave_shift(uzfktemp, icell, 1, 0, 1, nwv, 1)
2016  this%nwavst(icell) = nwv
2017  this%etact(icell) = dzero
2018  else
2019  itest = 1
2020  end if
2021  !
2022  ! -- end aet-pet loop for head dependent et
2023  fmp = fm
2024  if (k > 100) then
2025  itest = 1
2026  elseif (ietflag < 2) then
2027  fmp = petsub
2028  itest = 1
2029  end if
2030  end do
2031 500 continue
2032  !
2033  ! -- deallocate temporary worker
2034  call uzfktemp%dealloc()
2035  !
2036  ! -- Return
2037  return
Here is the call graph for this function:

◆ uzflow()

subroutine uzfcellgroupmodule::uzflow ( class(uzfcellgrouptype)  this,
real(dp), intent(inout)  thick,
real(dp), intent(inout)  thickold,
real(dp), intent(in)  delt,
integer(i4b), intent(in)  ietflag,
integer(i4b), intent(in)  icell,
integer(i4b), intent(inout)  ierr 
)
private

Definition at line 1125 of file UzfCellGroup.f90.

1126  ! -- dummy
1127  class(UzfCellGroupType) :: this
1128  real(DP), intent(inout) :: thickold
1129  real(DP), intent(inout) :: thick
1130  real(DP), intent(in) :: delt
1131  integer(I4B), intent(in) :: ietflag
1132  integer(I4B), intent(in) :: icell
1133  integer(I4B), intent(inout) :: ierr
1134  ! -- local
1135  real(DP) :: ffcheck, time, feps1, feps2
1136  real(DP) :: thetadif, thetab, fluxb, oldsflx
1137  integer(I4B) :: itrailflg, itester
1138  !
1139  time = dzero
1140  this%totflux(icell) = dzero
1141  itrailflg = 0
1142  oldsflx = this%uzflst(this%nwavst(icell), icell)
1143  call factors(feps1, feps2)
1144  !
1145  ! -- check for falling or rising water table
1146  if ((thick - thickold) > feps1) then
1147  thetadif = abs(this%uzthst(1, icell) - this%thtr(icell))
1148  if (thetadif > dem6) then
1149  call this%wave_shift(this, icell, icell, -1, &
1150  this%nwavst(icell) + 1, 2, -1)
1151  if (this%uzdpst(2, icell) < dem30) &
1152  this%uzdpst(2, icell) = (this%ntrail(icell) + dtwo) * dem6
1153  if (this%uzthst(2, icell) > this%thtr(icell)) then
1154  this%uzspst(2, icell) = this%uzflst(2, icell) / &
1155  (this%uzthst(2, icell) - this%thtr(icell))
1156  else
1157  this%uzspst(2, icell) = dzero
1158  end if
1159  this%uzthst(1, icell) = this%thtr(icell)
1160  this%uzflst(1, icell) = dzero
1161  this%uzspst(1, icell) = dzero
1162  this%uzdpst(1, icell) = thick
1163  this%nwavst(icell) = this%nwavst(icell) + 1
1164  if (this%nwavst(icell) >= this%nwav(icell)) then
1165  ! -- too many waves error
1166  ierr = 1
1167  return
1168  end if
1169  else
1170  this%uzdpst(1, icell) = thick
1171  end if
1172  end if
1173  thetab = this%uzthst(1, icell)
1174  fluxb = this%uzflst(1, icell)
1175  this%totflux(icell) = dzero
1176  itester = 0
1177  ffcheck = (this%surflux(icell) - this%uzflst(this%nwavst(icell), icell))
1178  !
1179  ! -- increase new waves in infiltration changes
1180  if (ffcheck > feps2 .OR. ffcheck < -feps2) then
1181  this%nwavst(icell) = this%nwavst(icell) + 1
1182  if (this%nwavst(icell) >= this%nwav(icell)) then
1183  !
1184  ! -- too many waves error
1185  ierr = 1
1186  return
1187  end if
1188  else if (this%nwavst(icell) == 1) then
1189  itester = 1
1190  end if
1191  if (this%nwavst(icell) > 1) then
1192  if (ffcheck < -feps2) then
1193  call this%trailwav(icell, ierr)
1194  if (ierr > 0) return
1195  itrailflg = 1
1196  end if
1197  call this%leadwav(time, itester, itrailflg, thetab, fluxb, ffcheck, &
1198  feps2, delt, icell)
1199  end if
1200  if (itester == 1) then
1201  this%totflux(icell) = this%totflux(icell) + &
1202  (delt - time) * this%uzflst(1, icell)
1203  time = dzero
1204  itester = 0
1205  end if
1206  !
1207  ! -- simulate et
1208  if (ietflag > 0) call this%uzet(icell, delt, ietflag, ierr)
1209  if (ierr > 0) return
1210  !
1211  ! -- Return
1212  return
Here is the call graph for this function:

◆ wave_shift()

subroutine uzfcellgroupmodule::wave_shift ( class(uzfcellgrouptype)  this,
type(uzfcellgrouptype)  this2,
integer(i4b), intent(in)  icell,
integer(i4b), intent(in)  icell2,
integer(i4b), intent(in)  shft,
integer(i4b), intent(in)  strt,
integer(i4b), intent(in)  stp,
integer(i4b), intent(in)  cntr 
)
private

Definition at line 1097 of file UzfCellGroup.f90.

1098  ! -- dummy
1099  class(UzfCellGroupType) :: this
1100  type(UzfCellGroupType) :: this2
1101  integer(I4B), intent(in) :: icell
1102  integer(I4B), intent(in) :: icell2
1103  integer(I4B), intent(in) :: shft
1104  integer(I4B), intent(in) :: strt
1105  integer(I4B), intent(in) :: stp
1106  integer(I4B), intent(in) :: cntr
1107  ! -- local
1108  integer(I4B) :: j
1109  !
1110  ! -- copy waves from one uzf cell group to another
1111  do j = strt, stp, cntr
1112  this%uzthst(j, icell) = this2%uzthst(j + shft, icell2)
1113  this%uzdpst(j, icell) = this2%uzdpst(j + shft, icell2)
1114  this%uzflst(j, icell) = this2%uzflst(j + shft, icell2)
1115  this%uzspst(j, icell) = this2%uzspst(j + shft, icell2)
1116  end do
1117  this%nwavst(icell) = this2%nwavst(icell2)
1118  !
1119  ! -- Return
1120  return