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

This module contains the WEL package methods. More...

Data Types

type  weltype
 

Functions/Subroutines

subroutine, public wel_create (packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath)
 @ brief Create a new package object More...
 
subroutine wel_da (this)
 @ brief Deallocate package memory More...
 
subroutine wel_allocate_scalars (this)
 @ brief Allocate scalars More...
 
subroutine wel_allocate_arrays (this, nodelist, auxvar)
 @ brief Allocate arrays More...
 
subroutine wel_options (this)
 @ brief Source additional options for package More...
 
subroutine log_wel_options (this, found)
 @ brief Log WEL specific package options More...
 
subroutine wel_rp (this)
 @ brief WEL read and prepare More...
 
subroutine wel_cf (this)
 @ brief Formulate the package hcof and rhs terms. More...
 
subroutine wel_fc (this, rhs, ia, idxglo, matrix_sln)
 @ brief Copy hcof and rhs terms into solution. More...
 
subroutine wel_fn (this, rhs, ia, idxglo, matrix_sln)
 @ brief Add Newton-Raphson terms for package into solution. More...
 
subroutine wel_afr_csv_init (this, fname)
 Initialize the auto flow reduce csv output file. More...
 
subroutine wel_afr_csv_write (this)
 Write out auto flow reductions only when & where they occur. More...
 
subroutine define_listlabel (this)
 @ brief Define the list label for the package More...
 
logical function wel_obs_supported (this)
 Determine if observations are supported. More...
 
subroutine wel_df_obs (this)
 Define the observation types available in the package. More...
 
subroutine wel_bd_obs (this)
 Save observations for the package. More...
 
real(dp) function q_mult (this, row)
 
real(dp) function wel_bound_value (this, col, row)
 @ brief Return a bound value More...
 

Variables

character(len=lenftype) ftype = 'WEL'
 package ftype More...
 
character(len=16) text = ' WEL'
 package flow text string More...
 

Detailed Description

This module contains the overridden methods for the standard WEL package. Several methods need to be overridden because of the AUTO_FLOW_REDUCE option. Overridden methods include:

  • bnd_cf (AUTO_FLOW_REDUCE)
  • bnd_fc (AUTO_FLOW_REDUCE)
  • bnd_fn (AUTO_FLOW_REDUCE Newton-Raphson terms)
  • bnd_ot_package_flows (write AUTO_FLOW_REDUCE terms to csv file)
  • bnd_da (deallocate AUTO_FLOW_REDUCE variables)
  • bnd_bd_obs (wel-reduction observation added)

Function/Subroutine Documentation

◆ define_listlabel()

subroutine welmodule::define_listlabel ( class(weltype), intent(inout)  this)

Method defined the list label for the WEL package. The list label is the heading that is written to iout when PRINT_INPUT option is used.

Parameters
[in,out]thisWelType object

Definition at line 480 of file gwf-wel.f90.

481  ! -- dummy variables
482  class(WelType), intent(inout) :: this !< WelType object
483  !
484  ! -- create the header list label
485  this%listlabel = trim(this%filtyp)//' NO.'
486  if (this%dis%ndim == 3) then
487  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'LAYER'
488  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'ROW'
489  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'COL'
490  elseif (this%dis%ndim == 2) then
491  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'LAYER'
492  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'CELL2D'
493  else
494  write (this%listlabel, '(a, a7)') trim(this%listlabel), 'NODE'
495  end if
496  write (this%listlabel, '(a, a16)') trim(this%listlabel), 'STRESS RATE'
497  if (this%inamedbound == 1) then
498  write (this%listlabel, '(a, a16)') trim(this%listlabel), 'BOUNDARY NAME'
499  end if

◆ log_wel_options()

subroutine welmodule::log_wel_options ( class(weltype), intent(inout)  this,
type(gwfwelparamfoundtype), intent(in)  found 
)
Parameters
[in,out]thisBndExtType object

Definition at line 234 of file gwf-wel.f90.

