MODFLOW 6  version 6.5.0.dev2
MODFLOW 6 Code Documentation
numericalmodelmodule Module Reference

Data Types

type  numericalmodeltype
 

Functions/Subroutines

subroutine model_df (this)
 
subroutine model_ac (this, sparse)
 
subroutine model_mc (this, matrix_sln)
 
subroutine model_ar (this)
 
subroutine model_rp (this)
 
subroutine model_ad (this)
 
subroutine model_reset (this)
 
subroutine model_solve (this)
 
subroutine model_cf (this, kiter)
 
subroutine model_fc (this, kiter, matrix_sln, inwtflag)
 
subroutine model_ptcchk (this, iptc)
 
subroutine model_ptc (this, vec_residual, iptc, ptcf)
 
subroutine model_nr (this, kiter, matrix, inwtflag)
 
subroutine model_cc (this, innertot, kiter, iend, icnvgmod, cpak, ipak, dpak)
 
subroutine model_nur (this, neqmod, x, xtemp, dx, inewtonur, dxmax, locmax)
 
subroutine model_cq (this, icnvg, isuppress_output)
 
subroutine model_bd (this, icnvg, isuppress_output)
 
subroutine model_bdcalc (this, icnvg)
 
subroutine model_bdsave (this, icnvg)
 
subroutine model_ot (this)
 
subroutine model_bdentry (this, budterm, budtxt, rowlabel)
 
subroutine model_fp (this)
 
subroutine model_da (this)
 
subroutine set_moffset (this, moffset)
 
subroutine get_mrange (this, mstart, mend)
 
subroutine set_idsoln (this, id)
 
subroutine allocate_scalars (this, modelname)
 
subroutine allocate_arrays (this)
 
subroutine set_xptr (this, xsln, sln_offset, varNameTgt, memPathTgt)
 
subroutine set_rhsptr (this, rhssln, sln_offset, varNameTgt, memPathTgt)
 
subroutine set_iboundptr (this, iboundsln, sln_offset, varNameTgt, memPathTgt)
 
subroutine get_mcellid (this, node, mcellid)
 
subroutine get_mnodeu (this, node, nodeu)
 
integer(i4b) function get_iasym (this)
 
class(numericalmodeltype) function, pointer castasnumericalmodelclass (obj)
 
subroutine, public addnumericalmodeltolist (list, model)
 
class(numericalmodeltype) function, pointer, public getnumericalmodelfromlist (list, idx)
 
subroutine create_lstfile (this, lst_fname, model_fname, defined, headertxt)
 

Function/Subroutine Documentation

◆ addnumericalmodeltolist()

subroutine, public numericalmodelmodule::addnumericalmodeltolist ( type(listtype), intent(inout)  list,
class(numericalmodeltype), intent(inout), pointer  model 
)

Definition at line 441 of file NumericalModel.f90.

442  implicit none
443  ! -- dummy
444  type(ListType), intent(inout) :: list
445  class(NumericalModelType), pointer, intent(inout) :: model
446  ! -- local
447  class(*), pointer :: obj
448  !
449  obj => model
450  call list%Add(obj)
451  !
452  return
Here is the caller graph for this function:

◆ allocate_arrays()

