MODFLOW 6  version 6.6.0.dev0
USGS Modular Hydrologic Model
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_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 700 of file UzfCellGroup.f90.

701  ! -- dummy
702  class(UzfCellGroupType) :: this
703  integer(I4B), intent(in) :: icell
704  integer(I4B), intent(in) :: jbelow
705  real(DP), intent(inout) :: trhs
706  real(DP), intent(inout) :: thcof
707  real(DP), intent(inout) :: deriv
708  real(DP), intent(in) :: delt
709  real(DP), intent(in) :: hgwf
710  ! -- local
711  real(DP) :: fcheck
712  real(DP) :: x, scale, range
713  !
714  ! -- initialize
715  range = dem5
716  deriv = dzero
717  thcof = dzero
718  trhs = this%uzfarea(icell) * this%totflux(icell) / delt
719  if (this%totflux(icell) < dem14) return
720  scale = done
721  !
722  ! -- smoothly reduce flow between cells when head close to cell top
723  x = hgwf - (this%celbot(icell) - range)
724  call sscurve(x, range, deriv, scale)
725  deriv = this%uzfarea(icell) * deriv * this%totflux(icell) / delt
726  this%finf(jbelow) = (done - scale) * this%totflux(icell) / delt
727  fcheck = this%finf(jbelow) - this%vks(jbelow)
728  !
729  ! -- reduce flow between cells when vks is too small
730  if (fcheck < dem14) fcheck = dzero
731  this%finf(jbelow) = this%finf(jbelow) - fcheck
732  this%surfluxbelow(icell) = this%finf(jbelow)
733  this%totflux(icell) = scale * this%totflux(icell) + fcheck * delt
734  trhs = this%uzfarea(icell) * this%totflux(icell) / delt
Here is the call graph for this function:

◆ advance()

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

Definition at line 551 of file UzfCellGroup.f90.

552  ! -- dummy
553  class(UzfCellGroupType) :: this
554  integer(I4B), intent(in) :: icell
555  !
556  ! -- set variables
557  this%surfseep(icell) = dzero

◆ caph()

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

Definition at line 1898 of file UzfCellGroup.f90.

1899  ! -- dummy
1900  class(UzfCellGroupType) :: this
1901  integer(I4B), intent(in) :: icell
1902  real(DP), intent(in) :: tho
1903  ! -- local
1904  real(DP) :: caph, lambda, star
1905  !
1906  caph = -dem6
1907  star = (tho - this%thtr(icell)) / (this%thts(icell) - this%thtr(icell))
1908  if (star < dem15) star = dem15
1909  lambda = dtwo / (this%eps(icell) - dthree)
1910  if (star > dem15) then
1911  if (tho - this%thts(icell) < dem15) then
1912  caph = this%ha(icell) * star**(-done / lambda)
1913  else
1914  caph = dzero
1915  end if
1916  end if

◆ dealloc()

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

Definition at line 241 of file UzfCellGroup.f90.

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

◆ 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 844 of file UzfCellGroup.f90.

845  ! -- Return
846  real(DP) :: etfunc_nlin
847  ! -- dummy
848  real(DP), intent(in) :: s
849  real(DP), intent(in) :: x
850  real(DP), intent(in) :: c
851  real(DP), intent(inout) :: det
852  real(DP), intent(inout) :: trhs
853  real(DP), intent(inout) :: thcof
854  real(DP), intent(in) :: hgwf
855  ! -- local
856  real(DP) :: etgw
857  real(DP) :: range
858  real(DP) :: depth, scale
859  !
860  depth = hgwf - (s - x)
861  if (depth < dzero) depth = dzero
862  etgw = c
863  range = dem3 * x
864  call scubic(depth, range, det, scale)
865  etgw = etgw * scale
866  trhs = -etgw
867  det = -det * etgw
868  etfunc_nlin = etgw
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 1094 of file UzfCellGroup.f90.

1095  ! -- dummy
1096  real(DP), intent(out) :: feps1
1097  real(DP), intent(out) :: feps2
1098  real(DP) :: factor1
1099  real(DP) :: factor2
1100  !
1101  ! calculate constants for uzflow
1102  factor1 = done
1103  factor2 = done
1104  feps1 = dem9
1105  feps2 = dem9
1106  if (itmuni == 1) then
1107  factor1 = done / 86400.d0
1108  else if (itmuni == 2) then
1109  factor1 = done / 1440.d0
1110  else if (itmuni == 3) then
1111  factor1 = done / 24.0d0
1112  else if (itmuni == 5) then
1113  factor1 = 365.0d0
1114  end if
1115  factor2 = done / 0.3048
1116  feps1 = feps1 * factor1 * factor2
1117  feps2 = feps2 * factor1 * factor2
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 1937 of file UzfCellGroup.f90.

1938  ! -- dummy
1939  class(UzfCellGroupType) :: this
1940  integer(I4B), intent(in) :: icell !< uzf cell containing depth
1941  real(DP), intent(in) :: depth !< depth within the cell
1942  ! -- return
1943  real(DP) :: theta_at_depth
1944  ! -- local
1945  real(DP) :: d1
1946  real(DP) :: d2
1947  real(DP) :: f1
1948  real(DP) :: f2
1949  !
1950  if (this%watab(icell) < this%celtop(icell)) then
1951  if (this%celtop(icell) - depth > this%watab(icell)) then
1952  d1 = depth - dem3
1953  d2 = depth + dem3
1954  f1 = this%unsat_stor(icell, d1)
1955  f2 = this%unsat_stor(icell, d2)
1956  theta_at_depth = this%thtr(icell) + (f2 - f1) / (d2 - d1)
1957  else
1958  theta_at_depth = this%thts(icell)
1959  end if
1960  else
1961  theta_at_depth = this%thts(icell)
1962  end if

◆ 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 1967 of file UzfCellGroup.f90.

1968  ! -- dummy
1969  class(UzfCellGroupType) :: this
1970  integer(I4B), intent(in) :: icell !< uzf cell containing depth
1971  ! -- return
1972  real(DP) :: watercontent
1973  ! -- local
1974  real(DP) :: top
1975  real(DP) :: bot
1976  real(DP) :: theta_r
1977  real(DP) :: thk
1978  real(DP) :: hgwf
1979  real(DP) :: fm
1980  real(DP) :: d
1981  !
1982  hgwf = this%watab(icell)
1983  top = this%celtop(icell)
1984  bot = this%celbot(icell)
1985  thk = top - max(bot, hgwf)
1986  if (thk > dzero) then
1987  theta_r = this%thtr(icell)
1988  d = thk
1989  fm = this%unsat_stor(icell, d)
1990  watercontent = fm / thk
1991  watercontent = watercontent + theta_r
1992  else
1993  watercontent = dzero
1994  end if

◆ 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 767 of file UzfCellGroup.f90.

768  ! -- dummy
769  class(UzfCellGroupType) :: this
770  integer(I4B), intent(in) :: icell
771  real(DP), intent(inout) :: deriv
772  real(DP), intent(inout) :: trhs
773  real(DP), intent(inout) :: thcof
774  real(DP), intent(inout) :: seep
775  real(DP), intent(out) :: scale
776  real(DP), intent(in) :: hgwf
777  ! -- local
778  real(DP) :: x, range, y, deriv1, d1, d2, Q
779  !
780  seep = dzero
781  deriv = dzero
782  deriv1 = dzero
783  d1 = dzero
784  d2 = dzero
785  scale = dzero
786  q = this%uzfarea(icell) * this%vks(icell)
787  range = this%surfdep(icell)
788  x = hgwf - this%celtop(icell)
789  call scubiclinear(x, range, deriv1, y)
790  scale = y
791  seep = scale * q * (hgwf - this%celtop(icell)) / range
792  trhs = scale * q * this%celtop(icell) / range
793  thcof = -scale * q / range
794  d1 = -deriv1 * q * x / range
795  d2 = -scale * q / range
796  deriv = d1 + d2
797  if (seep < dzero) then
798  seep = dzero
799  deriv = dzero
800  trhs = dzero
801  thcof = dzero
802  end if
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 98 of file UzfCellGroup.f90.

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