235  ! -- modules
237  ! -- dummy variables
238  class(WelType), intent(inout) :: this !< BndExtType object
239  type(GwfWelParamFoundType), intent(in) :: found
240  ! -- local variables
241  ! -- format
242  character(len=*), parameter :: fmtflowred = &
243  &"(4x, 'AUTOMATIC FLOW REDUCTION OF WELLS IMPLEMENTED.')"
244  character(len=*), parameter :: fmtflowredv = &
245  &"(4x, 'AUTOMATIC FLOW REDUCTION FRACTION (',g15.7,').')"
246  !
247  ! -- log found options
248  write (this%iout, '(/1x,a)') 'PROCESSING '//trim(adjustl(this%text)) &
249  //' OPTIONS'
250  !
251  if (found%flowred) then
252  if (this%iflowred > 0) &
253  write (this%iout, fmtflowred)
254  write (this%iout, fmtflowredv) this%flowred
255  end if
256  !
257  if (found%afrcsvfile) then
258  ! -- currently no-op
259  end if
260  !
261  if (found%mover) then
262  write (this%iout, '(4x,A)') 'MOVER OPTION ENABLED'
263  end if
264  !
265  ! -- close logging block
266  write (this%iout, '(1x,a)') &
267  'END OF '//trim(adjustl(this%text))//' OPTIONS'

◆ q_mult()

real(dp) function welmodule::q_mult ( class(weltype), intent(inout)  this,
integer(i4b), intent(in)  row 
)
private
Parameters
[in,out]thisBndExtType object

Definition at line 602 of file gwf-wel.f90.

603  ! -- modules
604  use constantsmodule, only: dzero
605  ! -- dummy variables
606  class(WelType), intent(inout) :: this !< BndExtType object
607  integer(I4B), intent(in) :: row
608  ! -- result
609  real(DP) :: q
610  !
611  if (this%iauxmultcol > 0) then
612  q = this%q(row) * this%auxvar(this%iauxmultcol, row)
613  else
614  q = this%q(row)
615  end if
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65

◆ wel_afr_csv_init()

subroutine welmodule::wel_afr_csv_init ( class(weltype), intent(inout)  this,
character(len=*), intent(in)  fname 
)
private
Parameters
[in,out]thisWelType object

Definition at line 427 of file gwf-wel.f90.

428  ! -- dummy variables
429  class(WelType), intent(inout) :: this !< WelType object
430  character(len=*), intent(in) :: fname
431  ! -- format
432  character(len=*), parameter :: fmtafrcsv = &
433  "(4x, 'AUTO FLOW REDUCE INFORMATION WILL BE SAVED TO FILE: ', a, /4x, &
434  &'OPENED ON UNIT: ', I0)"
435 
436  this%ioutafrcsv = getunit()
437  call openfile(this%ioutafrcsv, this%iout, fname, 'CSV', &
438  filstat_opt='REPLACE')
439  write (this%iout, fmtafrcsv) trim(adjustl(fname)), &
440  this%ioutafrcsv
441  write (this%ioutafrcsv, '(a)') &
442  'time,period,step,boundnumber,cellnumber,rate-requested,&
443  &rate-actual,wel-reduction'
Here is the call graph for this function:

◆ wel_afr_csv_write()

subroutine welmodule::wel_afr_csv_write ( class(weltype), intent(inout)  this)
private
Parameters
[in,out]thisWelType object

Definition at line 447 of file gwf-wel.f90.

448  ! -- modules
449  use tdismodule, only: totim, kstp, kper
450  ! -- dummy variables
451  class(WelType), intent(inout) :: this !< WelType object
452  ! -- local
453  integer(I4B) :: i
454  integer(I4B) :: nodereduced
455  integer(I4B) :: nodeuser
456  real(DP) :: v
457  ! -- format
458  do i = 1, this%nbound
459  nodereduced = this%nodelist(i)
460  !
461  ! -- test if node is constant or inactive
462  if (this%ibound(nodereduced) <= 0) then
463  cycle
464  end if
465  v = this%q_mult(i) + this%rhs(i)
466  if (v < dzero) then
467  nodeuser = this%dis%get_nodeuser(nodereduced)
468  write (this%ioutafrcsv, '(*(G0,:,","))') &
469  totim, kper, kstp, i, nodeuser, this%q_mult(i), this%simvals(i), v
470  end if
471  end do
real(dp), pointer, public totim
time relative to start of simulation
Definition: tdis.f90:32
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23

◆ wel_allocate_arrays()

