MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
prtmodule Module Reference

Data Types

type  prtmodeltype
 Particle tracking (PRT) model. More...
 

Functions/Subroutines

subroutine, public prt_cr (filename, id, modelname)
 Create a new particle tracking model object. More...
 
subroutine prt_df (this)
 Define packages. More...
 
subroutine prt_ar (this)
 Allocate and read. More...
 
subroutine prt_rp (this)
 Read and prepare (calls package read and prepare routines) More...
 
subroutine prt_ad (this)
 Time step advance (calls package advance subroutines) More...
 
subroutine prt_cq (this, icnvg, isuppress_output)
 Calculate intercell flow (flowja) More...
 
subroutine prt_cq_sto (this)
 Calculate particle mass storage. More...
 
subroutine prt_bd (this, icnvg, isuppress_output)
 Calculate flows and budget. More...
 
subroutine prt_ot (this)
 Print and/or save model output. More...
 
subroutine prt_ot_flow (this, icbcfl, ibudfl, icbcun)
 Save flows. More...
 
subroutine prt_ot_saveflow (this, nja, flowja, icbcfl, icbcun)
 Save intercell flows. More...
 
subroutine prt_ot_printflow (this, ibudfl, flowja)
 Print intercell flows. More...
 
subroutine prt_ot_dv (this, idvsave, idvprint, ipflag)
 Print dependent variables. More...
 
subroutine prt_ot_bdsummary (this, ibudfl, ipflag)
 Print budget summary. More...
 
subroutine prt_da (this)
 Deallocate. More...
 
subroutine allocate_scalars (this, modelname)
 Allocate memory for scalars. More...
 
subroutine allocate_arrays (this)
 Allocate arrays. More...
 
subroutine package_create (this, filtyp, ipakid, ipaknum, pakname, mempath, inunit, iout)
 Create boundary condition packages for this model. More...
 
subroutine ftype_check (this, indis)
 Check to make sure required input files have been specified. More...
 
subroutine prt_solve (this)
 Solve the model. More...
 
subroutine create_bndpkgs (this, bndpkgs, pkgtypes, pkgnames, mempaths, inunits)
 Source package info and begin to process. More...
 
subroutine create_packages (this)
 Source package info and begin to process. More...
 
subroutine log_namfile_options (this, found)
 Write model namfile options to list file. More...
 

Variables

integer(i4b), parameter nbditems = 1
 