◆ 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 1409 of file UzfCellGroup.f90.

1410  ! -- Return
1411  real(DP) :: leadspeed
1412  ! -- dummy
1413  real(DP), intent(in) :: theta1
1414  real(DP), intent(in) :: theta2
1415  real(DP), intent(in) :: flux1
1416  real(DP), intent(inout) :: flux2
1417  real(DP), intent(in) :: thts
1418  real(DP), intent(in) :: thtr
1419  real(DP), intent(in) :: eps
1420  real(DP), intent(in) :: vks
1421  ! -- local
1422  real(DP) :: comp1, comp2, thsrinv, epsfksths
1423  real(DP) :: eps_m1, fhold, comp3
1424  !
1425  eps_m1 = eps - done
1426  thsrinv = done / (thts - thtr)
1427  epsfksths = eps * vks * thsrinv
1428  comp1 = theta2 - theta1
1429  comp2 = abs(flux2 - flux1)
1430  comp3 = theta1 - thtr
1431  if (comp2 < dem15) flux2 = flux1 + dem15
1432  if (abs(comp1) < dem30) then
1433  if (comp3 > dem30) fhold = (comp3 * thsrinv)**eps
1434  if (fhold < dem30) fhold = dem30
1435  leadspeed = epsfksths * (fhold**eps_m1)
1436  else
1437  leadspeed = (flux2 - flux1) / (theta2 - theta1)
1438  end if
1439  if (leadspeed < dem30) leadspeed = dem30
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 1208 of file UzfCellGroup.f90.

1210  ! -- dummy
1211  class(UzfCellGroupType) :: this
1212  real(DP), intent(inout) :: thetab
1213  real(DP), intent(inout) :: fluxb
1214  real(DP), intent(in) :: feps2
1215  real(DP), intent(inout) :: time
1216  integer(I4B), intent(inout) :: itester
1217  integer(I4B), intent(inout) :: itrailflg
1218  real(DP), intent(inout) :: ffcheck
1219  real(DP), intent(in) :: delt
1220  integer(I4B), intent(in) :: icell
1221  ! -- local
1222  real(DP) :: bottomtime, shortest, fcheck
1223  real(DP) :: eps_m1, timenew, bottom, timedt
1224  real(DP) :: thtsrinv, diff, fluxhld2
1225  real(DP) :: flux1, flux2, theta1, theta2, ftest
1226  real(DP), allocatable, dimension(:) :: checktime
1227  integer(I4B) :: iflx, iremove, j, l
1228  integer(I4B) :: nwavp1, jshort
1229  integer(I4B), allocatable, dimension(:) :: more
1230  !
1231  allocate (checktime(this%nwavst(icell)))
1232  allocate (more(this%nwavst(icell)))
1233  ftest = dzero
1234  eps_m1 = dble(this%eps(icell)) - done
1235  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1236  !
1237  ! -- initialize new wave
1238  if (itrailflg == 0) then
1239  if (ffcheck > feps2) then
1240  this%uzflst(this%nwavst(icell), icell) = this%surflux(icell)
1241  if (this%uzflst(this%nwavst(icell), icell) < dem30) &
1242  this%uzflst(this%nwavst(icell), icell) = dzero
1243  this%uzthst(this%nwavst(icell), icell) = &
1244  (((this%uzflst(this%nwavst(icell), icell) / this%vks(icell))** &
1245  (done / this%eps(icell))) * (this%thts(icell) - this%thtr(icell))) &
1246  + this%thtr(icell)
1247  theta2 = this%uzthst(this%nwavst(icell), icell)
1248  flux2 = this%uzflst(this%nwavst(icell), icell)
1249  flux1 = this%uzflst(this%nwavst(icell) - 1, icell)
1250  theta1 = this%uzthst(this%nwavst(icell) - 1, icell)
1251  this%uzspst(this%nwavst(icell), icell) = &
1252  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1253  this%thtr(icell), this%eps(icell), this%vks(icell))
1254  this%uzdpst(this%nwavst(icell), icell) = dzero
1255  end if
1256  end if
1257  !
1258  ! -- route all waves and interception of waves over times step
1259  diff = done
1260  timedt = dzero
1261  iflx = 0
1262  fluxhld2 = this%uzflst(1, icell)
1263  if (this%nwavst(icell) == 0) itester = 1
1264  if (itester /= 1) then
1265  do while (diff > dem6)
1266  nwavp1 = this%nwavst(icell) + 1
1267  timedt = delt - time
1268  do j = 1, this%nwavst(icell)
1269  checktime(j) = dep20
1270  more(j) = 0
1271  end do
1272  shortest = timedt
1273  if (this%nwavst(icell) > 2) then
1274  j = 2
1275  !
1276  ! -- calculate time until wave overtakes wave ahead
1277  nwavp1 = this%nwavst(icell) + 1
1278  do while (j < nwavp1)
1279  ftest = this%uzspst(j - 1, icell) - this%uzspst(j, icell)
1280  if (abs(ftest) > dem30) then
1281  checktime(j) = (this%uzdpst(j, icell) - &
1282  this%uzdpst(j - 1, icell)) / ftest
1283  if (checktime(j) < dem30) checktime(j) = dep20
1284  end if
1285  j = j + 1
1286  end do
1287  end if
1288  !
1289  ! - calc time until wave reaches bottom of cell
1290  bottomtime = dep20
1291  if (this%nwavst(icell) > 1) then
1292  if (this%uzspst(2, icell) > dzero) then
1293  bottom = this%uzspst(2, icell)
1294  if (bottom < dem15) bottom = dem15
1295  bottomtime = (this%uzdpst(1, icell) - this%uzdpst(2, icell)) / bottom
1296  if (bottomtime < dzero) bottomtime = dem12
1297  end if
1298  end if
1299  !
1300  ! -- calc time for wave interception
1301  jshort = 0
1302  do j = this%nwavst(icell), 3, -1
1303  if (shortest - checktime(j) > -dem9) then
1304  more(j) = 1
1305  jshort = j
1306  shortest = checktime(j)
1307  end if
1308  end do
1309  do j = 3, this%nwavst(icell)
1310  if (shortest - checktime(j) < dem9) then
1311  if (j /= jshort) more(j) = 0
1312  end if
1313  end do
1314  !
1315  ! -- what happens first, waves hits bottom or interception
1316  iremove = 0
1317  timenew = time
1318  fcheck = (time + shortest) - delt
1319  if (shortest < dem7) fcheck = -done
1320  if (bottomtime < shortest .AND. time + bottomtime < delt) then
1321  j = 2
1322  do while (j < nwavp1)
1323  !
1324  ! -- route waves
1325  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1326  this%uzspst(j, icell) * bottomtime
1327  j = j + 1
1328  end do
1329  fluxb = this%uzflst(2, icell)
1330  thetab = this%uzthst(2, icell)
1331  iflx = 1
1332  call this%wave_shift(this, icell, icell, 1, 1, &
1333  this%nwavst(icell) - 1, 1)
1334  iremove = 1
1335  timenew = time + bottomtime
1336  this%uzspst(1, icell) = dzero
1337  !
1338  ! -- do waves intercept before end of time step
1339  else if (fcheck < dzero .AND. this%nwavst(icell) > 2) then
1340  j = 2
1341  do while (j < nwavp1)
1342  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1343  this%uzspst(j, icell) * shortest
1344  j = j + 1
1345  end do
1346  !
1347  ! -- combine waves that intercept, remove a wave
1348  j = 3
1349  l = j
1350  do while (j < this%nwavst(icell) + 1)
1351  if (more(j) == 1) then
1352  l = j
1353  theta2 = this%uzthst(j, icell)
1354  flux2 = this%uzflst(j, icell)
1355  if (j == 3) then
1356  flux1 = fluxb
1357  theta1 = thetab
1358  else
1359  flux1 = this%uzflst(j - 2, icell)
1360  theta1 = this%uzthst(j - 2, icell)
1361  end if
1362  this%uzspst(j, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1363  this%thts(icell), &
1364  this%thtr(icell), &
1365  this%eps(icell), this%vks(icell))
1366  !
1367  ! -- update waves
1368  call this%wave_shift(this, icell, icell, 1, l - 1, &
1369  this%nwavst(icell) - 1, 1)
1370  l = this%nwavst(icell) + 1
1371  iremove = iremove + 1
1372  end if
1373  j = j + 1
1374  end do
1375  timenew = timenew + shortest
1376  !
1377  ! -- calc. total flux to bottom during remaining time in step
1378  else
1379  j = 2
1380  do while (j < nwavp1)
1381  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1382  this%uzspst(j, icell) * timedt
1383  j = j + 1
1384  end do
1385  timenew = delt
1386  end if
1387  this%totflux(icell) = this%totflux(icell) + fluxhld2 * (timenew - time)
1388  if (iflx == 1) then
1389  fluxhld2 = this%uzflst(1, icell)
1390  iflx = 0
1391  end if
1392  !
1393  ! -- remove dead waves
1394  this%nwavst(icell) = this%nwavst(icell) - iremove
1395  time = timenew
1396  diff = delt - time
1397  if (this%nwavst(icell) == 1) then
1398  itester = 1
1399  exit
1400  end if
1401  end do
1402  end if
1403  deallocate (checktime)
1404  deallocate (more)
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 1920 of file UzfCellGroup.f90.