subroutine welmodule::wel_allocate_arrays ( class(weltype this,
integer(i4b), dimension(:), optional, pointer, contiguous  nodelist,
real(dp), dimension(:, :), optional, pointer, contiguous  auxvar 
)

Allocate and initialize arrays for the WEL package

Definition at line 160 of file gwf-wel.f90.

161  ! -- modules
163  ! -- dummy
164  class(WelType) :: this
165  integer(I4B), dimension(:), pointer, contiguous, optional :: nodelist
166  real(DP), dimension(:, :), pointer, contiguous, optional :: auxvar
167  ! -- local
168  !
169  ! -- call BndExtType allocate scalars
170  call this%BndExtType%allocate_arrays(nodelist, auxvar)
171  !
172  ! -- set constant head array input context pointer
173  call mem_setptr(this%q, 'Q', this%input_mempath)
174  !
175  ! -- checkin constant head array input context pointer
176  call mem_checkin(this%q, 'Q', this%memoryPath, &
177  'Q', this%input_mempath)

◆ wel_allocate_scalars()

subroutine welmodule::wel_allocate_scalars ( class(weltype this)

Allocate and initialize scalars for the WEL package. The base model allocate scalars method is also called.

Parameters
thisWelType object

Definition at line 135 of file gwf-wel.f90.

136  ! -- modules
138  ! -- dummy variables
139  class(WelType) :: this !< WelType object
140  !
141  ! -- call base type allocate scalars
142  call this%BndExtType%allocate_scalars()
143  !
144  ! -- allocate the object and assign values to object variables
145  call mem_allocate(this%iflowred, 'IFLOWRED', this%memoryPath)
146  call mem_allocate(this%flowred, 'FLOWRED', this%memoryPath)
147  call mem_allocate(this%ioutafrcsv, 'IOUTAFRCSV', this%memoryPath)
148  !
149  ! -- Set values
150  this%iflowred = 0
151  this%ioutafrcsv = 0
152  this%flowred = dzero

◆ wel_bd_obs()

subroutine welmodule::wel_bd_obs ( class(weltype this)
private

Method to save simulated values for the WEL package.

Parameters
thisWelType object

Definition at line 551 of file gwf-wel.f90.

552  ! -- dummy variables
553  class(WelType) :: this !< WelType object
554  ! -- local variables
555  integer(I4B) :: i
556  integer(I4B) :: n
557  integer(I4B) :: jj
558  real(DP) :: v
559  type(ObserveType), pointer :: obsrv => null()
560  !
561  ! -- clear the observations
562  call this%obs%obs_bd_clear()
563  !
564  ! -- Save simulated values for all of package's observations.
565  do i = 1, this%obs%npakobs
566  obsrv => this%obs%pakobs(i)%obsrv
567  if (obsrv%BndFound) then
568  do n = 1, obsrv%indxbnds_count
569  v = dnodata
570  jj = obsrv%indxbnds(n)
571  select case (obsrv%ObsTypeId)
572  case ('TO-MVR')
573  if (this%imover == 1) then
574  v = this%pakmvrobj%get_qtomvr(jj)
575  if (v > dzero) then
576  v = -v
577  end if
578  end if
579  case ('WEL')
580  v = this%simvals(jj)
581  case ('WEL-REDUCTION')
582  if (this%iflowred > 0) then
583  v = this%q_mult(jj) + this%rhs(jj)
584  end if
585  case default
586  errmsg = 'Unrecognized observation type: '//trim(obsrv%ObsTypeId)
587  call store_error(errmsg)
588  end select
589  call this%obs%SaveOneSimval(obsrv, v)
590  end do
591  else
592  call this%obs%SaveOneSimval(obsrv, dnodata)
593  end if
594  end do
595  !
596  ! -- Write the auto flow reduce csv file entries for this step
597  if (this%ioutafrcsv > 0) then
598  call this%wel_afr_csv_write()
599  end if
Here is the call graph for this function:

◆ wel_bound_value()

real(dp) function welmodule::wel_bound_value ( class(weltype), intent(inout)  this,
integer(i4b), intent(in)  col,
integer(i4b), intent(in)  row 
)

Return a bound value associated with an ncolbnd index and row.

Parameters
[in,out]thisBndExtType object

Definition at line 624 of file gwf-wel.f90.

625  ! -- modules
626  use constantsmodule, only: dzero
627  ! -- dummy variables
628  class(WelType), intent(inout) :: this !< BndExtType object
629  integer(I4B), intent(in) :: col
630  integer(I4B), intent(in) :: row
631  ! -- result
632  real(DP) :: bndval
633  !
634  select case (col)
635  case (1)
636  bndval = this%q_mult(row)
637  case default
638  errmsg = 'Programming error. WEL bound value requested column '&
639  &'outside range of ncolbnd (1).'
640  call store_error(errmsg)
641  call store_error_filename(this%input_fname)
642  end select
Here is the call graph for this function:

◆ wel_cf()

subroutine welmodule::wel_cf ( class(weltype this)

Formulate the hcof and rhs terms for the WEL package that will be added to the coefficient matrix and right-hand side vector.

Parameters
thisWelType object

Definition at line 296 of file gwf-wel.f90.

297  ! -- dummy variables
298  class(WelType) :: this !< WelType object
299  ! -- local variables
300  integer(I4B) :: i, node, ict
301  real(DP) :: qmult
302  real(DP) :: q
303  real(DP) :: tp
304  real(DP) :: bt
305  real(DP) :: thick
306  !
307  ! -- Return if no wells
308  if (this%nbound == 0) return
309  !
310  ! -- Calculate hcof and rhs for each well entry
311  do i = 1, this%nbound
312  node = this%nodelist(i)
313  this%hcof(i) = dzero
314  if (this%ibound(node) <= 0) then
315  this%rhs(i) = dzero
316  cycle
317  end if
318  q = this%q_mult(i)
319  if (this%iflowred /= 0 .and. q < dzero) then
320  ict = this%icelltype(node)
321  if (ict /= 0) then
322  tp = this%dis%top(node)
323  bt = this%dis%bot(node)
324  thick = tp - bt
325  tp = bt + this%flowred * thick
326  qmult = sqsaturation(tp, bt, this%xnew(node))
327  q = q * qmult
328  end if
329  end if
330  this%rhs(i) = -q
331  end do
Here is the call graph for this function:

◆ wel_create()

subroutine, public welmodule::wel_create ( class(bndtype), pointer  packobj,
integer(i4b), intent(in)  id,
integer(i4b), intent(in)  ibcnum,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout,
character(len=*), intent(in)  namemodel,
character(len=*), intent(in)  pakname,
character(len=*), intent(in)  mempath 
)

Create a new WEL Package object

Parameters
packobjpointer to default package type
[in]idpackage id
[in]ibcnumboundary condition number
[in]inunitunit number of WEL package input file
[in]ioutunit number of model listing file
[in]namemodelmodel name
[in]paknamepackage name
[in]mempathinput mempath

Definition at line 73 of file gwf-wel.f90.

75  ! -- dummy variables
76  class(BndType), pointer :: packobj !< pointer to default package type
77  integer(I4B), intent(in) :: id !< package id
78  integer(I4B), intent(in) :: ibcnum !< boundary condition number
79  integer(I4B), intent(in) :: inunit !< unit number of WEL package input file
80  integer(I4B), intent(in) :: iout !< unit number of model listing file
81  character(len=*), intent(in) :: namemodel !< model name
82  character(len=*), intent(in) :: pakname !< package name
83  character(len=*), intent(in) :: mempath !< input mempath
84  ! -- local variables
85  type(WelType), pointer :: welobj
86  !
87  ! -- allocate the object and assign values to object variables
88  allocate (welobj)
89  packobj => welobj
90  !
91  ! -- create name and memory path
92  call packobj%set_names(ibcnum, namemodel, pakname, ftype, mempath)
93  packobj%text = text
94  !
95  ! -- allocate scalars
96  call welobj%allocate_scalars()
97  !
98  ! -- initialize package
99  call packobj%pack_initialize()
100 
101  packobj%inunit = inunit
102  packobj%iout = iout
103  packobj%id = id
104  packobj%ibcnum = ibcnum
105  packobj%ictMemPath = create_mem_path(namemodel, 'NPF')
Here is the call graph for this function:
Here is the caller graph for this function:

◆ wel_da()

subroutine welmodule::wel_da ( class(weltype this)
private

Deallocate WEL package scalars and arrays.

Parameters
thisWelType object

Definition at line 113 of file gwf-wel.f90.

114  ! -- modules
116  ! -- dummy variables
117  class(WelType) :: this !< WelType object
118  !
119  ! -- Deallocate parent package
120  call this%BndExtType%bnd_da()
121  !
122  ! -- scalars
123  call mem_deallocate(this%iflowred)
124  call mem_deallocate(this%flowred)
125  call mem_deallocate(this%ioutafrcsv)
126  call mem_deallocate(this%q, 'Q', this%memoryPath)

◆ wel_df_obs()

subroutine welmodule::wel_df_obs ( class(weltype this)
private

Method to define the observation types available in the WEL package.

Parameters
thisWelType object

Definition at line 525 of file gwf-wel.f90.

526  ! -- dummy variables
527  class(WelType) :: this !< WelType object
528  ! -- local variables
529  integer(I4B) :: indx
530  !
531  ! -- initialize observations
532  call this%obs%StoreObsType('wel', .true., indx)
533  this%obs%obsData(indx)%ProcessIdPtr => defaultobsidprocessor
534  !
535  ! -- Store obs type and assign procedure pointer
536  ! for to-mvr observation type.
537  call this%obs%StoreObsType('to-mvr', .true., indx)
538  this%obs%obsData(indx)%ProcessIdPtr => defaultobsidprocessor
539  !
540  ! -- Store obs type and assign procedure pointer
541  ! for wel-reduction observation type.
542  call this%obs%StoreObsType('wel-reduction', .true., indx)
543  this%obs%obsData(indx)%ProcessIdPtr => defaultobsidprocessor
Here is the call graph for this function:

◆ wel_fc()

subroutine welmodule::wel_fc ( class(weltype this,
real(dp), dimension(:), intent(inout)  rhs,
integer(i4b), dimension(:), intent(in)  ia,
integer(i4b), dimension(:), intent(in)  idxglo,
class(matrixbasetype), pointer  matrix_sln 
)
private

Add the hcof and rhs terms for the WEL package to the coefficient matrix and right-hand side vector.

Parameters
thisWelType object
[in,out]rhsright-hand side vector for model
[in]iasolution CRS row pointers
[in]idxglomapping vector for model (local) to solution (global)
matrix_slnsolution coefficient matrix

Definition at line 340 of file gwf-wel.f90.

341  ! -- dummy variables
342  class(WelType) :: this !< WelType object
343  real(DP), dimension(:), intent(inout) :: rhs !< right-hand side vector for model
344  integer(I4B), dimension(:), intent(in) :: ia !< solution CRS row pointers
345  integer(I4B), dimension(:), intent(in) :: idxglo !< mapping vector for model (local) to solution (global)
346  class(MatrixBaseType), pointer :: matrix_sln !< solution coefficient matrix
347  ! -- local variables
348  integer(I4B) :: i
349  integer(I4B) :: n
350  integer(I4B) :: ipos
351  !
352  ! -- pakmvrobj fc
353  if (this%imover == 1) then
354  call this%pakmvrobj%fc()
355  end if
356  !
357  ! -- Copy package rhs and hcof into solution rhs and amat
358  do i = 1, this%nbound
359  n = this%nodelist(i)
360  rhs(n) = rhs(n) + this%rhs(i)
361  ipos = ia(n)
362  call matrix_sln%add_value_pos(idxglo(ipos), this%hcof(i))
363  !
364  ! -- If mover is active and this well is discharging,
365  ! store available water (as positive value).
366  if (this%imover == 1 .and. this%rhs(i) > dzero) then
367  call this%pakmvrobj%accumulate_qformvr(i, this%rhs(i))
368  end if
369  end do

◆ wel_fn()

subroutine welmodule::wel_fn ( class(weltype this,
real(dp), dimension(:), intent(inout)  rhs,
integer(i4b), dimension(:), intent(in)  ia,
integer(i4b), dimension(:), intent(in)  idxglo,
class(matrixbasetype), pointer  matrix_sln 
)
private

Calculate and add the Newton-Raphson terms for the WEL package to the coefficient matrix and right-hand side vector.

Parameters
thisWelType object
[in,out]rhsright-hand side vector for model
[in]iasolution CRS row pointers
[in]idxglomapping vector for model (local) to solution (global)
matrix_slnsolution coefficient matrix

Definition at line 378 of file gwf-wel.f90.

379  ! -- dummy variables
380  class(WelType) :: this !< WelType object
381  real(DP), dimension(:), intent(inout) :: rhs !< right-hand side vector for model
382  integer(I4B), dimension(:), intent(in) :: ia !< solution CRS row pointers
383  integer(I4B), dimension(:), intent(in) :: idxglo !< mapping vector for model (local) to solution (global)
384  class(MatrixBaseType), pointer :: matrix_sln !< solution coefficient matrix
385  ! -- local variables
386  integer(I4B) :: i
387  integer(I4B) :: node
388  integer(I4B) :: ipos
389  integer(I4B) :: ict
390  real(DP) :: drterm
391  real(DP) :: q
392  real(DP) :: tp
393  real(DP) :: bt
394  real(DP) :: thick
395  !
396  ! -- Copy package rhs and hcof into solution rhs and amat
397  do i = 1, this%nbound
398  node = this%nodelist(i)
399  !
400  ! -- test if node is constant or inactive
401  if (this%ibound(node) <= 0) then
402  cycle
403  end if
404  !
405  ! -- well rate is possibly head dependent
406  ict = this%icelltype(node)
407  if (this%iflowred /= 0 .and. ict /= 0) then
408  ipos = ia(node)
409  q = -this%rhs(i)
410  if (q < dzero) then
411  ! -- calculate derivative for well
412  tp = this%dis%top(node)
413  bt = this%dis%bot(node)
414  thick = tp - bt
415  tp = bt + this%flowred * thick
416  drterm = sqsaturationderivative(tp, bt, this%xnew(node))
417  drterm = drterm * this%q_mult(i)
418  !--fill amat and rhs with newton-raphson terms
419  call matrix_sln%add_value_pos(idxglo(ipos), drterm)
420  rhs(node) = rhs(node) + drterm * this%xnew(node)
421  end if
422  end if
423  end do
Here is the call graph for this function:

◆ wel_obs_supported()

logical function welmodule::wel_obs_supported ( class(weltype this)
private

Function to determine if observations are supported by the WEL package. Observations are supported by the WEL package.

Returns
wel_obs_supported boolean indicating if observations are supported
Parameters
thisWelType object

Definition at line 512 of file gwf-wel.f90.

513  ! -- dummy variables
514  class(WelType) :: this !< WelType object
515  !
516  ! -- set boolean
517  wel_obs_supported = .true.

◆ wel_options()

subroutine welmodule::wel_options ( class(weltype), intent(inout)  this)

Source additional options for WEL package.

Parameters
[in,out]thisWelType object

Definition at line 185 of file gwf-wel.f90.

186  ! -- modules
187  use inputoutputmodule, only: urword
190  ! -- dummy variables
191  class(WelType), intent(inout) :: this !< WelType object
192  ! -- local variables
193  character(len=LINELENGTH) :: fname
194  type(GwfWelParamFoundType) :: found
195  ! -- formats
196  character(len=*), parameter :: fmtflowred = &
197  &"(4x, 'AUTOMATIC FLOW REDUCTION OF WELLS IMPLEMENTED.')"
198  character(len=*), parameter :: fmtflowredv = &
199  &"(4x, 'AUTOMATIC FLOW REDUCTION FRACTION (',g15.7,').')"
200  !
201  ! -- source base BndExtType options
202  call this%BndExtType%source_options()
203  !
204  ! -- source well options from input context
205  call mem_set_value(this%flowred, 'FLOWRED', this%input_mempath, found%flowred)
206  call mem_set_value(fname, 'AFRCSVFILE', this%input_mempath, found%afrcsvfile)
207  call mem_set_value(this%imover, 'MOVER', this%input_mempath, found%mover)
208  !
209  if (found%flowred) then
210  !
211  this%iflowred = 1
212  !
213  if (this%flowred <= dzero) then
214  this%flowred = dem1
215  else if (this%flowred > done) then
216  this%flowred = done
217  end if
218  end if
219  !
220  if (found%afrcsvfile) then
221  call this%wel_afr_csv_init(fname)
222  end if
223  !
224  if (found%mover) then
225  this%imover = 1
226  end if
227  !
228  ! -- log WEL specific options
229  call this%log_wel_options(found)
subroutine, public urword(line, icol, istart, istop, ncode, n, r, iout, in)
Extract a word from a string.
Here is the call graph for this function:

◆ wel_rp()

subroutine welmodule::wel_rp ( class(weltype), intent(inout)  this)

Definition at line 273 of file gwf-wel.f90.

274  use tdismodule, only: kper
275  ! -- dummy
276  class(WelType), intent(inout) :: this
277  ! -- local
278  !
279  if (this%iper /= kper) return
280  !
281  ! -- Call the parent class read and prepare
282  call this%BndExtType%bnd_rp()
283  !
284  ! -- Write the list to iout if requested
285  if (this%iprpak /= 0) then
286  call this%write_list()
287  end if

Variable Documentation

◆ ftype

character(len=lenftype) welmodule::ftype = 'WEL'
private

Definition at line 36 of file gwf-wel.f90.

36  character(len=LENFTYPE) :: ftype = 'WEL' !< package ftype

◆ text

character(len=16) welmodule::text = ' WEL'
private

Definition at line 37 of file gwf-wel.f90.

37  character(len=16) :: text = ' WEL' !< package flow text string