character(len=lenbudtxt), dimension(nbditemsbudtxt
 
integer(i4b), parameter, public prt_nbasepkg = 50
 PRT base package array descriptors. More...
 
character(len=lenpackagetype), dimension(prt_nbasepkg), public prt_basepkg
 
integer(i4b), parameter, public prt_nmultipkg = 50
 PRT multi package array descriptors. More...
 
character(len=lenpackagetype), dimension(prt_nmultipkg), public prt_multipkg
 
integer(i4b), parameter niunit_prt = PRT_NBASEPKG + PRT_NMULTIPKG
 

Function/Subroutine Documentation

◆ allocate_arrays()

subroutine prtmodule::allocate_arrays ( class(prtmodeltype this)
private

Definition at line 811 of file prt.f90.

813  class(PrtModelType) :: this
814  integer(I4B) :: n
815 
816  ! Allocate arrays in parent type
817  this%nja = this%dis%nja
818  call this%NumericalModelType%allocate_arrays()
819 
820  ! Allocate and initialize arrays
821  call mem_allocate(this%masssto, this%dis%nodes, &
822  'MASSSTO', this%memoryPath)
823  call mem_allocate(this%massstoold, this%dis%nodes, &
824  'MASSSTOOLD', this%memoryPath)
825  call mem_allocate(this%ratesto, this%dis%nodes, &
826  'RATESTO', this%memoryPath)
827  ! explicit model, so these must be manually allocated
828  call mem_allocate(this%x, this%dis%nodes, 'X', this%memoryPath)
829  call mem_allocate(this%rhs, this%dis%nodes, 'RHS', this%memoryPath)
830  call mem_allocate(this%ibound, this%dis%nodes, 'IBOUND', this%memoryPath)
831  do n = 1, this%dis%nodes
832  this%masssto(n) = dzero
833  this%massstoold(n) = dzero
834  this%ratesto(n) = dzero
835  this%x(n) = dzero
836  this%rhs(n) = dzero
837  this%ibound(n) = 1
838  end do

◆ allocate_scalars()

subroutine prtmodule::allocate_scalars ( class(prtmodeltype this,
character(len=*), intent(in)  modelname 
)

Definition at line 782 of file prt.f90.

783  ! dummy
784  class(PrtModelType) :: this
785  character(len=*), intent(in) :: modelname
786 
787  ! allocate members from parent class
788  call this%NumericalModelType%allocate_scalars(modelname)
789 
790  ! allocate members that are part of model class
791  call mem_allocate(this%infmi, 'INFMI', this%memoryPath)
792  call mem_allocate(this%inmip, 'INMIP', this%memoryPath)
793  call mem_allocate(this%inmvt, 'INMVT', this%memoryPath)
794  call mem_allocate(this%inmst, 'INMST', this%memoryPath)
795  call mem_allocate(this%inadv, 'INADV', this%memoryPath)
796  call mem_allocate(this%indsp, 'INDSP', this%memoryPath)
797  call mem_allocate(this%inssm, 'INSSM', this%memoryPath)
798  call mem_allocate(this%inoc, 'INOC ', this%memoryPath)
799 
800  this%infmi = 0
801  this%inmip = 0
802  this%inmvt = 0
803  this%inmst = 0
804  this%inadv = 0
805  this%indsp = 0
806  this%inssm = 0
807  this%inoc = 0

◆ create_bndpkgs()

subroutine prtmodule::create_bndpkgs ( class(prtmodeltype this,
integer(i4b), dimension(:), intent(inout), allocatable  bndpkgs,
type(characterstringtype), dimension(:), intent(inout), pointer, contiguous  pkgtypes,
type(characterstringtype), dimension(:), intent(inout), pointer, contiguous  pkgnames,
type(characterstringtype), dimension(:), intent(inout), pointer, contiguous  mempaths,
integer(i4b), dimension(:), intent(inout), pointer, contiguous  inunits 
)

Definition at line 994 of file prt.f90.

996  ! modules
999  ! dummy
1000  class(PrtModelType) :: this
1001  integer(I4B), dimension(:), allocatable, intent(inout) :: bndpkgs
1002  type(CharacterStringType), dimension(:), contiguous, &
1003  pointer, intent(inout) :: pkgtypes
1004  type(CharacterStringType), dimension(:), contiguous, &
1005  pointer, intent(inout) :: pkgnames
1006  type(CharacterStringType), dimension(:), contiguous, &
1007  pointer, intent(inout) :: mempaths
1008  integer(I4B), dimension(:), contiguous, &
1009  pointer, intent(inout) :: inunits
1010  ! local
1011  integer(I4B) :: ipakid, ipaknum
1012  character(len=LENFTYPE) :: pkgtype, bndptype
1013  character(len=LENPACKAGENAME) :: pkgname
1014  character(len=LENMEMPATH) :: mempath
1015  integer(I4B), pointer :: inunit
1016  integer(I4B) :: n
1017 
1018  if (allocated(bndpkgs)) then
1019  ! create stress packages
1020  ipakid = 1
1021  bndptype = ''
1022  do n = 1, size(bndpkgs)
1023  pkgtype = pkgtypes(bndpkgs(n))
1024  pkgname = pkgnames(bndpkgs(n))
1025  mempath = mempaths(bndpkgs(n))
1026  inunit => inunits(bndpkgs(n))
1027 
1028  if (bndptype /= pkgtype) then
1029  ipaknum = 1
1030  bndptype = pkgtype
1031  end if
1032 
1033  call this%package_create(pkgtype, ipakid, ipaknum, pkgname, mempath, &
1034  inunit, this%iout)
1035  ipakid = ipakid + 1
1036  ipaknum = ipaknum + 1
1037  end do
1038 
1039  ! cleanup
1040  deallocate (bndpkgs)
1041  end if
1042 
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
integer(i4b), parameter lenpackagename
maximum length of the package name
Definition: Constants.f90:23
This class is used to store a single deferred-length character string. It was designed to work in an ...
Definition: CharString.f90:23

◆ create_packages()

subroutine prtmodule::create_packages ( class(prtmodeltype this)

Definition at line 1046 of file prt.f90.

1047  ! modules
1050  use arrayhandlersmodule, only: expandarray
1051  use memorymanagermodule, only: mem_setptr
1053  use simvariablesmodule, only: idm_context
1054  use budgetmodule, only: budget_cr
1058  use prtmipmodule, only: mip_cr
1059  use prtfmimodule, only: fmi_cr
1060  use prtocmodule, only: oc_cr
1061  ! dummy
1062  class(PrtModelType) :: this
1063  ! local
1064  type(CharacterStringType), dimension(:), contiguous, &
1065  pointer :: pkgtypes => null()
1066  type(CharacterStringType), dimension(:), contiguous, &
1067  pointer :: pkgnames => null()
1068  type(CharacterStringType), dimension(:), contiguous, &
1069  pointer :: mempaths => null()
1070  integer(I4B), dimension(:), contiguous, &
1071  pointer :: inunits => null()
1072  character(len=LENMEMPATH) :: model_mempath
1073  character(len=LENFTYPE) :: pkgtype
1074  character(len=LENPACKAGENAME) :: pkgname
1075  character(len=LENMEMPATH) :: mempath
1076  integer(I4B), pointer :: inunit
1077  integer(I4B), dimension(:), allocatable :: bndpkgs
1078  integer(I4B) :: n
1079  integer(I4B) :: indis = 0 ! DIS enabled flag
1080  character(len=LENMEMPATH) :: mempathmip = ''
1081  character(len=LENMEMPATH) :: mempathfmi = ''
1082 
1083  ! set input memory paths, input/model and input/model/namfile
1084  model_mempath = create_mem_path(component=this%name, context=idm_context)
1085 
1086  ! set pointers to model path package info
1087  call mem_setptr(pkgtypes, 'PKGTYPES', model_mempath)
1088  call mem_setptr(pkgnames, 'PKGNAMES', model_mempath)
1089  call mem_setptr(mempaths, 'MEMPATHS', model_mempath)
1090  call mem_setptr(inunits, 'INUNITS', model_mempath)
1091 
1092  do n = 1, size(pkgtypes)
1093  ! attributes for this input package
1094  pkgtype = pkgtypes(n)
1095  pkgname = pkgnames(n)
1096  mempath = mempaths(n)
1097  inunit => inunits(n)
1098 
1099  ! create dis package first as it is a prerequisite for other packages
1100  select case (pkgtype)
1101  case ('DIS6')
1102  indis = 1
1103  call dis_cr(this%dis, this%name, mempath, indis, this%iout)
1104  case ('DISV6')
1105  indis = 1
1106  call disv_cr(this%dis, this%name, mempath, indis, this%iout)
1107  case ('DISU6')
1108  indis = 1
1109  call disu_cr(this%dis, this%name, mempath, indis, this%iout)
1110  case ('MIP6')
1111  this%inmip = 1
1112  mempathmip = mempath
1113  case ('FMI6')
1114  this%infmi = 1
1115  mempathfmi = mempath
1116  case ('OC6')
1117  this%inoc = inunit
1118  case ('PRP6')
1119  call expandarray(bndpkgs)
1120  bndpkgs(size(bndpkgs)) = n
1121  case default
1122  call pstop(1, "Unrecognized package type: "//pkgtype)
1123  end select
1124  end do
1125 
1126  ! Create budget manager
1127  call budget_cr(this%budget, this%name)
1128 
1129  ! Create tracking method pools
1130  call create_method_pool()
1133 
1134  ! Create packages that are tied directly to model
1135  call mip_cr(this%mip, this%name, mempathmip, this%inmip, this%iout, this%dis)
1136  call fmi_cr(this%fmi, this%name, mempathfmi, this%infmi, this%iout)
1137  call oc_cr(this%oc, this%name, this%inoc, this%iout)
1138 
1139  ! Check to make sure that required ftype's have been specified
1140  call this%ftype_check(indis)
1141 
1142  ! Create boundary packages
1143  call this%create_bndpkgs(bndpkgs, pkgtypes, pkgnames, mempaths, inunits)
This module contains the BudgetModule.
Definition: Budget.f90:20
subroutine, public budget_cr(this, name_model)
@ brief Create a new budget object
Definition: Budget.f90:84
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
Cell-level tracking methods.
subroutine, public create_method_cell_pool()
Create the cell method pool.
Model-level tracking methods.
Definition: MethodPool.f90:2
subroutine, public create_method_pool()
Create the method pool.
Definition: MethodPool.f90:18
Subcell-level tracking methods.
subroutine, public create_method_subcell_pool()
Create the subcell method pool.
subroutine, public fmi_cr(fmiobj, name_model, input_mempath, inunit, iout)
Create a new PrtFmi object.
Definition: prt-fmi.f90:44
subroutine, public mip_cr(mip, name_model, input_mempath, inunit, iout, dis)
Create a model input object.
Definition: prt-mip.f90:34
subroutine, public oc_cr(ocobj, name_model, inunit, iout)
@ brief Create an output control object
Definition: prt-oc.f90:53
This module contains simulation variables.
Definition: SimVariables.f90:9
character(len=linelength) idm_context
Here is the call graph for this function:

◆ ftype_check()

subroutine prtmodule::ftype_check ( class(prtmodeltype this,
integer(i4b), intent(in)  indis 
)

Definition at line 891 of file prt.f90.

892  ! dummy
893  class(PrtModelType) :: this
894  integer(I4B), intent(in) :: indis
895  ! local
896  character(len=LINELENGTH) :: errmsg
897 
898  ! Check for DIS(u) and MIP. Stop if not present.
899  if (indis == 0) then
900  write (errmsg, '(1x,a)') &
901  'Discretization (DIS6, DISV6, or DISU6) package not specified.'
902  call store_error(errmsg)
903  end if
904  if (this%inmip == 0) then
905  write (errmsg, '(1x,a)') &
906  'Model input (MIP6) package not specified.'
907  call store_error(errmsg)
908  end if
909 
910  if (count_errors() > 0) then
911  write (errmsg, '(1x,a)') 'One or more required package(s) not specified.'
912  call store_error(errmsg)
913  call store_error_filename(this%filename)
914  end if
Here is the call graph for this function:

◆ log_namfile_options()

subroutine prtmodule::log_namfile_options ( class(prtmodeltype this,
type(prtnamparamfoundtype), intent(in)  found 
)

Definition at line 1147 of file prt.f90.

1149  class(PrtModelType) :: this
1150  type(PrtNamParamFoundType), intent(in) :: found
1151 
1152  write (this%iout, '(1x,a)') 'NAMEFILE OPTIONS:'
1153 
1154  if (found%print_input) then
1155  write (this%iout, '(4x,a)') 'STRESS PACKAGE INPUT WILL BE PRINTED '// &
1156  'FOR ALL MODEL STRESS PACKAGES'
1157  end if
1158 
1159  if (found%print_flows) then
1160  write (this%iout, '(4x,a)') 'PACKAGE FLOWS WILL BE PRINTED '// &
1161  'FOR ALL MODEL PACKAGES'
1162  end if
1163 
1164  if (found%save_flows) then
1165  write (this%iout, '(4x,a)') &
1166  'FLOWS WILL BE SAVED TO BUDGET FILE SPECIFIED IN OUTPUT CONTROL'
1167  end if
1168 
1169  write (this%iout, '(1x,a)') 'END NAMEFILE OPTIONS:'

◆ package_create()

subroutine prtmodule::package_create ( class(prtmodeltype this,
character(len=*), intent(in)  filtyp,
integer(i4b), intent(in)  ipakid,
integer(i4b), intent(in)  ipaknum,
character(len=*), intent(in)  pakname,
character(len=*), intent(in)  mempath,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout 
)

Definition at line 842 of file prt.f90.

844  ! modules
845  use constantsmodule, only: linelength
846  use prtprpmodule, only: prp_create
847  use apimodule, only: api_create
848  ! dummy
849  class(PrtModelType) :: this
850  character(len=*), intent(in) :: filtyp
851  character(len=LINELENGTH) :: errmsg
852  integer(I4B), intent(in) :: ipakid
853  integer(I4B), intent(in) :: ipaknum
854  character(len=*), intent(in) :: pakname
855  character(len=*), intent(in) :: mempath
856  integer(I4B), intent(in) :: inunit
857  integer(I4B), intent(in) :: iout
858  ! local
859  class(BndType), pointer :: packobj
860  class(BndType), pointer :: packobj2
861  integer(I4B) :: ip
862 
863  ! This part creates the package object
864  select case (filtyp)
865  case ('PRP6')
866  call prp_create(packobj, ipakid, ipaknum, inunit, iout, &
867  this%name, pakname, mempath, this%fmi)
868  case ('API6')
869  call api_create(packobj, ipakid, ipaknum, inunit, iout, &
870  this%name, pakname, mempath)
871  case default
872  write (errmsg, *) 'Invalid package type: ', filtyp
873  call store_error(errmsg, terminate=.true.)
874  end select
875 
876  ! Packages is the bndlist that is associated with the parent model
877  ! The following statement puts a pointer to this package in the ipakid
878  ! position of packages.
879  do ip = 1, this%bndlist%Count()
880  packobj2 => getbndfromlist(this%bndlist, ip)
881  if (packobj2%packName == pakname) then
882  write (errmsg, '(a,a)') 'Cannot create package. Package name '// &
883  'already exists: ', trim(pakname)
884  call store_error(errmsg, terminate=.true.)
885  end if
886  end do
887  call addbndtolist(this%bndlist, packobj)
This module contains the API package methods.
Definition: gwf-api.f90:12
subroutine, public api_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath)
@ brief Create a new package object
Definition: gwf-api.f90:51
subroutine, public prp_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, input_mempath, fmi)
Create a new particle release point package.
Definition: prt-prp.f90:103
Here is the call graph for this function:

◆ prt_ad()

subroutine prtmodule::prt_ad ( class(prtmodeltype this)

Definition at line 347 of file prt.f90.

348  ! modules
350  ! dummy
351  class(PrtModelType) :: this
352  class(BndType), pointer :: packobj
353  ! local
354  integer(I4B) :: irestore
355  integer(I4B) :: ip, n, i
356 
357  ! Reset state variable
358  irestore = 0
359  if (ifailedstepretry > 0) irestore = 1
360 
361  ! Copy masssto into massstoold
362  do n = 1, this%dis%nodes
363  this%massstoold(n) = this%masssto(n)
364  end do
365 
366  ! Advance fmi
367  call this%fmi%fmi_ad()
368 
369  ! Advance
370  do ip = 1, this%bndlist%Count()
371  packobj => getbndfromlist(this%bndlist, ip)
372  call packobj%bnd_ad()
373  if (isimcheck > 0) then
374  call packobj%bnd_ck()
375  end if
376  end do
377  !
378  ! Initialize the flowja array. Flowja is calculated each time,
379  ! even if output is suppressed. (Flowja represents flow of particle
380  ! mass and is positive into a cell. Currently, each particle is assigned
381  ! unit mass.) Flowja is updated continually as particles are tracked
382  ! over the time step and at the end of the time step. The diagonal
383  ! position of the flowja array will contain the flow residual.
384  do i = 1, this%dis%nja
385  this%flowja(i) = dzero
386  end do
integer(i4b) isimcheck
simulation input check flag (1) to check input, (0) to ignore checks
integer(i4b) ifailedstepretry
current retry for this time step
Here is the call graph for this function:

◆ prt_ar()

subroutine prtmodule::prt_ar ( class(prtmodeltype this)

(1) allocates and reads packages part of this model, (2) allocates memory for arrays part of this model object

Definition at line 232 of file prt.f90.

233  ! modules
234  use constantsmodule, only: dhnoflo
235  use prtprpmodule, only: prtprptype
236  use prtmipmodule, only: prtmiptype
238  ! dummy
239  class(PrtModelType) :: this
240  ! locals
241  integer(I4B) :: ip, nprp
242  class(BndType), pointer :: packobj
243 
244  ! Set up basic packages
245  call this%fmi%fmi_ar(this%ibound)
246  if (this%inmip > 0) call this%mip%mip_ar()
247 
248  ! Set up output control and budget
249  call this%oc%oc_ar(this%dis, dhnoflo)
250  call this%budget%set_ibudcsv(this%oc%ibudcsv)
251 
252  ! Select tracking events
253  call this%tracks%select_events( &
254  this%oc%trackrelease, &
255  this%oc%trackfeatexit, &
256  this%oc%tracktimestep, &
257  this%oc%trackterminate, &
258  this%oc%trackweaksink, &
259  this%oc%trackusertime, &
260  this%oc%tracksubfexit, &
261  this%oc%trackdropped)
262 
263  ! Set up boundary pkgs and pkg-scoped track files
264  nprp = 0
265  do ip = 1, this%bndlist%Count()
266  packobj => getbndfromlist(this%bndlist, ip)
267  select type (packobj)
268  type is (prtprptype)
269  nprp = nprp + 1
270  call packobj%prp_set_pointers(this%ibound, this%mip%izone)
271  call packobj%bnd_ar()
272  call packobj%bnd_ar()
273  if (packobj%itrkout > 0) then
274  call this%tracks%init_file( &
275  packobj%itrkout, &
276  iprp=nprp)
277  end if
278  if (packobj%itrkcsv > 0) then
279  call this%tracks%init_file( &
280  packobj%itrkcsv, &
281  csv=.true., &
282  iprp=nprp)
283  end if
284  class default
285  call packobj%bnd_ar()
286  end select
287  end do
288 
289  ! Set up model-scoped track files
290  if (this%oc%itrkout > 0) &
291  call this%tracks%init_file(this%oc%itrkout)
292  if (this%oc%itrkcsv > 0) &
293  call this%tracks%init_file(this%oc%itrkcsv, csv=.true.)
294 
295  ! Set up the tracking method
296  select type (dis => this%dis)
297  type is (distype)
298  call method_dis%init( &
299  fmi=this%fmi, &
300  events=this%events, &
301  izone=this%mip%izone, &
302  flowja=this%flowja, &
303  porosity=this%mip%porosity, &
304  retfactor=this%mip%retfactor, &
305  tracktimes=this%oc%tracktimes)
306  this%method => method_dis
307  type is (disvtype)
308  call method_disv%init( &
309  fmi=this%fmi, &
310  events=this%events, &
311  izone=this%mip%izone, &
312  flowja=this%flowja, &
313  porosity=this%mip%porosity, &
314  retfactor=this%mip%retfactor, &
315  tracktimes=this%oc%tracktimes)
316  this%method => method_disv
317  end select
318 
319  ! Subscribe track output manager to events
320  call this%events%subscribe(this%tracks)
321 
322  ! Set verbose tracing if requested
323  if (this%oc%dump_event_trace) this%tracks%iout = 0
real(dp), parameter dhnoflo
real no flow constant
Definition: Constants.f90:93
type(methoddisvtype), pointer, public method_disv
Definition: MethodPool.f90:12
type(methoddistype), pointer, public method_dis
Definition: MethodPool.f90:11
Particle release point (PRP) package.
Definition: prt-prp.f90:38
Here is the call graph for this function:

◆ prt_bd()

subroutine prtmodule::prt_bd ( class(prtmodeltype this,
integer(i4b), intent(in)  icnvg,
integer(i4b), intent(in)  isuppress_output 
)

(1) Calculate intercell flows (flowja) (2) Calculate package contributions to model budget

Definition at line 489 of file prt.f90.

490  ! modules
491  use tdismodule, only: delt
492  use budgetmodule, only: rate_accumulator
493  ! dummy
494  class(PrtModelType) :: this
495  integer(I4B), intent(in) :: icnvg
496  integer(I4B), intent(in) :: isuppress_output
497  ! local
498  integer(I4B) :: ip
499  class(BndType), pointer :: packobj
500  real(DP) :: rin
501  real(DP) :: rout
502 
503  ! Budget routines (start by resetting). Sole purpose of this section
504  ! is to add in and outs to model budget. All ins and out for a model
505  ! should be added here to this%budget. In a subsequent exchange call,
506  ! exchange flows might also be added.
507  call this%budget%reset()
508  call rate_accumulator(this%ratesto, rin, rout)
509  call this%budget%addentry(rin, rout, delt, budtxt(1), &
510  isuppress_output, ' PRT')
511  do ip = 1, this%bndlist%Count()
512  packobj => getbndfromlist(this%bndlist, ip)
513  call packobj%bnd_bd(this%budget)
514  end do
subroutine, public rate_accumulator(flow, rin, rout)
@ brief Rate accumulator subroutine
Definition: Budget.f90:632
real(dp), pointer, public delt
length of the current time step
Definition: tdis.f90:29
Here is the call graph for this function:

◆ prt_cq()

subroutine prtmodule::prt_cq ( class(prtmodeltype this,
integer(i4b), intent(in)  icnvg,
integer(i4b), intent(in)  isuppress_output 
)

Definition at line 390 of file prt.f90.

391  ! modules
392  use sparsemodule, only: csr_diagsum
393  use tdismodule, only: delt
394  use prtprpmodule, only: prtprptype
395  ! dummy
396  class(PrtModelType) :: this
397  integer(I4B), intent(in) :: icnvg
398  integer(I4B), intent(in) :: isuppress_output
399  ! local
400  integer(I4B) :: i
401  integer(I4B) :: ip
402  class(BndType), pointer :: packobj
403  real(DP) :: tled
404 
405  ! Flowja is calculated each time, even if output is suppressed.
406  ! Flowja represents flow of particle mass and is positive into a cell.
407  ! Currently, each particle is assigned unit mass.
408  !
409  ! Reciprocal of time step size.
410  tled = done / delt
411  !
412  ! Flowja was updated continually as particles were tracked over the
413  ! time step. At this point, flowja contains the net particle mass
414  ! exchanged between cells during the time step. To convert these to
415  ! flow rates (particle mass per time), divide by the time step size.
416  do i = 1, this%dis%nja
417  this%flowja(i) = this%flowja(i) * tled
418  end do
419 
420  ! Particle mass storage
421  call this%prt_cq_sto()
422 
423  ! Go through packages and call cq routines. Just a formality.
424  do ip = 1, this%bndlist%Count()
425  packobj => getbndfromlist(this%bndlist, ip)
426  call packobj%bnd_cq(this%masssto, this%flowja)
427  end do
428 
429  ! Finalize calculation of flowja by adding face flows to the diagonal.
430  ! This results in the flow residual being stored in the diagonal
431  ! position for each cell.
432  call csr_diagsum(this%dis%con%ia, this%flowja)
subroutine csr_diagsum(ia, flowja)
Definition: Sparse.f90:263
Here is the call graph for this function:

◆ prt_cq_sto()

subroutine prtmodule::prt_cq_sto ( class(prtmodeltype this)

Definition at line 436 of file prt.f90.

437  ! modules
438  use tdismodule, only: delt
439  use prtprpmodule, only: prtprptype
440  ! dummy
441  class(PrtModelType) :: this
442  ! local
443  integer(I4B) :: ip
444  class(BndType), pointer :: packobj
445  integer(I4B) :: n
446  integer(I4B) :: np
447  integer(I4B) :: idiag
448  integer(I4B) :: istatus
449  real(DP) :: tled
450  real(DP) :: rate
451 
452  ! Reciprocal of time step size.
453  tled = done / delt
454 
455  ! Particle mass storage rate
456  do n = 1, this%dis%nodes
457  this%masssto(n) = dzero
458  this%ratesto(n) = dzero
459  end do
460  do ip = 1, this%bndlist%Count()
461  packobj => getbndfromlist(this%bndlist, ip)
462  select type (packobj)
463  type is (prtprptype)
464  do np = 1, packobj%nparticles
465  istatus = packobj%particles%istatus(np)
466  ! this may need to change if istatus flags change
467  if ((istatus > 0) .and. (istatus /= term_unreleased)) then
468  n = packobj%particles%itrdomain(np, level_feature)
469  ! Each particle currently assigned unit mass
470  this%masssto(n) = this%masssto(n) + done
471  end if
472  end do
473  end select
474  end do
475  do n = 1, this%dis%nodes
476  rate = -(this%masssto(n) - this%massstoold(n)) * tled
477  this%ratesto(n) = rate
478  idiag = this%dis%con%ia(n)
479  this%flowja(idiag) = this%flowja(idiag) + rate
480  end do
Here is the call graph for this function:

◆ prt_cr()

subroutine, public prtmodule::prt_cr ( character(len=*), intent(in)  filename,
integer(i4b), intent(in)  id,
character(len=*), intent(in)  modelname 
)

Definition at line 122 of file prt.f90.

123  ! modules
124  use listsmodule, only: basemodellist
127  use compilerversion
132  ! dummy
133  character(len=*), intent(in) :: filename
134  integer(I4B), intent(in) :: id
135  character(len=*), intent(in) :: modelname
136  ! local
137  type(PrtModelType), pointer :: this
138  class(BaseModelType), pointer :: model
139  character(len=LENMEMPATH) :: input_mempath
140  character(len=LINELENGTH) :: lst_fname
141  type(PrtNamParamFoundType) :: found
142 
143  ! Allocate a new PRT Model (this)
144  allocate (this)
145 
146  ! Set this before any allocs in the memory manager can be done
147  this%memoryPath = create_mem_path(modelname)
148 
149  ! Allocate event system and track output manager
150  allocate (this%events)
151  allocate (this%tracks)
152 
153  ! Allocate scalars and add model to basemodellist
154  call this%allocate_scalars(modelname)
155  model => this
156  call addbasemodeltolist(basemodellist, model)
157 
158  ! Assign variables
159  this%filename = filename
160  this%name = modelname
161  this%macronym = 'PRT'
162  this%id = id
163 
164  ! Set input model namfile memory path
165  input_mempath = create_mem_path(modelname, 'NAM', idm_context)
166 
167  ! Copy options from input context
168  call mem_set_value(this%iprpak, 'PRINT_INPUT', input_mempath, &
169  found%print_input)
170  call mem_set_value(this%iprflow, 'PRINT_FLOWS', input_mempath, &
171  found%print_flows)
172  call mem_set_value(this%ipakcb, 'SAVE_FLOWS', input_mempath, &
173  found%save_flows)
174 
175  ! Create the list file
176  call this%create_lstfile(lst_fname, filename, found%list, &
177  'PARTICLE TRACKING MODEL (PRT)')
178 
179  ! Activate save_flows if found
180  if (found%save_flows) then
181  this%ipakcb = -1
182  end if
183 
184  ! Create model packages
185  call this%create_packages()
186 
187  ! Log options
188  if (this%iout > 0) then
189  call this%log_namfile_options(found)
190  end if
191 
subroutine, public addbasemodeltolist(list, model)
Definition: BaseModel.f90:161
type(listtype), public basemodellist
Definition: mf6lists.f90:16
Here is the call graph for this function:
Here is the caller graph for this function:

◆ prt_da()

subroutine prtmodule::prt_da ( class(prtmodeltype this)

Definition at line 717 of file prt.f90.

718  ! modules
725  ! dummy
726  class(PrtModelType) :: this
727  ! local
728  integer(I4B) :: ip
729  class(BndType), pointer :: packobj
730 
731  ! Deallocate idm memory
732  call memorystore_remove(this%name, 'NAM', idm_context)
733  call memorystore_remove(component=this%name, context=idm_context)
734 
735  ! Internal packages
736  call this%dis%dis_da()
737  call this%fmi%fmi_da()
738  call this%mip%mip_da()
739  call this%budget%budget_da()
740  call this%oc%oc_da()
741  deallocate (this%dis)
742  deallocate (this%fmi)
743  deallocate (this%mip)
744  deallocate (this%budget)
745  deallocate (this%oc)
746 
747  ! Method objects
750  call destroy_method_pool()
751 
752  ! Boundary packages
753  do ip = 1, this%bndlist%Count()
754  packobj => getbndfromlist(this%bndlist, ip)
755  call packobj%bnd_da()
756  deallocate (packobj)
757  end do
758 
759  ! Scalars
760  call mem_deallocate(this%infmi)
761  call mem_deallocate(this%inmip)
762  call mem_deallocate(this%inadv)
763  call mem_deallocate(this%indsp)
764  call mem_deallocate(this%inssm)
765  call mem_deallocate(this%inmst)
766  call mem_deallocate(this%inmvt)
767  call mem_deallocate(this%inoc)
768 
769  ! Arrays
770  call mem_deallocate(this%masssto)
771  call mem_deallocate(this%massstoold)
772  call mem_deallocate(this%ratesto)
773 
774  call this%tracks%destroy()
775  deallocate (this%events)
776  deallocate (this%tracks)
777 
778  call this%NumericalModelType%model_da()
subroutine, public memorystore_remove(component, subcomponent, context)
subroutine, public destroy_method_cell_pool()
Destroy the cell method pool.
subroutine, public destroy_method_pool()
Destroy the method pool.
Definition: MethodPool.f90:24
subroutine, public destroy_method_subcell_pool()
Destroy the subcell method pool.
Here is the call graph for this function:

◆ prt_df()

subroutine prtmodule::prt_df ( class(prtmodeltype this)

(1) call df routines for each package (2) set variables and pointers

Definition at line 199 of file prt.f90.

200  ! modules
201  use prtprpmodule, only: prtprptype
202  ! dummy
203  class(PrtModelType) :: this
204  ! local
205  integer(I4B) :: ip
206  class(BndType), pointer :: packobj
207 
208  ! Define packages and utility objects
209  call this%dis%dis_df()
210  call this%fmi%fmi_df(this%dis, 1)
211  call this%oc%oc_df()
212  call this%budget%budget_df(niunit_prt, 'MASS', 'M')
213 
214  ! Define packages and assign iout for time series managers
215  do ip = 1, this%bndlist%Count()
216  packobj => getbndfromlist(this%bndlist, ip)
217  call packobj%bnd_df(this%dis%nodes, this%dis)
218  packobj%TsManager%iout = this%iout
219  packobj%TasManager%iout = this%iout
220  end do
221 
222  ! Allocate model arrays
223  call this%allocate_arrays()
224 
Here is the call graph for this function:

◆ prt_ot()

subroutine prtmodule::prt_ot ( class(prtmodeltype this)

Definition at line 518 of file prt.f90.

519  use tdismodule, only: tdis_ot, endofperiod
520  ! dummy
521  class(PrtModelType) :: this
522  ! local
523  integer(I4B) :: idvsave
524  integer(I4B) :: idvprint
525  integer(I4B) :: icbcfl
526  integer(I4B) :: icbcun
527  integer(I4B) :: ibudfl
528  integer(I4B) :: ipflag
529 
530  ! Note: particle tracking output is handled elsewhere
531 
532  ! Set write and print flags
533  idvsave = 0
534  idvprint = 0
535  icbcfl = 0
536  ibudfl = 0
537  if (this%oc%oc_save('CONCENTRATION')) idvsave = 1
538  if (this%oc%oc_print('CONCENTRATION')) idvprint = 1
539  if (this%oc%oc_save('BUDGET')) icbcfl = 1
540  if (this%oc%oc_print('BUDGET')) ibudfl = 1
541  icbcun = this%oc%oc_save_unit('BUDGET')
542 
543  ! Override ibudfl and idvprint flags for nonconvergence
544  ! and end of period
545  ibudfl = this%oc%set_print_flag('BUDGET', 1, endofperiod)
546  idvprint = this%oc%set_print_flag('CONCENTRATION', 1, endofperiod)
547 
548  ! Save and print flows
549  call this%prt_ot_flow(icbcfl, ibudfl, icbcun)
550 
551  ! Save and print dependent variables
552  call this%prt_ot_dv(idvsave, idvprint, ipflag)
553 
554  ! Print budget summaries
555  call this%prt_ot_bdsummary(ibudfl, ipflag)
556 
557  ! Timing Output; if any dependent variables or budgets
558  ! are printed, then ipflag is set to 1.
559  if (ipflag == 1) call tdis_ot(this%iout)
logical(lgp), pointer, public endofperiod
flag indicating end of stress period
Definition: tdis.f90:27
subroutine, public tdis_ot(iout)
Print simulation time.
Definition: tdis.f90:274
Here is the call graph for this function:

◆ prt_ot_bdsummary()

subroutine prtmodule::prt_ot_bdsummary ( class(prtmodeltype this,
integer(i4b), intent(in)  ibudfl,
integer(i4b), intent(inout)  ipflag 
)
private

Definition at line 687 of file prt.f90.

688  ! modules
689  use tdismodule, only: kstp, kper, totim, delt
690  ! dummy
691  class(PrtModelType) :: this
692  integer(I4B), intent(in) :: ibudfl
693  integer(I4B), intent(inout) :: ipflag
694  ! local
695  class(BndType), pointer :: packobj
696  integer(I4B) :: ip
697 
698  ! Package budget summary
699  do ip = 1, this%bndlist%Count()
700  packobj => getbndfromlist(this%bndlist, ip)
701  call packobj%bnd_ot_bdsummary(kstp, kper, this%iout, ibudfl)
702  end do
703 
704  ! model budget summary
705  call this%budget%finalize_step(delt)
706  if (ibudfl /= 0) then
707  ipflag = 1
708  ! model budget summary
709  call this%budget%budget_ot(kstp, kper, this%iout)
710  end if
711 
712  ! Write to budget csv
713  call this%budget%writecsv(totim)
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
Here is the call graph for this function:

◆ prt_ot_dv()

subroutine prtmodule::prt_ot_dv ( class(prtmodeltype this,
integer(i4b), intent(in)  idvsave,
integer(i4b), intent(in)  idvprint,
integer(i4b), intent(inout)  ipflag 
)

Definition at line 666 of file prt.f90.

667  ! dummy
668  class(PrtModelType) :: this
669  integer(I4B), intent(in) :: idvsave
670  integer(I4B), intent(in) :: idvprint
671  integer(I4B), intent(inout) :: ipflag
672  ! local
673  class(BndType), pointer :: packobj
674  integer(I4B) :: ip
675 
676  ! Print advanced package dependent variables
677  do ip = 1, this%bndlist%Count()
678  packobj => getbndfromlist(this%bndlist, ip)
679  call packobj%bnd_ot_dv(idvsave, idvprint)
680  end do
681 
682  ! save head and print head
683  call this%oc%oc_ot(ipflag)
Here is the call graph for this function:

◆ prt_ot_flow()

subroutine prtmodule::prt_ot_flow ( class(prtmodeltype this,
integer(i4b), intent(in)  icbcfl,
integer(i4b), intent(in)  ibudfl,
integer(i4b), intent(in)  icbcun 
)

Definition at line 563 of file prt.f90.

564  use prtprpmodule, only: prtprptype
565  class(PrtModelType) :: this
566  integer(I4B), intent(in) :: icbcfl
567  integer(I4B), intent(in) :: ibudfl
568  integer(I4B), intent(in) :: icbcun
569  class(BndType), pointer :: packobj
570  integer(I4B) :: ip
571 
572  ! Save PRT flows
573  call this%prt_ot_saveflow(this%dis%nja, this%flowja, icbcfl, icbcun)
574  do ip = 1, this%bndlist%Count()
575  packobj => getbndfromlist(this%bndlist, ip)
576  call packobj%bnd_ot_model_flows(icbcfl=icbcfl, ibudfl=0, icbcun=icbcun)
577  end do
578 
579  ! Save advanced package flows
580  do ip = 1, this%bndlist%Count()
581  packobj => getbndfromlist(this%bndlist, ip)
582  call packobj%bnd_ot_package_flows(icbcfl=icbcfl, ibudfl=0)
583  end do
584 
585  ! Print PRT flows
586  call this%prt_ot_printflow(ibudfl, this%flowja)
587  do ip = 1, this%bndlist%Count()
588  packobj => getbndfromlist(this%bndlist, ip)
589  call packobj%bnd_ot_model_flows(icbcfl=icbcfl, ibudfl=ibudfl, icbcun=0)
590  end do
591 
592  ! Print advanced package flows
593  do ip = 1, this%bndlist%Count()
594  packobj => getbndfromlist(this%bndlist, ip)
595  call packobj%bnd_ot_package_flows(icbcfl=0, ibudfl=ibudfl)
596  end do
Here is the call graph for this function:

◆ prt_ot_printflow()

subroutine prtmodule::prt_ot_printflow ( class(prtmodeltype this,
integer(i4b), intent(in)  ibudfl,
real(dp), dimension(:), intent(inout)  flowja 
)
private

Definition at line 627 of file prt.f90.

628  ! modules
629  use tdismodule, only: kper, kstp
630  use constantsmodule, only: lenbigline
631  ! dummy
632  class(PrtModelType) :: this
633  integer(I4B), intent(in) :: ibudfl
634  real(DP), intent(inout), dimension(:) :: flowja
635  ! local
636  character(len=LENBIGLINE) :: line
637  character(len=30) :: tempstr
638  integer(I4B) :: n, ipos, m
639  real(DP) :: qnm
640  ! formats
641  character(len=*), parameter :: fmtiprflow = &
642  "(/,4x,'CALCULATED INTERCELL FLOW &
643  &FOR PERIOD ', i0, ' STEP ', i0)"
644 
645  ! Write flowja to list file if requested
646  if (ibudfl /= 0 .and. this%iprflow > 0) then
647  write (this%iout, fmtiprflow) kper, kstp
648  do n = 1, this%dis%nodes
649  line = ''
650  call this%dis%noder_to_string(n, tempstr)
651  line = trim(tempstr)//':'
652  do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
653  m = this%dis%con%ja(ipos)
654  call this%dis%noder_to_string(m, tempstr)
655  line = trim(line)//' '//trim(tempstr)
656  qnm = flowja(ipos)
657  write (tempstr, '(1pg15.6)') qnm
658  line = trim(line)//' '//trim(adjustl(tempstr))
659  end do
660  write (this%iout, '(a)') trim(line)
661  end do
662  end if
integer(i4b), parameter lenbigline
maximum length of a big line
Definition: Constants.f90:15

◆ prt_ot_saveflow()

subroutine prtmodule::prt_ot_saveflow ( class(prtmodeltype this,
integer(i4b), intent(in)  nja,
real(dp), dimension(nja), intent(in)  flowja,
integer(i4b), intent(in)  icbcfl,
integer(i4b), intent(in)  icbcun 
)

Definition at line 600 of file prt.f90.

601  ! dummy
602  class(PrtModelType) :: this
603  integer(I4B), intent(in) :: nja
604  real(DP), dimension(nja), intent(in) :: flowja
605  integer(I4B), intent(in) :: icbcfl
606  integer(I4B), intent(in) :: icbcun
607  ! local
608  integer(I4B) :: ibinun
609 
610  ! Set unit number for binary output
611  if (this%ipakcb < 0) then
612  ibinun = icbcun
613  elseif (this%ipakcb == 0) then
614  ibinun = 0
615  else
616  ibinun = this%ipakcb
617  end if
618  if (icbcfl == 0) ibinun = 0
619 
620  ! Write the face flows if requested
621  if (ibinun /= 0) then
622  call this%dis%record_connection_array(flowja, ibinun, this%iout)
623  end if

◆ prt_rp()

subroutine prtmodule::prt_rp ( class(prtmodeltype this)

Definition at line 327 of file prt.f90.

328  use tdismodule, only: readnewdata
329  ! dummy
330  class(PrtModelType) :: this
331  ! local
332  class(BndType), pointer :: packobj
333  integer(I4B) :: ip
334 
335  ! Check with TDIS on whether or not it is time to RP
336  if (.not. readnewdata) return
337 
338  ! Read and prepare
339  if (this%inoc > 0) call this%oc%oc_rp()
340  do ip = 1, this%bndlist%Count()
341  packobj => getbndfromlist(this%bndlist, ip)
342  call packobj%bnd_rp()
343  end do
logical(lgp), pointer, public readnewdata
flag indicating time to read new data
Definition: tdis.f90:26
Here is the call graph for this function:

◆ prt_solve()

subroutine prtmodule::prt_solve ( class(prtmodeltype this)
private

Definition at line 918 of file prt.f90.

920  use prtprpmodule, only: prtprptype
923  ! dummy
924  class(PrtModelType) :: this
925  ! local
926  integer(I4B) :: np, ip
927  class(BndType), pointer :: packobj
928  type(ParticleType), pointer :: particle
929  real(DP) :: tmax
930  integer(I4B) :: iprp
931 
932  ! A single particle is reused in the tracking loops
933  ! to avoid allocating and deallocating it each time.
934  ! get() and put() retrieve and store particle state.
935  call create_particle(particle)
936  ! Loop over PRP packages and particles within them.
937  iprp = 0
938  do ip = 1, this%bndlist%Count()
939  packobj => getbndfromlist(this%bndlist, ip)
940  select type (packobj)
941  type is (prtprptype)
942  iprp = iprp + 1
943  do np = 1, packobj%nparticles
944  ! Get the particle from the store
945  call packobj%particles%get(particle, this%id, iprp, np)
946  ! If particle is permanently unreleased, cycle.
947  ! Raise a termination event if we haven't yet.
948  ! TODO: when we have generic dynamic vectors,
949  ! consider terminating permanently unreleased
950  ! in PRP instead of here. For now, status -8
951  ! indicates the permanently unreleased event
952  ! is not yet recorded, status 8 it has been.
953  if (particle%istatus == (-1 * term_unreleased)) then
954  call this%method%terminate(particle, status=term_unreleased)
955  call packobj%particles%put(particle, np)
956  end if
957  if (particle%istatus > active) cycle ! Skip terminated particles
958  particle%istatus = active ! Set active status in case of release
959  ! If the particle was released this time step, emit a release event
960  if (particle%trelease >= totimc) call this%method%release(particle)
961  ! Maximum time is the end of the time step or the particle
962  ! stop time, whichever comes first, unless it's the final
963  ! time step and the extend option is on, in which case
964  ! it's just the particle stop time.
965  if (endofsimulation .and. particle%extend) then
966  tmax = particle%tstop
967  else
968  tmax = min(totimc + delt, particle%tstop)
969  end if
970  ! Apply the tracking method until the maximum time.
971  call this%method%apply(particle, tmax)
972  ! If the particle timed out, terminate it.
973  ! "Timed out" means it's still active but
974  ! - it reached its stop time, or
975  ! - the simulation is over and not extended.
976  ! We can't detect timeout within the tracking
977  ! method because the method just receives the
978  ! maximum time with no context on what it is.
979  ! TODO maybe think about changing that?
980  if (particle%istatus <= active .and. &
981  (particle%ttrack == particle%tstop .or. &
982  (endofsimulation .and. .not. particle%extend))) &
983  call this%method%terminate(particle, status=term_timeout)
984  ! Return the particle to the store
985  call packobj%particles%put(particle, np)
986  end do
987  end select
988  end do
989  call particle%destroy()
990  deallocate (particle)
@, public release
particle was released
@, public terminate
particle terminated
@ term_timeout
terminated at stop time or end of simulation
Definition: Particle.f90:37
@ term_unreleased
terminated permanently unreleased
Definition: Particle.f90:35
logical(lgp), pointer, public endofsimulation
flag indicating end of simulation
Definition: tdis.f90:28
real(dp), pointer, public totimc
simulation time at start of time step
Definition: tdis.f90:33
Here is the call graph for this function:

Variable Documentation

◆ budtxt

character(len=lenbudtxt), dimension(nbditems) prtmodule::budtxt
private

Definition at line 39 of file prt.f90.

39  character(len=LENBUDTXT), dimension(NBDITEMS) :: budtxt

◆ nbditems

integer(i4b), parameter prtmodule::nbditems = 1
private

Definition at line 38 of file prt.f90.

38  integer(I4B), parameter :: NBDITEMS = 1

◆ niunit_prt

integer(i4b), parameter prtmodule::niunit_prt = PRT_NBASEPKG + PRT_NMULTIPKG
private

Definition at line 117 of file prt.f90.

117  integer(I4B), parameter :: NIUNIT_PRT = prt_nbasepkg + prt_nmultipkg

◆ prt_basepkg

character(len=lenpackagetype), dimension(prt_nbasepkg), public prtmodule::prt_basepkg

Definition at line 98 of file prt.f90.

98  character(len=LENPACKAGETYPE), dimension(PRT_NBASEPKG) :: PRT_BASEPKG

◆ prt_multipkg

character(len=lenpackagetype), dimension(prt_nmultipkg), public prtmodule::prt_multipkg

Definition at line 112 of file prt.f90.

112  character(len=LENPACKAGETYPE), dimension(PRT_NMULTIPKG) :: PRT_MULTIPKG

◆ prt_nbasepkg

integer(i4b), parameter, public prtmodule::prt_nbasepkg = 50

PRT6 model base package types. Only listed packages are candidates for input and these will be loaded in the order specified.

Definition at line 97 of file prt.f90.

97  integer(I4B), parameter :: PRT_NBASEPKG = 50

◆ prt_nmultipkg

integer(i4b), parameter, public prtmodule::prt_nmultipkg = 50

PRT6 model multi-instance package types. Only listed packages are candidates for input and these will be loaded in the order specified.

Definition at line 111 of file prt.f90.

111  integer(I4B), parameter :: PRT_NMULTIPKG = 50