1921  ! -- Return
1922  real(DP) :: rate_et_z
1923  ! -- dummy
1924  class(UzfCellGroupType) :: this
1925  integer(I4B), intent(in) :: icell
1926  real(DP), intent(in) :: factor, fktho, h
1927  !
1928  rate_et_z = factor * fktho * (h - this%hroot(icell))
1929  if (rate_et_z < dzero) rate_et_z = dzero

◆ 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 739 of file UzfCellGroup.f90.

740  ! -- dummy
741  class(UzfCellGroupType) :: this
742  integer(I4B), intent(in) :: icell
743  real(DP), intent(inout) :: deriv
744  real(DP), intent(inout) :: finfact
745  real(DP), intent(inout) :: thcof
746  real(DP), intent(inout) :: trhs
747  real(DP), intent(in) :: hgwf
748  ! -- local
749  real(DP) :: x, range, scale, q
750  !
751  range = this%surfdep(icell)
752  q = this%surflux(icell)
753  finfact = q
754  trhs = finfact * this%uzfarea(icell)
755  x = this%celtop(icell) - hgwf
756  call slinear(x, range, deriv, scale)
757  deriv = -q * deriv * this%uzfarea(icell) * scale
758  if (scale < done) then
759  finfact = q * scale
760  trhs = finfact * this%uzfarea(icell) * this%celtop(icell) / range
761  thcof = finfact * this%uzfarea(icell) / range
762  end if
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 942 of file UzfCellGroup.f90.

943  ! -- dummy
944  class(UzfCellGroupType) :: this
945  real(DP), intent(inout) :: totfluxtot
946  real(DP), intent(in) :: delt
947  integer(I4B), intent(in) :: ietflag
948  integer(I4B), intent(in) :: icell
949  integer(I4B), intent(inout) :: ierr
950  ! -- local
951  real(DP) :: thick, thickold
952  integer(I4B) :: idelt, iwav, ik
953  !
954  ! -- initialize
955  this%totflux(icell) = dzero
956  this%etact(icell) = dzero
957  thick = this%celtop(icell) - this%watab(icell)
958  thickold = this%celtop(icell) - this%watabold(icell)
959  !
960  ! -- no uz, clear waves
961  if (thickold < dzero) then
962  do iwav = 1, 6
963  this%uzthst(iwav, icell) = this%thtr(icell)
964  this%uzdpst(iwav, icell) = dzero
965  this%uzspst(iwav, icell) = dzero
966  this%uzflst(iwav, icell) = dzero
967  this%nwavst(icell) = 1
968  end do
969  end if
970  idelt = 1
971  do ik = 1, idelt
972  call this%uzflow(thick, thickold, delt, ietflag, icell, ierr)
973  if (ierr > 0) return
974  totfluxtot = totfluxtot + this%totflux(icell)
975  end do

◆ setbelowpet()

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

Definition at line 491 of file UzfCellGroup.f90.

492  ! -- modules
493  use tdismodule, only: delt
494  ! -- dummy
495  class(UzfCellGroupType) :: this
496  integer(I4B), intent(in) :: icell
497  integer(I4B), intent(in) :: jbelow
498  ! -- dummy
499  real(DP) :: pet
500  !
501  pet = dzero
502  !
503  ! -- transfer unmet pet to lower cell
504  !
505  if (this%extdpuz(jbelow) > dem3) then
506  pet = this%pet(icell) - this%etact(icell) / delt - &
507  this%gwet(icell) / this%uzfarea(icell)
508  if (pet < dzero) pet = dzero
509  end if
510  this%pet(jbelow) = pet
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 335 of file UzfCellGroup.f90.

337  ! -- dummy
338  class(UzfCellGroupType) :: this
339  integer(I4B), intent(in) :: icell
340  real(DP), intent(in) :: area
341  real(DP), intent(in) :: top
342  real(DP), intent(in) :: bot
343  real(DP), intent(in) :: surfdep
344  real(DP), intent(in) :: vks
345  real(DP), intent(in) :: thtr
346  real(DP), intent(in) :: thts
347  real(DP), intent(in) :: thti
348  real(DP), intent(in) :: eps
349  integer(I4B), intent(in) :: ntrail
350  integer(I4B), intent(in) :: landflag
351  integer(I4B), intent(in) :: ivertcon
352  !
353  ! -- set the values for uzf cell icell
354  this%landflag(icell) = landflag
355  this%ivertcon(icell) = ivertcon
356  this%surfdep(icell) = surfdep
357  this%uzfarea(icell) = area
358  this%cellarea(icell) = area
359  if (this%landflag(icell) == 1) then
360  this%celtop(icell) = top - dhalf * this%surfdep(icell)
361  else
362  this%celtop(icell) = top
363  end if
364  this%celbot(icell) = bot
365  this%vks(icell) = vks
366  this%thtr(icell) = thtr
367  this%thts(icell) = thts
368  this%thti(icell) = thti
369  this%eps(icell) = eps
370  this%ntrail(icell) = ntrail
371  this%pet(icell) = dzero
372  this%extdp(icell) = dzero
373  this%extwc(icell) = dzero
374  this%ha(icell) = dzero
375  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 428 of file UzfCellGroup.f90.