subroutine numericalmodelmodule::allocate_arrays ( class(numericalmodeltype this)

Definition at line 302 of file NumericalModel.f90.

303  use constantsmodule, only: dzero
305  class(NumericalModelType) :: this
306  integer(I4B) :: i
307  !
308  call mem_allocate(this%xold, this%neq, 'XOLD', this%memoryPath)
309  call mem_allocate(this%flowja, this%nja, 'FLOWJA', this%memoryPath)
310  call mem_allocate(this%idxglo, this%nja, 'IDXGLO', this%memoryPath)
311  !
312  ! -- initialize
313  do i = 1, size(this%flowja)
314  this%flowja(i) = dzero
315  end do
316  !
317  ! -- return
318  return
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:64

◆ allocate_scalars()

subroutine numericalmodelmodule::allocate_scalars ( class(numericalmodeltype this,
character(len=*), intent(in)  modelname 
)
private

Definition at line 276 of file NumericalModel.f90.

278  class(NumericalModelType) :: this
279  character(len=*), intent(in) :: modelname
280  !
281  ! -- allocate basetype members
282  call this%BaseModelType%allocate_scalars(modelname)
283  !
284  ! -- allocate members from this type
285  call mem_allocate(this%neq, 'NEQ', this%memoryPath)
286  call mem_allocate(this%nja, 'NJA', this%memoryPath)
287  call mem_allocate(this%icnvg, 'ICNVG', this%memoryPath)
288  call mem_allocate(this%moffset, 'MOFFSET', this%memoryPath)
289  allocate (this%filename)
290  allocate (this%bndlist)
291  !
292  this%filename = ''
293  this%neq = 0
294  this%nja = 0
295  this%icnvg = 0
296  this%moffset = 0
297  !
298  ! -- return
299  return

◆ castasnumericalmodelclass()

class(numericalmodeltype) function, pointer numericalmodelmodule::castasnumericalmodelclass ( class(*), intent(inout), pointer  obj)
private

Definition at line 426 of file NumericalModel.f90.

427  implicit none
428  class(*), pointer, intent(inout) :: obj
429  class(NumericalModelType), pointer :: res
430  !
431  res => null()
432  if (.not. associated(obj)) return
433  !
434  select type (obj)
435  class is (numericalmodeltype)
436  res => obj
437  end select
438  return
Here is the caller graph for this function:

◆ create_lstfile()

subroutine numericalmodelmodule::create_lstfile ( class(numericalmodeltype this,
character(len=*), intent(inout)  lst_fname,
character(len=*), intent(in)  model_fname,
logical(lgp), intent(in)  defined,
character(len=*), intent(in)  headertxt 
)
private

Definition at line 470 of file NumericalModel.f90.

471  ! -- modules
472  use kindmodule, only: lgp
474  ! -- dummy
475  class(NumericalModelType) :: this
476  character(len=*), intent(inout) :: lst_fname
477  character(len=*), intent(in) :: model_fname
478  logical(LGP), intent(in) :: defined
479  character(len=*), intent(in) :: headertxt
480  ! -- local
481  integer(I4B) :: i, istart, istop
482  !
483  ! -- set list file name if not provided
484  if (.not. defined) then
485  !
486  ! -- initialize
487  lst_fname = ' '
488  istart = 0
489  istop = len_trim(model_fname)
490  !
491  ! -- identify '.' character position from back of string
492  do i = istop, 1, -1
493  if (model_fname(i:i) == '.') then
494  istart = i
495  exit
496  end if
497  end do
498  !
499  ! -- if not found start from string end
500  if (istart == 0) istart = istop + 1
501  !
502  ! -- set list file name
503  lst_fname = model_fname(1:istart)
504  istop = istart + 3
505  lst_fname(istart:istop) = '.lst'
506  end if
507  !
508  ! -- create the list file
509  this%iout = getunit()
510  call openfile(this%iout, 0, lst_fname, 'LIST', filstat_opt='REPLACE')
511  !
512  ! -- write list file header
513  call write_listfile_header(this%iout, headertxt)
514  !
515  ! -- return
516  return
integer(i4b) function, public getunit()
Get a free unit number.
subroutine, public openfile(iu, iout, fname, ftype, fmtarg_opt, accarg_opt, filstat_opt, mode_opt)
Open a file.
Definition: InputOutput.f90:30
This module defines variable data types.
Definition: kind.f90:8
Here is the call graph for this function:

◆ get_iasym()

integer(i4b) function numericalmodelmodule::get_iasym ( class(numericalmodeltype this)

Definition at line 420 of file NumericalModel.f90.

421  class(NumericalModelType) :: this
422  integer(I4B) :: iasym
423  iasym = 0

◆ get_mcellid()

subroutine numericalmodelmodule::get_mcellid ( class(numericalmodeltype this,
integer(i4b), intent(in)  node,
character(len=*), intent(inout)  mcellid 
)

Definition at line 370 of file NumericalModel.f90.

371  use bndmodule, only: bndtype, getbndfromlist
372  class(NumericalModelType) :: this
373  integer(I4B), intent(in) :: node
374  character(len=*), intent(inout) :: mcellid
375  ! -- local
376  character(len=20) :: cellid
377  integer(I4B) :: ip, ipaknode, istart, istop
378  class(BndType), pointer :: packobj
379 
380  if (node < 1) then
381  cellid = ''
382  else if (node <= this%dis%nodes) then
383  call this%dis%noder_to_string(node, cellid)
384  else
385  cellid = '***ERROR***'
386  ipaknode = node - this%dis%nodes
387  istart = 1
388  do ip = 1, this%bndlist%Count()
389  packobj => getbndfromlist(this%bndlist, ip)
390  if (packobj%npakeq == 0) cycle
391  istop = istart + packobj%npakeq - 1
392  if (istart <= ipaknode .and. ipaknode <= istop) then
393  write (cellid, '(a, a, a, i0, a, i0, a)') '(', &
394  trim(packobj%filtyp), '_', &
395  packobj%ibcnum, '-', ipaknode - packobj%ioffset, ')'
396  exit
397  end if
398  istart = istop + 1
399  end do
400  end if
401  write (mcellid, '(i0, a, a, a, a)') this%id, '_', this%macronym, '-', &
402  trim(adjustl(cellid))
403  return
This module contains the base boundary package.
class(bndtype) function, pointer, public getbndfromlist(list, idx)
Get boundary from package list.
@ brief BndType
Here is the call graph for this function:

◆ get_mnodeu()

subroutine numericalmodelmodule::get_mnodeu ( class(numericalmodeltype this,
integer(i4b), intent(in)  node,
integer(i4b), intent(inout)  nodeu 
)

Definition at line 406 of file NumericalModel.f90.

407  use bndmodule, only: bndtype, getbndfromlist
408  class(NumericalModelType) :: this
409  integer(I4B), intent(in) :: node
410  integer(I4B), intent(inout) :: nodeu
411  ! -- local
412  if (node <= this%dis%nodes) then
413  nodeu = this%dis%get_nodeuser(node)
414  else
415  nodeu = -(node - this%dis%nodes)
416  end if
417  return
Here is the call graph for this function:

◆ get_mrange()

subroutine numericalmodelmodule::get_mrange ( class(numericalmodeltype this,
integer(i4b), intent(inout)  mstart,
integer(i4b), intent(inout)  mend 
)
private

Definition at line 262 of file NumericalModel.f90.

263  class(NumericalModelType) :: this
264  integer(I4B), intent(inout) :: mstart
265  integer(I4B), intent(inout) :: mend
266  mstart = this%moffset + 1
267  mend = mstart + this%neq - 1

◆ getnumericalmodelfromlist()

class(numericalmodeltype) function, pointer, public numericalmodelmodule::getnumericalmodelfromlist ( type(listtype), intent(inout)  list,
integer(i4b), intent(in)  idx 
)

Definition at line 455 of file NumericalModel.f90.

456  implicit none
457  ! -- dummy
458  type(ListType), intent(inout) :: list
459  integer(I4B), intent(in) :: idx
460  class(NumericalModelType), pointer :: res
461  ! -- local
462  class(*), pointer :: obj
463  !
464  obj => list%GetItem(idx)
465  res => castasnumericalmodelclass(obj)
466  !
467  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ model_ac()

subroutine numericalmodelmodule::model_ac ( class(numericalmodeltype this,
type(sparsematrix), intent(inout)  sparse 
)
private

Definition at line 91 of file NumericalModel.f90.

92  class(NumericalModelType) :: this
93  type(sparsematrix), intent(inout) :: sparse

◆ model_ad()

subroutine numericalmodelmodule::model_ad ( class(numericalmodeltype this)
private

Definition at line 109 of file NumericalModel.f90.

110  class(NumericalModelType) :: this

◆ model_ar()

subroutine numericalmodelmodule::model_ar ( class(numericalmodeltype this)
private

Definition at line 101 of file NumericalModel.f90.

102  class(NumericalModelType) :: this

◆ model_bd()

subroutine numericalmodelmodule::model_bd ( class(numericalmodeltype this,
integer(i4b), intent(in)  icnvg,
integer(i4b), intent(in)  isuppress_output 
)
private

Definition at line 191 of file NumericalModel.f90.

192  class(NumericalModelType) :: this
193  integer(I4B), intent(in) :: icnvg
194  integer(I4B), intent(in) :: isuppress_output

◆ model_bdcalc()

subroutine numericalmodelmodule::model_bdcalc ( class(numericalmodeltype this,
integer(i4b), intent(in)  icnvg 
)
private

Definition at line 197 of file NumericalModel.f90.

198  class(NumericalModelType) :: this
199  integer(I4B), intent(in) :: icnvg

◆ model_bdentry()

subroutine numericalmodelmodule::model_bdentry ( class(numericalmodeltype this,
real(dp), dimension(:, :), intent(in)  budterm,
character(len=lenbudtxt), dimension(:), intent(in)  budtxt,
character(len=*), intent(in)  rowlabel 
)
private

Definition at line 211 of file NumericalModel.f90.

212  class(NumericalModelType) :: this
213  real(DP), dimension(:, :), intent(in) :: budterm
214  character(len=LENBUDTXT), dimension(:), intent(in) :: budtxt
215  character(len=*), intent(in) :: rowlabel

◆ model_bdsave()

subroutine numericalmodelmodule::model_bdsave ( class(numericalmodeltype this,
integer(i4b), intent(in)  icnvg 
)
private

Definition at line 202 of file NumericalModel.f90.

203  class(NumericalModelType) :: this
204  integer(I4B), intent(in) :: icnvg

◆ model_cc()

subroutine numericalmodelmodule::model_cc ( class(numericalmodeltype this,
integer(i4b), intent(in)  innertot,
integer(i4b), intent(in)  kiter,
integer(i4b), intent(in)  iend,
integer(i4b), intent(in)  icnvgmod,
character(len=lenpakloc), intent(inout)  cpak,
integer(i4b), intent(inout)  ipak,
real(dp), intent(inout)  dpak 
)
private

Definition at line 163 of file NumericalModel.f90.

164  class(NumericalModelType) :: this
165  integer(I4B), intent(in) :: innertot
166  integer(I4B), intent(in) :: kiter
167  integer(I4B), intent(in) :: iend
168  integer(I4B), intent(in) :: icnvgmod
169  character(len=LENPAKLOC), intent(inout) :: cpak
170  integer(I4B), intent(inout) :: ipak
171  real(DP), intent(inout) :: dpak

◆ model_cf()

subroutine numericalmodelmodule::model_cf ( class(numericalmodeltype this,
integer(i4b), intent(in)  kiter 
)
private

Definition at line 131 of file NumericalModel.f90.

132  class(NumericalModelType) :: this
133  integer(I4B), intent(in) :: kiter

◆ model_cq()

subroutine numericalmodelmodule::model_cq ( class(numericalmodeltype this,
integer(i4b), intent(in)  icnvg,
integer(i4b), intent(in)  isuppress_output 
)
private

Definition at line 185 of file NumericalModel.f90.

186  class(NumericalModelType) :: this
187  integer(I4B), intent(in) :: icnvg
188  integer(I4B), intent(in) :: isuppress_output

◆ model_da()

subroutine numericalmodelmodule::model_da ( class(numericalmodeltype this)
private

Definition at line 222 of file NumericalModel.f90.

223  ! -- modules
225  class(NumericalModelType) :: this
226 
227  ! -- Scalars
228  call mem_deallocate(this%neq)
229  call mem_deallocate(this%nja)
230  call mem_deallocate(this%icnvg)
231  call mem_deallocate(this%moffset)
232  deallocate (this%filename)
233  !
234  ! -- Arrays
235  call mem_deallocate(this%xold)
236  call mem_deallocate(this%flowja)
237  call mem_deallocate(this%idxglo)
238  !
239  ! -- derived types
240  call this%bndlist%Clear()
241  deallocate (this%bndlist)
242  !
243  ! -- nullify pointers
244  call mem_deallocate(this%x, 'X', this%memoryPath)
245  call mem_deallocate(this%rhs, 'RHS', this%memoryPath)
246  call mem_deallocate(this%ibound, 'IBOUND', this%memoryPath)
247  !
248  ! -- BaseModelType
249  call this%BaseModelType%model_da()
250  !
251  !
252  ! -- Return
253  return

◆ model_df()

subroutine numericalmodelmodule::model_df ( class(numericalmodeltype this)
private

Definition at line 87 of file NumericalModel.f90.

88  class(NumericalModelType) :: this

◆ model_fc()

subroutine numericalmodelmodule::model_fc ( class(numericalmodeltype this,
integer(i4b), intent(in)  kiter,
class(matrixbasetype), pointer  matrix_sln,
integer(i4b), intent(in)  inwtflag 
)
private

Definition at line 136 of file NumericalModel.f90.

137  class(NumericalModelType) :: this
138  integer(I4B), intent(in) :: kiter
139  class(MatrixBaseType), pointer :: matrix_sln
140  integer(I4B), intent(in) :: inwtflag

◆ model_fp()

subroutine numericalmodelmodule::model_fp ( class(numericalmodeltype this)
private

Definition at line 218 of file NumericalModel.f90.

219  class(NumericalModelType) :: this

◆ model_mc()

subroutine numericalmodelmodule::model_mc ( class(numericalmodeltype this,
class(matrixbasetype), pointer  matrix_sln 
)
private

Definition at line 96 of file NumericalModel.f90.

97  class(NumericalModelType) :: this
98  class(MatrixBaseType), pointer :: matrix_sln

◆ model_nr()

subroutine numericalmodelmodule::model_nr ( class(numericalmodeltype this,
integer(i4b), intent(in)  kiter,
class(matrixbasetype), pointer  matrix,
integer(i4b), intent(in)  inwtflag 
)
private

Definition at line 156 of file NumericalModel.f90.

157  class(NumericalModelType) :: this
158  integer(I4B), intent(in) :: kiter
159  class(MatrixBaseType), pointer :: matrix
160  integer(I4B), intent(in) :: inwtflag

◆ model_nur()

subroutine numericalmodelmodule::model_nur ( class(numericalmodeltype this,
integer(i4b), intent(in)  neqmod,
real(dp), dimension(neqmod), intent(inout)  x,
real(dp), dimension(neqmod), intent(in)  xtemp,
real(dp), dimension(neqmod), intent(inout)  dx,
integer(i4b), intent(inout)  inewtonur,
real(dp), intent(inout)  dxmax,
integer(i4b), intent(inout)  locmax 
)
private

Definition at line 174 of file NumericalModel.f90.

175  class(NumericalModelType) :: this
176  integer(I4B), intent(in) :: neqmod
177  real(DP), dimension(neqmod), intent(inout) :: x
178  real(DP), dimension(neqmod), intent(in) :: xtemp
179  real(DP), dimension(neqmod), intent(inout) :: dx
180  integer(I4B), intent(inout) :: inewtonur
181  real(DP), intent(inout) :: dxmax
182  integer(I4B), intent(inout) :: locmax

◆ model_ot()

subroutine numericalmodelmodule::model_ot ( class(numericalmodeltype this)
private

Definition at line 207 of file NumericalModel.f90.

208  class(NumericalModelType) :: this

◆ model_ptc()

subroutine numericalmodelmodule::model_ptc ( class(numericalmodeltype this,
class(vectorbasetype), pointer  vec_residual,
integer(i4b), intent(inout)  iptc,
real(dp), intent(inout)  ptcf 
)
private

Definition at line 149 of file NumericalModel.f90.

150  class(NumericalModelType) :: this
151  class(VectorBaseType), pointer :: vec_residual
152  integer(I4B), intent(inout) :: iptc
153  real(DP), intent(inout) :: ptcf

◆ model_ptcchk()

subroutine numericalmodelmodule::model_ptcchk ( class(numericalmodeltype this,
integer(i4b), intent(inout)  iptc 
)
private

Definition at line 143 of file NumericalModel.f90.

144  class(NumericalModelType) :: this
145  integer(I4B), intent(inout) :: iptc
146  iptc = 0

◆ model_reset()

subroutine numericalmodelmodule::model_reset ( class(numericalmodeltype this)
private

Definition at line 113 of file NumericalModel.f90.

114  use bndmodule, only: bndtype, getbndfromlist
115  class(NumericalModelType) :: this
116  ! local
117  class(BndType), pointer :: packobj
118  integer(I4B) :: ip
119 
120  do ip = 1, this%bndlist%Count()
121  packobj => getbndfromlist(this%bndlist, ip)
122  call packobj%bnd_reset()
123  end do
124 
Here is the call graph for this function:

◆ model_rp()

subroutine numericalmodelmodule::model_rp ( class(numericalmodeltype this)
private

Definition at line 105 of file NumericalModel.f90.

106  class(NumericalModelType) :: this

◆ model_solve()

subroutine numericalmodelmodule::model_solve ( class(numericalmodeltype this)

Definition at line 127 of file NumericalModel.f90.

128  class(NumericalModelType) :: this

◆ set_iboundptr()

subroutine numericalmodelmodule::set_iboundptr ( class(numericalmodeltype this,
integer(i4b), dimension(:), intent(in), pointer, contiguous  iboundsln,
integer(i4b)  sln_offset,
character(len=*), intent(in)  varNameTgt,
character(len=*), intent(in)  memPathTgt 
)

Definition at line 353 of file NumericalModel.f90.

355  ! -- dummy
356  class(NumericalModelType) :: this
357  integer(I4B), dimension(:), pointer, contiguous, intent(in) :: iboundsln
358  integer(I4B) :: sln_offset
359  character(len=*), intent(in) :: varNameTgt
360  character(len=*), intent(in) :: memPathTgt
361  ! -- local
362  integer(I4B) :: offset
363  ! -- code
364  offset = this%moffset - sln_offset
365  this%ibound => iboundsln(offset + 1:offset + this%neq)
366  call mem_checkin(this%ibound, 'IBOUND', this%memoryPath, varnametgt, &
367  mempathtgt)

◆ set_idsoln()

subroutine numericalmodelmodule::set_idsoln ( class(numericalmodeltype this,
integer(i4b), intent(in)  id 
)
private

Definition at line 270 of file NumericalModel.f90.

271  class(NumericalModelType) :: this
272  integer(I4B), intent(in) :: id
273  this%idsoln = id

◆ set_moffset()

subroutine numericalmodelmodule::set_moffset ( class(numericalmodeltype this,
integer(i4b), intent(in)  moffset 
)

Definition at line 256 of file NumericalModel.f90.

257  class(NumericalModelType) :: this
258  integer(I4B), intent(in) :: moffset
259  this%moffset = moffset

◆ set_rhsptr()

subroutine numericalmodelmodule::set_rhsptr ( class(numericalmodeltype this,
real(dp), dimension(:), intent(in), pointer, contiguous  rhssln,
integer(i4b)  sln_offset,
character(len=*), intent(in)  varNameTgt,
character(len=*), intent(in)  memPathTgt 
)

Definition at line 337 of file NumericalModel.f90.

339  ! -- dummy
340  class(NumericalModelType) :: this
341  real(DP), dimension(:), pointer, contiguous, intent(in) :: rhssln
342  integer(I4B) :: sln_offset
343  character(len=*), intent(in) :: varNameTgt
344  character(len=*), intent(in) :: memPathTgt
345  ! -- local
346  integer(I4B) :: offset
347  ! -- code
348  offset = this%moffset - sln_offset
349  this%rhs => rhssln(offset + 1:offset + this%neq)
350  call mem_checkin(this%rhs, 'RHS', this%memoryPath, varnametgt, mempathtgt)

◆ set_xptr()

subroutine numericalmodelmodule::set_xptr ( class(numericalmodeltype this,
real(dp), dimension(:), intent(in), pointer, contiguous  xsln,
integer(i4b)  sln_offset,
character(len=*), intent(in)  varNameTgt,
character(len=*), intent(in)  memPathTgt 
)

Definition at line 321 of file NumericalModel.f90.

323  ! -- dummy
324  class(NumericalModelType) :: this
325  real(DP), dimension(:), pointer, contiguous, intent(in) :: xsln
326  integer(I4B) :: sln_offset
327  character(len=*), intent(in) :: varNameTgt
328  character(len=*), intent(in) :: memPathTgt
329  ! -- local
330  integer(I4B) :: offset
331  ! -- code
332  offset = this%moffset - sln_offset
333  this%x => xsln(offset + 1:offset + this%neq)
334  call mem_checkin(this%x, 'X', this%memoryPath, varnametgt, mempathtgt)