30 integer(I4B),
pointer :: iadvwt => null()
31 real(dp),
pointer :: ats_percel => null()
32 integer(I4B),
dimension(:),
pointer,
contiguous :: ibound => null()
34 real(dp),
pointer :: eqnsclfac => null()
58 subroutine adv_cr(advobj, name_model, input_mempath, inunit, iout, fmi, &
62 character(len=*),
intent(in) :: name_model
63 character(len=*),
intent(in) :: input_mempath
64 integer(I4B),
intent(in) :: inunit
65 integer(I4B),
intent(in) :: iout
67 real(dp),
intent(in),
pointer :: eqnsclfac
73 call advobj%set_names(1, name_model,
'ADV',
'ADV', input_mempath)
76 call advobj%allocate_scalars()
79 advobj%inunit = inunit
82 advobj%eqnsclfac => eqnsclfac
94 character(len=*),
parameter :: fmtadv = &
95 "(1x,/1x,'ADV -- ADVECTION PACKAGE, VERSION 1, 8/25/2017', &
96 &' INPUT READ FROM MEMPATH: ', A, //)"
99 if (.not.
present(adv_options))
then
102 write (this%iout, fmtadv) this%input_mempath
105 call this%source_options()
109 this%iadvwt = adv_options%iAdvScheme
123 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(in) :: ibound
125 integer(I4B) :: iadvwt_value
129 this%ibound => ibound
132 iadvwt_value = this%iadvwt
133 select case (iadvwt_value)
135 this%face_interpolation = &
138 this%face_interpolation = &
141 this%face_interpolation = &
146 this%face_interpolation = &
149 call store_error(
"Unknown advection scheme", terminate=.true.)
157 subroutine adv_dt(this, dtmax, msg, thetam)
160 real(DP),
intent(out) :: dtmax
161 character(len=*),
intent(inout) :: msg
162 real(DP),
dimension(:),
intent(in) :: thetam
167 integer(I4B) :: nrmax
168 character(len=LINELENGTH) :: cellstr
171 real(DP) :: flowsumpos
172 real(DP) :: flowsumneg
174 real(DP) :: cell_volume
181 if (this%ats_percel ==
dnodata)
then
187 do n = 1, this%dis%nodes
188 if (this%ibound(n) == 0) cycle
191 do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
192 if (this%dis%con%mask(ipos) == 0) cycle
193 m = this%dis%con%ja(ipos)
194 if (this%ibound(m) == 0) cycle
195 flownm = this%fmi%gwfflowja(ipos)
196 if (flownm <
dzero)
then
197 flowsumneg = flowsumneg - flownm
199 flowsumpos = flowsumpos + flownm
202 flowmax = max(flowsumneg, flowsumpos)
203 if (flowmax <
dprec) cycle
204 cell_volume = this%dis%get_cell_volume(n, this%dis%top(n))
205 dt = cell_volume * this%fmi%gwfsat(n) * thetam(n) / flowmax
206 dt = dt * this%ats_percel
213 call this%dis%noder_to_string(nrmax, cellstr)
214 write (msg, *) adjustl(trim(this%memoryPath))//
'-'//trim(cellstr)
222 subroutine adv_fc(this, nodes, matrix_sln, idxglo, cnew, rhs)
226 integer(I4B),
intent(in) :: nodes
228 integer(I4B),
intent(in),
dimension(:) :: idxglo
229 real(DP),
intent(in),
dimension(:),
target :: cnew
230 real(DP),
dimension(:),
intent(inout) :: rhs
232 integer(I4B) :: n, m, idiag, ipos
237 call this%face_interpolation%set_field(cnew)
239 if (this%ibound(n) == 0) cycle
240 idiag = this%dis%con%ia(n)
241 do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
242 m = this%dis%con%ja(ipos)
243 if (this%ibound(m) == 0) cycle
244 if (this%dis%con%mask(ipos) == 0) cycle
246 qnm = this%fmi%gwfflowja(ipos) * this%eqnsclfac
247 coefficients = this%face_interpolation%compute(n, m, ipos)
249 call matrix_sln%add_value_pos(idxglo(idiag), qnm * coefficients%c_n)
250 call matrix_sln%add_value_pos(idxglo(ipos), qnm * coefficients%c_m)
251 rhs(n) = rhs(n) + qnm * coefficients%rhs
262 real(DP),
intent(in),
dimension(:),
target :: cnew
263 real(DP),
intent(inout),
dimension(:) :: flowja
265 integer(I4B) :: nodes
266 integer(I4B) :: n, m, ipos
272 nodes = this%dis%nodes
274 call this%face_interpolation%set_field(cnew)
276 if (this%ibound(n) == 0) cycle
277 do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
278 m = this%dis%con%ja(ipos)
279 if (this%ibound(m) == 0) cycle
280 qnm = this%fmi%gwfflowja(ipos) * this%eqnsclfac
282 coefficients = this%face_interpolation%compute(n, m, ipos)
283 flowja(ipos) = flowja(ipos) &
284 + qnm * coefficients%c_n * cnew(n) &
285 + qnm * coefficients%c_m * cnew(m) &
286 - qnm * coefficients%rhs
301 if (this%inunit > 0)
then
305 this%ibound => null()
312 call this%NumericalPackageType%da()
326 call this%NumericalPackageType%allocate_scalars()
329 call mem_allocate(this%iadvwt,
'IADVWT', this%memoryPath)
330 call mem_allocate(this%ats_percel,
'ATS_PERCEL', this%memoryPath)
352 character(len=LENVARNAME),
dimension(4) :: scheme = &
353 &[character(len=LENVARNAME) ::
'UPSTREAM',
'CENTRAL',
'TVD',
'UTVD']
354 logical(lgp) :: found_scheme, found_atspercel
356 character(len=*),
parameter :: fmtiadvwt = &
357 &
"(4x,'ADVECTION WEIGHTING SCHEME HAS BEEN SET TO: ', a)"
360 call mem_set_value(this%iadvwt,
'SCHEME', this%input_mempath, &
361 scheme, found_scheme)
362 call mem_set_value(this%ats_percel,
'ATS_PERCEL', this%input_mempath, &
365 if (found_scheme)
then
367 if (this%iadvwt == 0)
then
368 write (
errmsg,
'(a, a)') &
369 'Unknown scheme, must be "UPSTREAM", "CENTRAL", "TVD" or "UTVD"'
374 this%iadvwt = this%iadvwt - 1
379 call dev_feature(
'UTVD is still under development, install the &
380 &nightly build or compile from source with IDEVELOPMODE = 1.')
383 if (found_atspercel)
then
384 if (this%ats_percel == dzero) this%ats_percel = dnodata
388 write (this%iout,
'(1x,a)')
'PROCESSING ADVECTION OPTIONS'
389 if (found_scheme)
then
390 write (this%iout, fmtiadvwt) trim(scheme(this%iadvwt + 1))
392 if (found_atspercel)
then
393 write (this%iout,
'(4x,a,1pg15.6)') &
394 'User-specified fractional cell distance for adaptive time &
395 &steps: ', this%ats_percel
397 write (this%iout,
'(1x,a)')
'END OF ADVECTION OPTIONS'
integer(i4b), parameter adv_scheme_tvd
integer(i4b), parameter adv_scheme_upstream
integer(i4b), parameter adv_scheme_utvd
integer(i4b), parameter adv_scheme_central
This module contains simulation constants.
integer(i4b), parameter linelength
maximum length of a standard line
real(dp), parameter dnodata
real no data constant
integer(i4b), parameter lenvarname
maximum length of a variable name
real(dp), parameter dzero
real constant zero
real(dp), parameter dprec
real constant machine precision
real(dp), parameter done
real constant 1
Disable development features in release mode.
subroutine, public dev_feature(errmsg, iunit)
Terminate if in release mode (guard development features)
This module defines variable data types.
This module contains the base numerical package type.
This module contains simulation methods.
subroutine, public store_error(msg, terminate)
Store an error message.
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
subroutine, public adv_cr(advobj, name_model, input_mempath, inunit, iout, fmi, eqnsclfac)
@ brief Create a new ADV object
subroutine source_options(this)
Source input options.
subroutine adv_df(this, adv_options)
Define ADV object.
subroutine adv_fc(this, nodes, matrix_sln, idxglo, cnew, rhs)
Fill coefficient method for ADV package.
subroutine adv_cq(this, cnew, flowja)
Calculate advection contribution to flowja.
subroutine adv_da(this)
Deallocate memory.
subroutine adv_ar(this, dis, ibound)
Allocate and read method for package.
subroutine adv_dt(this, dtmax, msg, thetam)
Calculate maximum time step length.
subroutine allocate_scalars(this)
Allocate scalars specific to the streamflow energy transport (SFE) package.
Decorator that adds caching to any gradient computation implementation.
Abstract interface for cell-based gradient computation.
Weighted least-squares gradient method for structured and unstructured grids.
Total Variation Diminishing (TVD) interpolation scheme.