429  ! -- dummy
430  class(UzfCellGroupType) :: this
431  integer(I4B), intent(in) :: icell
432  integer(I4B), intent(in) :: jbelow
433  real(DP), intent(in) :: pet
434  real(DP), intent(in) :: extdp
435  ! -- local
436  real(DP) :: thick
437  !
438  if (this%landflag(icell) == 1) then
439  this%pet(icell) = pet
440  this%gwpet(icell) = pet
441  else
442  this%pet(icell) = dzero
443  this%gwpet(icell) = dzero
444  end if
445  thick = this%celtop(icell) - this%celbot(icell)
446  this%extdp(icell) = extdp
447  if (this%landflag(icell) > 0) then
448  this%landtop(icell) = this%celtop(icell)
449  this%petmax(icell) = this%pet(icell)
450  end if
451  !
452  ! -- set uz extinction depth
453  if (this%landtop(icell) - this%extdp(icell) < this%celbot(icell)) then
454  this%extdpuz(icell) = thick
455  else
456  this%extdpuz(icell) = this%celtop(icell) - &
457  (this%landtop(icell) - this%extdp(icell))
458  end if
459  if (this%extdpuz(icell) < dzero) this%extdpuz(icell) = dzero
460  if (this%extdpuz(icell) > dem7 .and. this%extdp(icell) < dem7) &
461  this%extdp(icell) = this%extdpuz(icell)
462  !
463  ! -- set pet for underlying cell
464  if (jbelow > 0) then
465  this%landtop(jbelow) = this%landtop(icell)
466  this%petmax(jbelow) = this%petmax(icell)
467  end if

◆ 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 529 of file UzfCellGroup.f90.

530  ! -- dummy
531  class(UzfCellGroupType) :: this
532  integer(I4B), intent(in) :: icell
533  integer(I4B), intent(in) :: jbelow
534  real(DP), intent(in) :: ha
535  real(DP), intent(in) :: hroot
536  real(DP), intent(in) :: rootact
537  !
538  ! -- set variables
539  this%ha(icell) = ha
540  this%hroot(icell) = hroot
541  this%rootact(icell) = rootact
542  if (jbelow > 0) then
543  this%ha(jbelow) = ha
544  this%hroot(jbelow) = hroot
545  this%rootact(jbelow) = rootact
546  end if

◆ 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 515 of file UzfCellGroup.f90.

516  ! -- dummy
517  class(UzfCellGroupType) :: this
518  integer(I4B), intent(in) :: icell
519  integer(I4B), intent(in) :: jbelow
520  real(DP), intent(in) :: extwc
521  !
522  ! -- set extinction water content
523  this%extwc(icell) = extwc
524  if (jbelow > 0) this%extwc(jbelow) = extwc

◆ setdatafinf()

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

Definition at line 396 of file UzfCellGroup.f90.

397  ! -- dummy
398  class(UzfCellGroupType) :: this
399  integer(I4B), intent(in) :: icell
400  real(DP), intent(in) :: finf
401  !
402  if (this%landflag(icell) == 1) then
403  this%sinf(icell) = finf
404  this%finf(icell) = finf
405  else
406  this%sinf(icell) = dzero
407  this%finf(icell) = dzero
408  end if
409  this%finf_rej(icell) = dzero
410  this%surflux(icell) = dzero
411  this%surfluxbelow(icell) = dzero

◆ setdatauzfarea()

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

Definition at line 416 of file UzfCellGroup.f90.

417  ! -- dummy
418  class(UzfCellGroupType) :: this
419  integer(I4B), intent(in) :: icell
420  real(DP), intent(in) :: areamult
421  !
422  ! -- set uzf area
423  this%uzfarea(icell) = this%cellarea(icell) * areamult

◆ setgwpet()

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

Definition at line 472 of file UzfCellGroup.f90.

473  ! -- modules
474  use tdismodule, only: delt
475  ! -- dummy
476  class(UzfCellGroupType) :: this
477  integer(I4B), intent(in) :: icell
478  ! -- dummy
479  real(DP) :: pet
480  !
481  pet = dzero
482  !
483  ! -- reduce pet for gw by uzet
484  pet = this%pet(icell) - this%etact(icell) / delt
485  if (pet < dzero) pet = dzero
486  this%gwpet(icell) = pet

◆ sethead()

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

Definition at line 380 of file UzfCellGroup.f90.

381  ! -- dummy
382  class(UzfCellGroupType) :: this
383  integer(I4B), intent(in) :: icell
384  real(DP), intent(in) :: hgwf
385  !
386  ! -- set initial head
387  this%watab(icell) = this%celbot(icell)
388  if (hgwf > this%celbot(icell)) this%watab(icell) = hgwf
389  if (this%watab(icell) > this%celtop(icell)) &
390  this%watab(icell) = this%celtop(icell)
391  this%watabold(icell) = this%watab(icell)

◆ setwaves()

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

Definition at line 893 of file UzfCellGroup.f90.

894  ! -- dummy
895  class(UzfCellGroupType) :: this
896  ! -- local
897  integer(I4B), intent(in) :: icell
898  real(DP) :: bottom, top
899  integer(I4B) :: jk
900  real(DP) :: thick
901  !
902  ! -- initialize
903  this%totflux(icell) = dzero
904  this%nwavst(icell) = 1
905  this%uzdpst(:, icell) = dzero
906  thick = this%celtop(icell) - this%watab(icell)
907  do jk = 1, this%nwav(icell)
908  this%uzthst(jk, icell) = this%thtr(icell)
909  end do
910  !
911  ! -- initialize waves for first stress period
912  if (thick > dzero) then
913  this%uzdpst(1, icell) = thick
914  this%uzthst(1, icell) = this%thti(icell)
915  top = this%uzthst(1, icell) - this%thtr(icell)
916  if (top < dzero) top = dzero
917  bottom = this%thts(icell) - this%thtr(icell)
918  if (bottom < dzero) bottom = dzero
919  this%uzflst(1, icell) = this%vks(icell) * (top / bottom)**this%eps(icell)
920  if (this%uzthst(1, icell) < this%thtr(icell)) &
921  this%uzthst(1, icell) = this%thtr(icell)
922  !
923  ! -- calculate water stored in the unsaturated zone
924  if (top > dzero) then
925  this%uzspst(1, icell) = dzero
926  else
927  this%uzflst(1, icell) = dzero
928  this%uzspst(1, icell) = dzero
929  end if
930  !
931  ! no unsaturated zone
932  else
933  this%uzflst(1, icell) = dzero
934  this%uzdpst(1, icell) = dzero
935  this%uzspst(1, icell) = dzero
936  this%uzthst(1, icell) = this%thtr(icell)
937  end if

◆ 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 807 of file UzfCellGroup.f90.

808  ! -- dummy
809  class(UzfCellGroupType) :: this
810  integer(I4B), intent(in) :: igwetflag
811  integer(I4B), intent(in) :: icell
812  real(DP), intent(in) :: hgwf
813  real(DP), intent(inout) :: trhs
814  real(DP), intent(inout) :: thcof
815  real(DP), intent(inout) :: det
816  ! -- local
817  real(DP) :: s, x, c, b, et
818  !
819  this%gwet(icell) = dzero
820  trhs = dzero
821  thcof = dzero
822  det = dzero
823  s = this%landtop(icell)
824  x = this%extdp(icell)
825  c = this%gwpet(icell)
826  b = this%celbot(icell)
827  if (b > hgwf) return
828  if (x < dem6) return
829  if (igwetflag == 1) then
830  et = etfunc_lin(s, x, c, det, trhs, thcof, hgwf, &
831  this%celtop(icell), this%celbot(icell))
832  else if (igwetflag == 2) then
833  et = etfunc_nlin(s, x, c, det, trhs, thcof, hgwf)
834  end if
835  ! this%gwet(icell) = et * this%uzfarea(icell)
836  trhs = trhs * this%uzfarea(icell)
837  thcof = thcof * this%uzfarea(icell)
838  this%gwet(icell) = trhs - (thcof * hgwf)
839  ! write(99,*)'in group', icell, this%gwet(icell)
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 563 of file UzfCellGroup.f90.

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

◆ trailwav()

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

Definition at line 1122 of file UzfCellGroup.f90.

1123  ! -- dummy
1124  class(UzfCellGroupType) :: this
1125  integer(I4B), intent(in) :: icell
1126  integer(I4B), intent(inout) :: ierr
1127  ! -- local
1128  real(DP) :: smoist, smoistinc, ftrail, eps_m1
1129  real(DP) :: thtsrinv
1130  real(DP) :: flux1, flux2, theta1, theta2
1131  real(DP) :: fnuminc
1132  integer(I4B) :: j, jj, jk, nwavstm1
1133  !
1134  ! -- initialize
1135  eps_m1 = dble(this%eps(icell)) - done
1136  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1137  nwavstm1 = this%nwavst(icell) - 1
1138  !
1139  ! -- initialize trailwaves
1140  smoist = (((this%surflux(icell) / this%vks(icell))** &
1141  (done / this%eps(icell))) * &
1142  (this%thts(icell) - this%thtr(icell))) + this%thtr(icell)
1143  if (this%uzthst(nwavstm1, icell) - smoist > dem9) then
1144  fnuminc = dzero
1145  do jk = 1, this%ntrail(icell)
1146  fnuminc = fnuminc + float(jk)
1147  end do
1148  smoistinc = (this%uzthst(nwavstm1, icell) - smoist) / (fnuminc - done)
1149  jj = this%ntrail(icell)
1150  ftrail = dble(this%ntrail(icell)) + done
1151  do j = this%nwavst(icell), this%nwavst(icell) + this%ntrail(icell) - 1
1152  if (j > this%nwav(icell)) then
1153  ! -- too many waves error
1154  ierr = 1
1155  return
1156  end if
1157  if (j > this%nwavst(icell)) then
1158  this%uzthst(j, icell) = this%uzthst(j - 1, icell) &
1159  - ((ftrail - float(jj)) * smoistinc)
1160  else
1161  this%uzthst(j, icell) = this%uzthst(j - 1, icell) - dem9
1162  end if
1163  jj = jj - 1
1164  if (this%uzthst(j, icell) <= this%thtr(icell) + dem9) &
1165  this%uzthst(j, icell) = this%thtr(icell) + dem9
1166  this%uzflst(j, icell) = &
1167  this%vks(icell) * (((this%uzthst(j, icell) - this%thtr(icell)) * &
1168  thtsrinv)**this%eps(icell))
1169  theta2 = this%uzthst(j - 1, icell)
1170  flux2 = this%uzflst(j - 1, icell)
1171  flux1 = this%uzflst(j, icell)
1172  theta1 = this%uzthst(j, icell)
1173  this%uzspst(j, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1174  this%thts(icell), this%thtr(icell), &
1175  this%eps(icell), this%vks(icell))
1176  this%uzdpst(j, icell) = dzero
1177  if (j == this%nwavst(icell)) then
1178  this%uzdpst(j, icell) = this%uzdpst(j, icell) + &
1179  (this%ntrail(icell) + 1) * dem9
1180  else
1181  this%uzdpst(j, icell) = this%uzdpst(j - 1, icell) - dem9
1182  end if
1183  end do
1184  this%nwavst(icell) = this%nwavst(icell) + this%ntrail(icell) - 1
1185  if (this%nwavst(icell) >= this%nwav(icell)) then
1186  ! -- too many waves error
1187  ierr = 1
1188  return
1189  end if
1190  else
1191  this%uzdpst(this%nwavst, icell) = dzero
1192  this%uzflst(this%nwavst, icell) = &
1193  this%vks(icell) * (((this%uzthst(this%nwavst, icell) - &
1194  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1195  this%uzthst(this%nwavst, icell) = smoist
1196  theta2 = this%uzthst(this%nwavst(icell) - 1, icell)
1197  flux2 = this%uzflst(this%nwavst(icell) - 1, icell)
1198  flux1 = this%uzflst(this%nwavst(icell), icell)
1199  theta1 = this%uzthst(this%nwavst(icell), icell)
1200  this%uzspst(this%nwavst(icell), icell) = &
1201  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1202  this%thtr(icell), this%eps(icell), this%vks(icell))
1203  end if
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 1444 of file UzfCellGroup.f90.

1445  ! -- Return
1446  real(DP) :: unsat_stor
1447  ! -- dummy
1448  class(UzfCellGroupType) :: this
1449  integer(I4B), intent(in) :: icell
1450  real(DP), intent(inout) :: d1
1451  ! -- local
1452  real(DP) :: fm
1453  integer(I4B) :: j, k, nwavm1, jj
1454  !
1455  fm = dzero
1456  j = this%nwavst(icell) + 1
1457  k = this%nwavst(icell)
1458  nwavm1 = k - 1
1459  if (d1 > this%uzdpst(1, icell)) d1 = this%uzdpst(1, icell)
1460  !
1461  ! -- find deepest wave above depth d1, counter held as j
1462  do while (k > 0)
1463  if (this%uzdpst(k, icell) - d1 < -dem30) j = k
1464  k = k - 1
1465  end do
1466  if (j > this%nwavst(icell)) then
1467  fm = fm + (this%uzthst(this%nwavst(icell), icell) - this%thtr(icell)) * d1
1468  elseif (this%nwavst(icell) > 1) then
1469  if (j > 1) then
1470  fm = fm + (this%uzthst(j - 1, icell) - this%thtr(icell)) &
1471  * (d1 - this%uzdpst(j, icell))
1472  end if
1473  do jj = j, nwavm1
1474  fm = fm + (this%uzthst(jj, icell) - this%thtr(icell)) &
1475  * (this%uzdpst(jj, icell) &
1476  - this%uzdpst(jj + 1, icell))
1477  end do
1478  fm = fm + (this%uzthst(this%nwavst(icell), icell) - this%thtr(icell)) &
1479  * (this%uzdpst(this%nwavst(icell), icell))
1480  else
1481  fm = fm + (this%uzthst(1, icell) - this%thtr(icell)) * d1
1482  end if
1483  unsat_stor = fm

◆ 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 1488 of file UzfCellGroup.f90.

1489  ! -- dummy
1490  class(UzfCellGroupType) :: this
1491  integer(I4B), intent(in) :: icell
1492  integer(I4B), intent(in) :: itest
1493  integer(I4B), intent(in) :: iss
1494  real(DP), intent(in) :: delt
1495  ! -- local
1496  real(DP) :: bot, depthsave, top
1497  real(DP) :: thick, thtsrinv
1498  integer(I4B) :: nwavhld, k, j
1499  !
1500  bot = this%watab(icell)
1501  top = this%celtop(icell)
1502  thick = top - bot
1503  nwavhld = this%nwavst(icell)
1504  if (itest == 1) then
1505  this%uzflst(1, icell) = dzero
1506  this%uzthst(1, icell) = this%thtr(icell)
1507  return
1508  end if
1509  if (iss == 1) then
1510  if (this%thts(icell) - this%thtr(icell) < dem7) then
1511  thtsrinv = done / dem7
1512  else
1513  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1514  end if
1515  this%totflux(icell) = this%surflux(icell) * delt
1516  this%watabold(icell) = this%watab(icell)
1517  this%uzthst(1, icell) = this%thti(icell)
1518  this%uzflst(1, icell) = &
1519  this%vks(icell) * (((this%uzthst(1, icell) - this%thtr(icell)) &
1520  * thtsrinv)**this%eps(icell))
1521  this%uzdpst(1, icell) = thick
1522  this%uzspst(1, icell) = thick
1523  this%nwavst(icell) = 1
1524  else
1525  !
1526  ! -- water table rises through waves
1527  if (this%watab(icell) - this%watabold(icell) > dem30) then
1528  depthsave = this%uzdpst(1, icell)
1529  j = 0
1530  k = this%nwavst(icell)
1531  do while (k > 0)
1532  if (this%uzdpst(k, icell) - thick < -dem30) j = k
1533  k = k - 1
1534  end do
1535  this%uzdpst(1, icell) = thick
1536  if (j > 1) then
1537  this%uzspst(1, icell) = dzero
1538  this%nwavst(icell) = this%nwavst(icell) - j + 2
1539  this%uzthst(1, icell) = this%uzthst(j - 1, icell)
1540  this%uzflst(1, icell) = this%uzflst(j - 1, icell)
1541  if (j > 2) call this%wave_shift(this, icell, icell, j - 2, 2, &
1542  nwavhld - (j - 2), 1)
1543  elseif (j == 0) then
1544  this%uzspst(1, icell) = dzero
1545  this%uzthst(1, icell) = this%uzthst(this%nwavst(icell), icell)
1546  this%uzflst(1, icell) = this%uzflst(this%nwavst(icell), icell)
1547  this%nwavst(icell) = 1
1548  end if
1549  end if
1550  !
1551  ! -- calculate new unsat. storage
1552  if (thick <= dzero) then
1553  this%uzspst(1, icell) = dzero
1554  this%nwavst(icell) = 1
1555  this%uzthst(1, icell) = this%thtr(icell)
1556  this%uzflst(1, icell) = dzero
1557  end if
1558  this%watabold(icell) = this%watab(icell)
1559  end if

◆ uz_rise()

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

Definition at line 873 of file UzfCellGroup.f90.

874  ! -- dummy
875  class(UzfCellGroupType) :: this
876  integer(I4B), intent(in) :: icell
877  real(DP), intent(inout) :: totfluxtot
878  ! -- local
879  real(DP) :: fm1, fm2, d1
880  !
881  ! -- additional recharge from a rising water table
882  if (this%watab(icell) - this%watabold(icell) > dem30) then
883  d1 = this%celtop(icell) - this%watabold(icell)
884  fm1 = this%unsat_stor(icell, d1)
885  d1 = this%celtop(icell) - this%watab(icell)
886  fm2 = this%unsat_stor(icell, d1)
887  totfluxtot = totfluxtot + (fm1 - fm2)
888  end if

◆ 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 1564 of file UzfCellGroup.f90.

1565  ! -- dummy
1566  class(UzfCellGroupType) :: this
1567  integer(I4B), intent(in) :: icell
1568  real(DP), intent(in) :: delt
1569  integer(I4B), intent(in) :: ietflag
1570  integer(I4B), intent(inout) :: ierr
1571  ! -- local
1572  type(UzfCellGroupType) :: uzfktemp
1573  real(DP) :: diff
1574  real(DP) :: thetaout
1575  real(DP) :: fm
1576  real(DP) :: st
1577  real(DP) :: thtsrinv
1578  real(DP) :: epsfksthts
1579  real(DP) :: fmp
1580  real(DP) :: fktho
1581  real(DP) :: theta1
1582  real(DP) :: theta2
1583  real(DP) :: flux1
1584  real(DP) :: flux2
1585  real(DP) :: hcap
1586  real(DP) :: factor
1587  real(DP) :: tho
1588  real(DP) :: depth
1589  real(DP) :: extwc1
1590  real(DP) :: petsub
1591  integer(I4B) :: i
1592  integer(I4B) :: j
1593  integer(I4B) :: jhold
1594  integer(I4B) :: jk
1595  integer(I4B) :: kj
1596  integer(I4B) :: kk
1597  integer(I4B) :: numadd
1598  integer(I4B) :: k
1599  integer(I4B) :: nwv
1600  integer(I4B) :: itest
1601  !
1602  ! -- initialize
1603  this%etact(icell) = dzero
1604  if (this%extdpuz(icell) < dem7) return
1605  petsub = this%rootact(icell) * this%pet(icell) * &
1606  this%extdpuz(icell) / this%extdp(icell)
1607  thetaout = delt * petsub / this%extdp(icell)
1608  if (ietflag == 1) thetaout = delt * this%pet(icell) / this%extdp(icell)
1609  if (thetaout < dem10) return
1610  depth = this%uzdpst(1, icell)
1611  st = this%unsat_stor(icell, depth)
1612  if (st < dem4) return
1613  !
1614  ! -- allocate temporary wave storage.
1615  nwv = this%nwavst(icell)
1616  itest = 0
1617  call uzfktemp%init(1, nwv)
1618  !
1619  ! store original wave characteristics
1620  call uzfktemp%wave_shift(this, 1, icell, 0, 1, nwv, 1)
1621  factor = done
1622  this%etact(icell) = dzero
1623  if (this%thts(icell) - this%thtr(icell) < dem7) then
1624  thtsrinv = 1.0 / dem7
1625  else
1626  thtsrinv = done / (this%thts(icell) - this%thtr(icell))
1627  end if
1628  epsfksthts = this%eps(icell) * this%vks(icell) * thtsrinv
1629  this%etact(icell) = dzero
1630  fmp = dzero
1631  extwc1 = this%extwc(icell) - this%thtr(icell)
1632  if (extwc1 < dem6) extwc1 = dem7
1633  numadd = 0
1634  fm = st
1635  k = 0
1636  !
1637  ! -- loop for reducing aet to pet when et is head dependent
1638  do while (itest == 0)
1639  k = k + 1
1640  if (k > 1 .AND. abs(fmp - petsub) > dem5 * petsub) then
1641  factor = factor / (fm / petsub)
1642  end if
1643  !
1644  ! -- one wave shallower than extdp
1645  if (this%nwavst(icell) == 1 .AND. &
1646  this%uzdpst(1, icell) <= this%extdpuz(icell)) then
1647  if (ietflag == 2) then
1648  tho = this%uzthst(1, icell)
1649  fktho = this%uzflst(1, icell)
1650  hcap = this%caph(icell, tho)
1651  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1652  end if
1653  if ((this%uzthst(1, icell) - thetaout) > this%thtr(icell) + extwc1) then
1654  this%uzthst(1, icell) = this%uzthst(1, icell) - thetaout
1655  this%uzflst(1, icell) = &
1656  this%vks(icell) * (((this%uzthst(1, icell) - &
1657  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1658  else if (this%uzthst(1, icell) > this%thtr(icell) + extwc1) then
1659  this%uzthst(1, icell) = this%thtr(icell) + extwc1
1660  this%uzflst(1, icell) = &
1661  this%vks(icell) * (((this%uzthst(1, icell) - &
1662  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1663  end if
1664  !
1665  ! -- all waves shallower than extinction depth
1666  else if (this%nwavst(icell) > 1 .AND. &
1667  this%uzdpst(this%nwavst(icell), icell) > this%extdpuz(icell)) then
1668  if (ietflag == 2) then
1669  tho = this%uzthst(this%nwavst(icell), icell)
1670  fktho = this%uzflst(this%nwavst(icell), icell)
1671  hcap = this%caph(icell, tho)
1672  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1673  end if
1674  if (this%uzthst(this%nwavst(icell), icell) - thetaout > &
1675  this%thtr(icell) + extwc1) then
1676  this%uzthst(this%nwavst(icell) + 1, icell) = &
1677  this%uzthst(this%nwavst(icell), icell) - thetaout
1678  numadd = 1
1679  else if (this%uzthst(this%nwavst(icell), icell) > &
1680  this%thtr(icell) + extwc1) then
1681  this%uzthst(this%nwavst(icell) + 1, icell) = this%thtr(icell) + extwc1
1682  numadd = 1
1683  end if
1684  if (numadd == 1) then
1685  this%uzflst(this%nwavst(icell) + 1, icell) = &
1686  this%vks(icell) * &
1687  (((this%uzthst(this%nwavst(icell) + 1, icell) - &
1688  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1689  theta2 = this%uzthst(this%nwavst(icell) + 1, icell)
1690  flux2 = this%uzflst(this%nwavst(icell) + 1, icell)
1691  flux1 = this%uzflst(this%nwavst(icell), icell)
1692  theta1 = this%uzthst(this%nwavst(icell), icell)
1693  this%uzspst(this%nwavst(icell) + 1, icell) = &
1694  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1695  this%thtr(icell), this%eps(icell), this%vks(icell))
1696  this%uzdpst(this%nwavst(icell) + 1, icell) = this%extdpuz(icell)
1697  this%nwavst(icell) = this%nwavst(icell) + 1
1698  if (this%nwavst(icell) > this%nwav(icell)) then
1699  !
1700  ! -- too many waves error, deallocate temp arrays and return
1701  ierr = 1
1702  goto 500
1703  end if
1704  else
1705  numadd = 0
1706  end if
1707  !
1708  ! -- one wave below extinction depth
1709  else if (this%nwavst(icell) == 1) then
1710  if (ietflag == 2) then
1711  tho = this%uzthst(1, icell)
1712  fktho = this%uzflst(1, icell)
1713  hcap = this%caph(icell, tho)
1714  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1715  end if
1716  if ((this%uzthst(1, icell) - thetaout) > this%thtr(icell) + extwc1) then
1717  if (thetaout > dem30) then
1718  this%uzthst(2, icell) = this%uzthst(1, icell) - thetaout
1719  this%uzflst(2, icell) = &
1720  this%vks(icell) * (((this%uzthst(2, icell) - this%thtr(icell)) * &
1721  thtsrinv)**this%eps(icell))
1722  this%uzdpst(2, icell) = this%extdpuz(icell)
1723  theta2 = this%uzthst(2, icell)
1724  flux2 = this%uzflst(2, icell)
1725  flux1 = this%uzflst(1, icell)
1726  theta1 = this%uzthst(1, icell)
1727  this%uzspst(2, icell) = &
1728  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1729  this%thtr(icell), this%eps(icell), this%vks(icell))
1730  this%nwavst(icell) = this%nwavst(icell) + 1
1731  if (this%nwavst(icell) > this%nwav(icell)) then
1732  !
1733  ! -- too many waves error
1734  ierr = 1
1735  goto 500
1736  end if
1737  end if
1738  else if (this%uzthst(1, icell) > this%thtr(icell) + extwc1) then
1739  if (thetaout > dem30) then
1740  this%uzthst(2, icell) = this%thtr(icell) + extwc1
1741  this%uzflst(2, icell) = &
1742  this%vks(icell) * (((this%uzthst(2, icell) - &
1743  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1744  this%uzdpst(2, icell) = this%extdpuz(icell)
1745  theta2 = this%uzthst(2, icell)
1746  flux2 = this%uzflst(2, icell)
1747  flux1 = this%uzflst(1, icell)
1748  theta1 = this%uzthst(1, icell)
1749  this%uzspst(2, icell) = &
1750  leadspeed(theta1, theta2, flux1, flux2, this%thts(icell), &
1751  this%thtr(icell), this%eps(icell), this%vks(icell))
1752  this%nwavst(icell) = this%nwavst(icell) + 1
1753  if (this%nwavst(icell) > this%nwav(icell)) then
1754  !
1755  ! -- too many waves error
1756  ierr = 1
1757  goto 500
1758  end if
1759  end if
1760  end if
1761  else
1762  !
1763  ! -- extinction depth splits waves
1764  if (this%uzdpst(1, icell) - this%extdpuz(icell) > dem7) then
1765  j = 2
1766  jk = 0
1767  !
1768  ! -- locate extinction depth between waves
1769  do while (jk == 0)
1770  diff = this%uzdpst(j, icell) - this%extdpuz(icell)
1771  if (diff > dzero) then
1772  j = j + 1
1773  else
1774  jk = 1
1775  end if
1776  end do
1777  kk = j
1778  if (this%uzthst(j, icell) > this%thtr(icell) + extwc1) then
1779  !
1780  ! -- create a wave at extinction depth
1781  if (abs(diff) > dem5) then
1782  call this%wave_shift(this, icell, icell, -1, &
1783  this%nwavst(icell) + 1, j, -1)
1784  this%uzdpst(j, icell) = this%extdpuz(icell)
1785  this%nwavst(icell) = this%nwavst(icell) + 1
1786  if (this%nwavst(icell) > this%nwav(icell)) then
1787  !
1788  ! -- too many waves error
1789  ierr = 1
1790  goto 500
1791  end if
1792  end if
1793  kk = j
1794  else
1795  jhold = this%nwavst(icell)
1796  i = j + 1
1797  do while (i < this%nwavst(icell))
1798  if (this%uzthst(i, icell) > this%thtr(icell) + extwc1) then
1799  jhold = i
1800  i = this%nwavst(icell) + 1
1801  end if
1802  i = i + 1
1803  end do
1804  j = jhold
1805  kk = jhold
1806  end if
1807  else
1808  kk = 1
1809  end if
1810  !
1811  ! -- all waves above extinction depth
1812  do while (kk <= this%nwavst(icell))
1813  if (ietflag == 2) then
1814  tho = this%uzthst(kk, icell)
1815  fktho = this%uzflst(kk, icell)
1816  hcap = this%caph(icell, tho)
1817  thetaout = this%rate_et_z(icell, factor, fktho, hcap)
1818  end if
1819  if (this%uzthst(kk, icell) > this%thtr(icell) + extwc1) then
1820  if (this%uzthst(kk, icell) - thetaout > &
1821  this%thtr(icell) + extwc1) then
1822  this%uzthst(kk, icell) = this%uzthst(kk, icell) - thetaout
1823  else if (this%uzthst(kk, icell) > this%thtr(icell) + extwc1) then
1824  this%uzthst(kk, icell) = this%thtr(icell) + extwc1
1825  end if
1826  if (kk == 1) then
1827  this%uzflst(kk, icell) = &
1828  this%vks(icell) * &
1829  (((this%uzthst(kk, icell) - &
1830  this%thtr(icell)) * thtsrinv)**this%eps(icell))
1831  end if
1832  if (kk > 1) then
1833  flux1 = &
1834  this%vks(icell) * ((this%uzthst(kk - 1, icell) - &
1835  this%thtr(icell)) * thtsrinv)**this%eps(icell)
1836  flux2 = &
1837  this%vks(icell) * ((this%uzthst(kk, icell) - &
1838  this%thtr(icell)) * thtsrinv)**this%eps(icell)
1839  this%uzflst(kk, icell) = flux2
1840  theta2 = this%uzthst(kk, icell)
1841  theta1 = this%uzthst(kk - 1, icell)
1842  this%uzspst(kk, icell) = leadspeed(theta1, theta2, flux1, flux2, &
1843  this%thts(icell), &
1844  this%thtr(icell), &
1845  this%eps(icell), this%vks(icell))
1846  end if
1847  end if
1848  kk = kk + 1
1849  end do
1850  end if
1851  !
1852  ! -- calculate aet
1853  kj = 1
1854  do while (kj <= this%nwavst(icell) - 1)
1855  if (abs(this%uzthst(kj, icell) - this%uzthst(kj + 1, icell)) < dem6) then
1856  call this%wave_shift(this, icell, icell, 1, kj + 1, &
1857  this%nwavst(icell) - 1, 1)
1858  kj = kj - 1
1859  this%nwavst(icell) = this%nwavst(icell) - 1
1860  end if
1861  kj = kj + 1
1862  end do
1863  depth = this%uzdpst(1, icell)
1864  fm = this%unsat_stor(icell, depth)
1865  this%etact(icell) = st - fm
1866  fm = this%etact(icell) / delt
1867  if (this%etact(icell) < dzero) then
1868  call this%wave_shift(uzfktemp, icell, 1, 0, 1, nwv, 1)
1869  this%nwavst(icell) = nwv
1870  this%etact(icell) = dzero
1871  elseif (petsub - fm < -dem15 .AND. ietflag == 2) then
1872  !
1873  ! -- aet greater than pet, reset and try again
1874  call this%wave_shift(uzfktemp, icell, 1, 0, 1, nwv, 1)
1875  this%nwavst(icell) = nwv
1876  this%etact(icell) = dzero
1877  else
1878  itest = 1
1879  end if
1880  !
1881  ! -- end aet-pet loop for head dependent et
1882  fmp = fm
1883  if (k > 100) then
1884  itest = 1
1885  elseif (ietflag < 2) then
1886  fmp = petsub
1887  itest = 1
1888  end if
1889  end do
1890 500 continue
1891  !
1892  ! -- deallocate temporary worker
1893  call uzfktemp%dealloc()
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 1005 of file UzfCellGroup.f90.

1006  ! -- dummy
1007  class(UzfCellGroupType) :: this
1008  real(DP), intent(inout) :: thickold
1009  real(DP), intent(inout) :: thick
1010  real(DP), intent(in) :: delt
1011  integer(I4B), intent(in) :: ietflag
1012  integer(I4B), intent(in) :: icell
1013  integer(I4B), intent(inout) :: ierr
1014  ! -- local
1015  real(DP) :: ffcheck, time, feps1, feps2
1016  real(DP) :: thetadif, thetab, fluxb, oldsflx
1017  integer(I4B) :: itrailflg, itester
1018  !
1019  time = dzero
1020  this%totflux(icell) = dzero
1021  itrailflg = 0
1022  oldsflx = this%uzflst(this%nwavst(icell), icell)
1023  call factors(feps1, feps2)
1024  !
1025  ! -- check for falling or rising water table
1026  if ((thick - thickold) > feps1) then
1027  thetadif = abs(this%uzthst(1, icell) - this%thtr(icell))
1028  if (thetadif > dem6) then
1029  call this%wave_shift(this, icell, icell, -1, &
1030  this%nwavst(icell) + 1, 2, -1)
1031  if (this%uzdpst(2, icell) < dem30) &
1032  this%uzdpst(2, icell) = (this%ntrail(icell) + dtwo) * dem6
1033  if (this%uzthst(2, icell) > this%thtr(icell)) then
1034  this%uzspst(2, icell) = this%uzflst(2, icell) / &
1035  (this%uzthst(2, icell) - this%thtr(icell))
1036  else
1037  this%uzspst(2, icell) = dzero
1038  end if
1039  this%uzthst(1, icell) = this%thtr(icell)
1040  this%uzflst(1, icell) = dzero
1041  this%uzspst(1, icell) = dzero
1042  this%uzdpst(1, icell) = thick
1043  this%nwavst(icell) = this%nwavst(icell) + 1
1044  if (this%nwavst(icell) >= this%nwav(icell)) then
1045  ! -- too many waves error
1046  ierr = 1
1047  return
1048  end if
1049  else
1050  this%uzdpst(1, icell) = thick
1051  end if
1052  end if
1053  thetab = this%uzthst(1, icell)
1054  fluxb = this%uzflst(1, icell)
1055  this%totflux(icell) = dzero
1056  itester = 0
1057  ffcheck = (this%surflux(icell) - this%uzflst(this%nwavst(icell), icell))
1058  !
1059  ! -- increase new waves in infiltration changes
1060  if (ffcheck > feps2 .OR. ffcheck < -feps2) then
1061  this%nwavst(icell) = this%nwavst(icell) + 1
1062  if (this%nwavst(icell) >= this%nwav(icell)) then
1063  !
1064  ! -- too many waves error
1065  ierr = 1
1066  return
1067  end if
1068  else if (this%nwavst(icell) == 1) then
1069  itester = 1
1070  end if
1071  if (this%nwavst(icell) > 1) then
1072  if (ffcheck < -feps2) then
1073  call this%trailwav(icell, ierr)
1074  if (ierr > 0) return
1075  itrailflg = 1
1076  end if
1077  call this%leadwav(time, itester, itrailflg, thetab, fluxb, ffcheck, &
1078  feps2, delt, icell)
1079  end if
1080  if (itester == 1) then
1081  this%totflux(icell) = this%totflux(icell) + &
1082  (delt - time) * this%uzflst(1, icell)
1083  time = dzero
1084  itester = 0
1085  end if
1086  !
1087  ! -- simulate et
1088  if (ietflag > 0) call this%uzet(icell, delt, ietflag, ierr)
1089  if (ierr > 0) 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 980 of file UzfCellGroup.f90.

981  ! -- dummy
982  class(UzfCellGroupType) :: this
983  type(UzfCellGroupType) :: this2
984  integer(I4B), intent(in) :: icell
985  integer(I4B), intent(in) :: icell2
986  integer(I4B), intent(in) :: shft
987  integer(I4B), intent(in) :: strt
988  integer(I4B), intent(in) :: stp
989  integer(I4B), intent(in) :: cntr
990  ! -- local
991  integer(I4B) :: j
992  !
993  ! -- copy waves from one uzf cell group to another
994  do j = strt, stp, cntr
995  this%uzthst(j, icell) = this2%uzthst(j + shft, icell2)
996  this%uzdpst(j, icell) = this2%uzdpst(j + shft, icell2)
997  this%uzflst(j, icell) = this2%uzflst(j + shft, icell2)
998  this%uzspst(j, icell) = this2%uzspst(j + shft, icell2)
999  end do
1000  this%nwavst(icell) = this2%nwavst(icell2)