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

This module contains the TspSsm Module. More...

Data Types

type  tspssmtype
 Derived type for the SSM Package. More...
 

Functions/Subroutines

subroutine, public ssm_cr (ssmobj, name_model, inunit, iout, fmi, eqnsclfac, depvartype)
 @ brief Create a new SSM package More...
 
subroutine ssm_df (this)
 @ brief Define SSM Package More...
 
subroutine ssm_ar (this, dis, ibound, cnew)
 @ brief Allocate and read SSM Package More...
 
subroutine ssm_rp (this)
 @ brief Read and prepare this SSM Package More...
 
subroutine ssm_ad (this)
 @ brief Advance the SSM Package More...
 
subroutine ssm_term (this, ipackage, ientry, rrate, rhsval, hcofval, cssm, qssm)
 @ brief Calculate the SSM mass flow rate and hcof and rhs values More...
 
subroutine get_ssm_conc (this, ipackage, ientry, nbound_flow, conc, lauxmixed)
 @ brief Provide bound concentration (or temperature) and mixed flag More...
 
subroutine ssm_fc (this, matrix_sln, idxglo, rhs)
 @ brief Fill coefficients More...
 
subroutine ssm_cq (this, flowja)
 @ brief Calculate flow More...
 
subroutine ssm_bd (this, isuppress_output, model_budget)
 @ brief Calculate the global SSM budget terms More...
 
subroutine ssm_ot_flow (this, icbcfl, ibudfl, icbcun)
 @ brief Output flows More...
 
subroutine ssm_da (this)
 @ brief Deallocate More...
 
subroutine allocate_scalars (this)
 @ brief Allocate scalars More...
 
subroutine allocate_arrays (this)
 @ brief Allocate arrays More...
 
subroutine read_options (this)
 @ brief Read package options More...
 
subroutine read_data (this)
 @ brief Read package data More...
 
subroutine read_sources_aux (this)
 @ brief Read SOURCES block More...
 
subroutine read_sources_fileinput (this)
 @ brief Read FILEINPUT block More...
 
subroutine set_iauxpak (this, ip, packname)
 @ brief Set iauxpak array value for package ip More...
 
subroutine set_ssmivec (this, ip, packname)
 @ brief Set ssmivec array value for package ip More...
 
subroutine pak_setup_outputtab (this)
 @ brief Setup the output table More...
 

Variables

character(len=lenftype) ftype = 'SSM'
 
character(len=lenpackagename) text = ' SOURCE-SINK MIX'
 

Detailed Description

This module contains the code for handling sources and sinks associated with groundwater flow model stress packages.

todo: need observations for SSM terms

Function/Subroutine Documentation

◆ allocate_arrays()

subroutine tspssmmodule::allocate_arrays ( class(tspssmtype this)

Allocate array variables for this derived type

Parameters
thisTspSsmType object

Definition at line 757 of file tsp-ssm.f90.

758  ! -- modules
760  ! -- dummy
761  class(TspSsmType) :: this !< TspSsmType object
762  ! -- local
763  integer(I4B) :: nflowpack
764  integer(I4B) :: i
765  !
766  ! -- Allocate
767  nflowpack = this%fmi%nflowpack
768  call mem_allocate(this%iauxpak, nflowpack, 'IAUXPAK', this%memoryPath)
769  call mem_allocate(this%isrctype, nflowpack, 'ISRCTYPE', this%memoryPath)
770  !
771  ! -- Initialize
772  do i = 1, nflowpack
773  this%iauxpak(i) = 0
774  this%isrctype(i) = 0
775  end do
776  !
777  ! -- Allocate the ssmivec array
778  allocate (this%ssmivec(nflowpack))
779  !
780  ! -- Return
781  return

◆ allocate_scalars()

subroutine tspssmmodule::allocate_scalars ( class(tspssmtype this)

Allocate scalar variables for this derived type

Parameters
thisTspSsmType object

Definition at line 733 of file tsp-ssm.f90.

734  ! -- modules
736  ! -- dummy
737  class(TspSsmType) :: this !< TspSsmType object
738  ! -- local
739  !
740  ! -- allocate scalars in NumericalPackageType
741  call this%NumericalPackageType%allocate_scalars()
742  !
743  ! -- Allocate
744  call mem_allocate(this%nbound, 'NBOUND', this%memoryPath)
745  !
746  ! -- Initialize
747  this%nbound = 0
748  !
749  ! -- Return
750  return

◆ get_ssm_conc()

subroutine tspssmmodule::get_ssm_conc ( class(tspssmtype this,
integer(i4b), intent(in)  ipackage,
integer(i4b), intent(in)  ientry,
integer(i4b), intent(in)  nbound_flow,
real(dp), intent(out)  conc,
logical(lgp), intent(out)  lauxmixed 
)

SSM concentrations and temperatures can be provided in auxiliary variables or through separate SPC files. If not provided, the default concentration (or temperature) is zero. This single routine provides the SSM bound concentration (or temperature) based on these different approaches. The mixed flag indicates whether or not the boundary as a mixed type.

Parameters
thisTspSsmType
[in]ipackagepackage number
[in]ientrybound number
[in]nbound_flowsize of flow package bound list
[out]concuser-specified concentration for this bound
[out]lauxmixeduser-specified flag for marking this as a mixed boundary

Definition at line 377 of file tsp-ssm.f90.

379  ! -- dummy
380  class(TspSsmType) :: this !< TspSsmType
381  integer(I4B), intent(in) :: ipackage !< package number
382  integer(I4B), intent(in) :: ientry !< bound number
383  integer(I4B), intent(in) :: nbound_flow !< size of flow package bound list
384  real(DP), intent(out) :: conc !< user-specified concentration for this bound
385  logical(LGP), intent(out) :: lauxmixed !< user-specified flag for marking this as a mixed boundary
386  ! -- local
387  integer(I4B) :: isrctype
388  integer(I4B) :: iauxpos
389 
390  conc = dzero
391  lauxmixed = .false.
392  isrctype = this%isrctype(ipackage)
393 
394  select case (isrctype)
395  case (1, 2)
396  iauxpos = this%iauxpak(ipackage)
397  conc = this%fmi%gwfpackages(ipackage)%auxvar(iauxpos, ientry)
398  if (isrctype == 2) lauxmixed = .true.
399  case (3, 4)
400  conc = this%ssmivec(ipackage)%get_value(ientry, nbound_flow)
401  if (isrctype == 4) lauxmixed = .true.
402  end select
403  !
404  ! -- Return
405  return

◆ pak_setup_outputtab()

subroutine tspssmmodule::pak_setup_outputtab ( class(tspssmtype), intent(inout)  this)

Setup the output table by creating the column headers.

Definition at line 1156 of file tsp-ssm.f90.

1157  ! -- dummy
1158  class(TspSsmType), intent(inout) :: this
1159  ! -- local
1160  character(len=LINELENGTH) :: title
1161  character(len=LINELENGTH) :: text
1162  integer(I4B) :: ntabcol
1163  !
1164  ! -- allocate and initialize the output table
1165  if (this%iprflow /= 0) then
1166  !
1167  ! -- dimension table
1168  ntabcol = 6
1169  !if (this%inamedbound > 0) then
1170  ! ntabcol = ntabcol + 1
1171  !end if
1172  !
1173  ! -- initialize the output table object
1174  title = 'SSM PACKAGE ('//trim(this%packName)// &
1175  ') FLOW RATES'
1176  call table_cr(this%outputtab, this%packName, title)
1177  call this%outputtab%table_df(1, ntabcol, this%iout, transient=.true.)
1178  text = 'NUMBER'
1179  call this%outputtab%initialize_column(text, 10, alignment=tabcenter)
1180  text = 'CELLID'
1181  call this%outputtab%initialize_column(text, 20, alignment=tableft)
1182  text = 'BOUND Q'
1183  call this%outputtab%initialize_column(text, 15, alignment=tabcenter)
1184  text = 'SSM CONC'
1185  call this%outputtab%initialize_column(text, 15, alignment=tabcenter)
1186  text = 'RATE'
1187  call this%outputtab%initialize_column(text, 15, alignment=tabcenter)
1188  text = 'PACKAGE NAME'
1189  call this%outputtab%initialize_column(text, 16, alignment=tabcenter)
1190  !if (this%inamedbound > 0) then
1191  ! text = 'NAME'
1192  ! call this%outputtab%initialize_column(text, 20, alignment=TABLEFT)
1193  !end if
1194  end if
1195  !
1196  ! -- Return
1197  return
Here is the call graph for this function:

◆ read_data()

subroutine tspssmmodule::read_data ( class(tspssmtype this)

Read and set the SSM Package data

Parameters
thisTspSsmType object

Definition at line 839 of file tsp-ssm.f90.

840  ! -- dummy
841  class(TspSsmType) :: this !< TspSsmType object
842  !
843  ! -- read and process required SOURCES block
844  call this%read_sources_aux()
845  !
846  ! -- read and process optional FILEINPUT block
847  call this%read_sources_fileinput()
848  !
849  ! -- Return
850  return

◆ read_options()

subroutine tspssmmodule::read_options ( class(tspssmtype this)

Read and set the SSM Package options

Parameters
thisTspSsmType object

Definition at line 788 of file tsp-ssm.f90.

789  ! -- modules
790  ! -- dummy
791  class(TspSsmType) :: this !< TspSsmType object
792  ! -- local
793  character(len=LINELENGTH) :: keyword
794  integer(I4B) :: ierr
795  logical :: isfound, endOfBlock
796  ! -- formats
797  character(len=*), parameter :: fmtiprflow = &
798  "(4x,'SSM FLOW INFORMATION WILL BE PRINTED TO LISTING FILE &
799  &WHENEVER ICBCFL IS NOT ZERO.')"
800  character(len=*), parameter :: fmtisvflow = &
801  "(4x,'CELL-BY-CELL FLOW INFORMATION WILL BE SAVED TO BINARY FILE &
802  &WHENEVER ICBCFL IS NOT ZERO.')"
803  !
804  ! -- get options block
805  call this%parser%GetBlock('OPTIONS', isfound, ierr, blockrequired=.false., &
806  supportopenclose=.true.)
807  !
808  ! -- parse options block if detected
809  if (isfound) then
810  write (this%iout, '(1x,a)') 'PROCESSING SSM OPTIONS'
811  do
812  call this%parser%GetNextLine(endofblock)
813  if (endofblock) exit
814  call this%parser%GetStringCaps(keyword)
815  select case (keyword)
816  case ('PRINT_FLOWS')
817  this%iprflow = 1
818  write (this%iout, fmtiprflow)
819  case ('SAVE_FLOWS')
820  this%ipakcb = -1
821  write (this%iout, fmtisvflow)
822  case default
823  write (errmsg, '(a,a)') 'Unknown SSM option: ', trim(keyword)
824  call store_error(errmsg)
825  call this%parser%StoreErrorUnit()
826  end select
827  end do
828  write (this%iout, '(1x,a)') 'END OF SSM OPTIONS'
829  end if
830  !
831  ! -- Return
832  return
Here is the call graph for this function:

◆ read_sources_aux()

subroutine tspssmmodule::read_sources_aux ( class(tspssmtype this)

Read SOURCES block and look for auxiliary columns in corresponding flow data.

Parameters
thisTspSsmType object

Definition at line 858 of file tsp-ssm.f90.

859  ! -- dummy
860  class(TspSsmType) :: this !< TspSsmType object
861  ! -- local
862  character(len=LINELENGTH) :: keyword
863  character(len=20) :: srctype
864  integer(I4B) :: ierr
865  integer(I4B) :: ip
866  integer(I4B) :: nflowpack
867  integer(I4B) :: isrctype
868  logical :: isfound, endOfBlock
869  logical :: pakfound
870  logical :: lauxmixed
871  ! -- formats
872  ! -- data
873  !
874  ! -- initialize
875  isfound = .false.
876  lauxmixed = .false.
877  nflowpack = this%fmi%nflowpack
878  !
879  ! -- get sources block
880  call this%parser%GetBlock('SOURCES', isfound, ierr, &
881  supportopenclose=.true., &
882  blockrequired=.true.)
883  if (isfound) then
884  write (this%iout, '(1x,a)') 'PROCESSING SOURCES'
885  do
886  call this%parser%GetNextLine(endofblock)
887  if (endofblock) exit
888  !
889  ! -- read package name and make sure it can be found
890  call this%parser%GetStringCaps(keyword)
891  pakfound = .false.
892  do ip = 1, nflowpack
893  if (trim(adjustl(this%fmi%gwfpackages(ip)%name)) == keyword) then
894  pakfound = .true.
895  exit
896  end if
897  end do
898  if (.not. pakfound) then
899  write (errmsg, '(a,a)') 'Flow package cannot be found: ', &
900  trim(keyword)
901  call store_error(errmsg)
902  call this%parser%StoreErrorUnit()
903  end if
904  !
905  ! -- Ensure package was not specified more than once in SOURCES block
906  if (this%isrctype(ip) /= 0) then
907  write (errmsg, '(a, a)') &
908  'A package cannot be specified more than once in the SSM SOURCES &
909  &block. The following package was specified more than once: ', &
910  trim(keyword)
911  call store_error(errmsg)
912  call this%parser%StoreErrorUnit()
913  end if
914  !
915  ! -- read the source type
916  call this%parser%GetStringCaps(srctype)
917  select case (srctype)
918  case ('AUX')
919  write (this%iout, '(1x,a)') 'AUX SOURCE DETECTED.'
920  isrctype = 1
921  case ('AUXMIXED')
922  write (this%iout, '(1x,a)') 'AUXMIXED SOURCE DETECTED.'
923  lauxmixed = .true.
924  isrctype = 2
925  case default
926  write (errmsg, '(a, a)') &
927  'SRCTYPE must be AUX or AUXMIXED. Found: ', trim(srctype)
928  call store_error(errmsg)
929  call this%parser%StoreErrorUnit()
930  end select
931  !
932  ! -- Store the source type (1 or 2)
933  this%isrctype(ip) = isrctype
934  !
935  ! -- Find and store the auxiliary column
936  call this%set_iauxpak(ip, trim(keyword))
937 
938  end do
939  write (this%iout, '(1x,a)') 'END PROCESSING SOURCES'
940  else
941  write (errmsg, '(a)') 'Required SOURCES block not found.'
942  call store_error(errmsg)
943  call this%parser%StoreErrorUnit()
944  end if
945  !
946  ! -- terminate if errors
947  if (count_errors() > 0) then
948  call this%parser%StoreErrorUnit()
949  end if
950  !
951  ! -- Return
952  return
Here is the call graph for this function:

◆ read_sources_fileinput()

subroutine tspssmmodule::read_sources_fileinput ( class(tspssmtype this)

Read optional FILEINPUT block and initialize an SPC input file reader for each entry.

Parameters
thisTspSsmType object

Definition at line 960 of file tsp-ssm.f90.

961  ! -- dummy
962  class(TspSsmType) :: this !< TspSsmType object
963  ! -- local
964  character(len=LINELENGTH) :: keyword
965  character(len=LINELENGTH) :: keyword2
966  character(len=20) :: srctype
967  integer(I4B) :: ierr
968  integer(I4B) :: ip
969  integer(I4B) :: nflowpack
970  integer(I4B) :: isrctype
971  logical :: isfound, endOfBlock
972  logical :: pakfound
973  logical :: lauxmixed
974  ! -- formats
975  ! -- data
976  !
977  ! -- initialize
978  isfound = .false.
979  lauxmixed = .false.
980  nflowpack = this%fmi%nflowpack
981  !
982  ! -- get sources_file block
983  call this%parser%GetBlock('FILEINPUT', isfound, ierr, &
984  supportopenclose=.true., &
985  blockrequired=.false.)
986  if (isfound) then
987  write (this%iout, '(1x,a)') 'PROCESSING FILEINPUT'
988  do
989  call this%parser%GetNextLine(endofblock)
990  if (endofblock) exit
991  !
992  ! -- read package name and make sure it can be found
993  call this%parser%GetStringCaps(keyword)
994  pakfound = .false.
995  do ip = 1, nflowpack
996  if (trim(adjustl(this%fmi%gwfpackages(ip)%name)) == keyword) then
997  pakfound = .true.
998  exit
999  end if
1000  end do
1001  if (.not. pakfound) then
1002  write (errmsg, '(a,a)') 'Flow package cannot be found: ', &
1003  trim(keyword)
1004  call store_error(errmsg)
1005  call this%parser%StoreErrorUnit()
1006  end if
1007  !
1008  ! -- Ensure package was not specified more than once in SOURCES block
1009  if (this%isrctype(ip) /= 0) then
1010  write (errmsg, '(a, a)') &
1011  'A package cannot be specified more than once in the SSM SOURCES &
1012  &and SOURCES_FILES blocks. The following package was specified &
1013  &more than once: ', &
1014  trim(keyword)
1015  call store_error(errmsg)
1016  call this%parser%StoreErrorUnit()
1017  end if
1018  !
1019  ! -- read the source type
1020  call this%parser%GetStringCaps(srctype)
1021  select case (srctype)
1022  case ('SPC6')
1023  write (this%iout, '(1x,a)') 'SPC6 SOURCE DETECTED.'
1024  isrctype = 3
1025  !
1026  ! verify filein is next
1027  call this%parser%GetStringCaps(keyword2)
1028  if (trim(adjustl(keyword2)) /= 'FILEIN') then
1029  errmsg = 'SPC6 keyword must be followed by "FILEIN" '// &
1030  'then by filename and optionally by <MIXED>.'
1031  call store_error(errmsg)
1032  call this%parser%StoreErrorUnit()
1033  end if
1034  !
1035  ! -- Use set_ssmivec to read file name and set up
1036  ! ssmi file object
1037  call this%set_ssmivec(ip, trim(keyword))
1038  !
1039  ! -- check for optional MIXED keyword and set isrctype to 4 if found
1040  call this%parser%GetStringCaps(keyword2)
1041  if (trim(keyword2) == 'MIXED') then
1042  isrctype = 4
1043  write (this%iout, '(1x,a,a)') 'ASSIGNED MIXED SSM TYPE TO PACKAGE ', &
1044  trim(keyword)
1045  end if
1046  case default
1047  write (errmsg, '(a,a)') &
1048  'SRCTYPE must be SPC6. Found: ', trim(srctype)
1049  call store_error(errmsg)
1050  call this%parser%StoreErrorUnit()
1051  end select
1052  !
1053  ! -- Store the source type (3 or 4)
1054  this%isrctype(ip) = isrctype
1055 
1056  end do
1057  write (this%iout, '(1x,a)') 'END PROCESSING FILEINPUT'
1058  else
1059  write (this%iout, '(1x,a)') &
1060  'OPTIONAL FILEINPUT BLOCK NOT FOUND. CONTINUING.'
1061  end if
1062  !
1063  ! -- terminate if errors
1064  if (count_errors() > 0) then
1065  call this%parser%StoreErrorUnit()
1066  end if
1067  !
1068  ! -- Return
1069  return
Here is the call graph for this function:

◆ set_iauxpak()

subroutine tspssmmodule::set_iauxpak ( class(tspssmtype), intent(inout)  this,
integer(i4b), intent(in)  ip,
character(len=*), intent(in)  packname 
)

The next call to parser will return the auxiliary name for package ip in the SSM SOURCES block. The routine searches through the auxiliary names in package ip and sets iauxpak to the column number corresponding to the correct auxiliary column.

Parameters
[in,out]thisTspSsmType
[in]ippackage number
[in]packnamename of package

Definition at line 1080 of file tsp-ssm.f90.

1081  ! -- dummy
1082  class(TspSsmType), intent(inout) :: this !< TspSsmType
1083  integer(I4B), intent(in) :: ip !< package number
1084  character(len=*), intent(in) :: packname !< name of package
1085  ! -- local
1086  character(len=LENAUXNAME) :: auxname
1087  logical :: auxfound
1088  integer(I4B) :: iaux
1089  !
1090  ! -- read name of auxiliary column
1091  call this%parser%GetStringCaps(auxname)
1092  auxfound = .false.
1093  do iaux = 1, this%fmi%gwfpackages(ip)%naux
1094  if (trim(this%fmi%gwfpackages(ip)%auxname(iaux)) == &
1095  trim(auxname)) then
1096  auxfound = .true.
1097  exit
1098  end if
1099  end do
1100  if (.not. auxfound) then
1101  write (errmsg, '(a, a)') &
1102  'Auxiliary name cannot be found: ', trim(auxname)
1103  call store_error(errmsg)
1104  call this%parser%StoreErrorUnit()
1105  end if
1106  !
1107  ! -- set iauxpak and write message
1108  this%iauxpak(ip) = iaux
1109  write (this%iout, '(4x, a, i0, a, a)') 'USING AUX COLUMN ', &
1110  iaux, ' IN PACKAGE ', trim(packname)
1111  !
1112  ! -- Return
1113  return
Here is the call graph for this function:

◆ set_ssmivec()

subroutine tspssmmodule::set_ssmivec ( class(tspssmtype), intent(inout)  this,
integer(i4b), intent(in)  ip,
character(len=*), intent(in)  packname 
)

The next call to parser will return the input file name for package ip in the SSM SOURCES block. The routine then initializes the SPC input file.

Parameters
[in,out]thisTspSsmType
[in]ippackage number
[in]packnamename of package

Definition at line 1122 of file tsp-ssm.f90.

1123  ! -- module
1124  use inputoutputmodule, only: openfile, getunit
1125  ! -- dummy
1126  class(TspSsmType), intent(inout) :: this !< TspSsmType
1127  integer(I4B), intent(in) :: ip !< package number
1128  character(len=*), intent(in) :: packname !< name of package
1129  ! -- local
1130  character(len=LINELENGTH) :: filename
1131  type(GwtSpcType), pointer :: ssmiptr
1132  integer(I4B) :: inunit
1133  !
1134  ! -- read file name
1135  call this%parser%GetString(filename)
1136  inunit = getunit()
1137  call openfile(inunit, this%iout, filename, 'SPC', filstat_opt='OLD')
1138 
1139  ! -- Create the SPC file object
1140  ssmiptr => this%ssmivec(ip)
1141  call ssmiptr%initialize(this%dis, ip, inunit, this%iout, this%name_model, &
1142  trim(packname))
1143 
1144  write (this%iout, '(4x, a, a, a, a, a)') 'USING SPC INPUT FILE ', &
1145  trim(filename), ' TO SET ', trim(this%depvartype), &
1146  'S FOR PACKAGE ', trim(packname)
1147  !
1148  ! -- Return
1149  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
Here is the call graph for this function:

◆ ssm_ad()

subroutine tspssmmodule::ssm_ad ( class(tspssmtype this)

This routine is called from gwt_ad(). It is called at the beginning of each time step. The total number of flow boundaries is counted and stored in thisnbound. Also, if any SPC input files are used to provide source and sink concentrations and time series are referenced in those files, then ssm concenrations must be interpolated for the time step.

Parameters
thisTspSsmType object

Definition at line 229 of file tsp-ssm.f90.

230  ! -- modules
231  ! -- dummy
232  class(TspSsmType) :: this !< TspSsmType object
233  ! -- local
234  integer(I4B) :: ip
235  type(GwtSpcType), pointer :: ssmiptr
236  integer(I4B) :: i
237  integer(I4B) :: node
238 ! ------------------------------------------------------------------------------
239  !
240  ! -- Calculate total number of existing flow boundaries. It is possible
241  ! that a node may equal zero. In this case, the bound should be
242  ! skipped and not written to ssm output.
243  this%nbound = 0
244  do ip = 1, this%fmi%nflowpack
245  if (this%fmi%iatp(ip) /= 0) cycle
246  do i = 1, this%fmi%gwfpackages(ip)%nbound
247  node = this%fmi%gwfpackages(ip)%nodelist(i)
248  if (node > 0) then
249  this%nbound = this%nbound + 1
250  end if
251  end do
252  end do
253  !
254  ! -- Call the ad method on any ssm input files so that values for
255  ! time-series are interpolated
256  do ip = 1, this%fmi%nflowpack
257  if (this%fmi%iatp(ip) /= 0) cycle
258  if (this%isrctype(ip) == 3 .or. this%isrctype(ip) == 4) then
259  ssmiptr => this%ssmivec(ip)
260  call ssmiptr%spc_ad(this%fmi%gwfpackages(ip)%nbound, &
261  this%fmi%gwfpackages(ip)%budtxt)
262  end if
263  end do
264  !
265  ! -- Return
266  return

◆ ssm_ar()

subroutine tspssmmodule::ssm_ar ( class(tspssmtype this,
class(disbasetype), intent(in), pointer  dis,
integer(i4b), dimension(:), pointer, contiguous  ibound,
real(dp), dimension(:), pointer, contiguous  cnew 
)

This routine is called from gwt_ar(). It allocates arrays, reads options and data, and sets up the output table.

Parameters
thisTspSsmType object
[in]disdiscretization package
iboundGWT model ibound
cnewGWT model dependent variable

Definition at line 142 of file tsp-ssm.f90.

143  ! -- modules
145  ! -- dummy
146  class(TspSsmType) :: this !< TspSsmType object
147  class(DisBaseType), pointer, intent(in) :: dis !< discretization package
148  integer(I4B), dimension(:), pointer, contiguous :: ibound !< GWT model ibound
149  real(DP), dimension(:), pointer, contiguous :: cnew !< GWT model dependent variable
150  ! -- local
151  ! -- formats
152  character(len=*), parameter :: fmtssm = &
153  "(1x,/1x,'SSM -- SOURCE-SINK MIXING PACKAGE, VERSION 1, 8/25/2017', &
154  &' INPUT READ FROM UNIT ', i0, //)"
155  !
156  ! --print a message identifying the storage package.
157  write (this%iout, fmtssm) this%inunit
158  !
159  ! -- store pointers to arguments that were passed in
160  this%dis => dis
161  this%ibound => ibound
162  this%cnew => cnew
163  !
164  ! -- Check to make sure that there are flow packages
165  if (this%fmi%nflowpack == 0) then
166  write (errmsg, '(a)') 'SSM package does not detect any boundary flows &
167  &that require SSM terms. Activate GWF-GWT &
168  &exchange or activate FMI package and provide a &
169  &budget file that contains boundary flows. If no &
170  &boundary flows are present in corresponding GWF &
171  &model then this SSM package should be removed.'
172  call store_error(errmsg)
173  call this%parser%StoreErrorUnit()
174  end if
175  !
176  ! -- Allocate arrays
177  call this%allocate_arrays()
178  !
179  ! -- Read ssm options
180  call this%read_options()
181  !
182  ! -- read the data block
183  call this%read_data()
184  !
185  ! -- setup the output table
186  call this%pak_setup_outputtab()
187  !
188  ! -- Return
189  return
Here is the call graph for this function:

◆ ssm_bd()

subroutine tspssmmodule::ssm_bd ( class(tspssmtype this,
integer(i4b), intent(in)  isuppress_output,
type(budgettype), intent(inout)  model_budget 
)

Calculate the global SSM budget terms using separate in and out entries for each flow package.

Parameters
thisTspSsmType object
[in]isuppress_outputflag to suppress output
[in,out]model_budgetbudget object for the GWT model

Definition at line 498 of file tsp-ssm.f90.

499  ! -- modules
500  use tdismodule, only: delt
501  use budgetmodule, only: budgettype
502  ! -- dummy
503  class(TspSsmType) :: this !< TspSsmType object
504  integer(I4B), intent(in) :: isuppress_output !< flag to suppress output
505  type(BudgetType), intent(inout) :: model_budget !< budget object for the GWT model
506  ! -- local
507  character(len=LENBUDROWLABEL) :: rowlabel
508  integer(I4B) :: ip
509  integer(I4B) :: i
510  integer(I4B) :: n
511  real(DP) :: rate
512  real(DP) :: rin
513  real(DP) :: rout
514  !
515  ! -- do for each flow package, unless it is being handled by an advanced
516  ! transport package
517  do ip = 1, this%fmi%nflowpack
518  !
519  ! -- cycle if package is being managed as an advanced package
520  if (this%fmi%iatp(ip) /= 0) cycle
521  !
522  ! -- Initialize the rate accumulators
523  rin = dzero
524  rout = dzero
525  !
526  ! -- do for each boundary
527  do i = 1, this%fmi%gwfpackages(ip)%nbound
528  n = this%fmi%gwfpackages(ip)%nodelist(i)
529  if (n <= 0) cycle
530  call this%ssm_term(ip, i, rrate=rate)
531  if (rate < dzero) then
532  rout = rout - rate
533  else
534  rin = rin + rate
535  end if
536  !
537  end do
538  !
539  rowlabel = 'SSM_'//adjustl(trim(this%fmi%flowpacknamearray(ip)))
540  call model_budget%addentry(rin, rout, delt, &
541  this%fmi%gwfpackages(ip)%budtxt, &
542  isuppress_output, rowlabel=rowlabel)
543  end do
544  !
545  ! -- Return
546  return
This module contains the BudgetModule.
Definition: Budget.f90:20
real(dp), pointer, public delt
length of the current time step
Definition: tdis.f90:29
Derived type for the Budget object.
Definition: Budget.f90:39

◆ ssm_cq()

subroutine tspssmmodule::ssm_cq ( class(tspssmtype this,
real(dp), dimension(:), intent(inout), contiguous  flowja 
)

Calculate the resulting mass flow between the boundary and the connected GWT model cell. Update the diagonal position of the flowja array so that it ultimately contains the solute balance residual.

Parameters
thisTspSsmType object
[in,out]flowjaflow across each face in the model grid

Definition at line 459 of file tsp-ssm.f90.

460  ! -- modules
461  ! -- dummy
462  class(TspSsmType) :: this !< TspSsmType object
463  real(DP), dimension(:), contiguous, intent(inout) :: flowja !< flow across each face in the model grid
464  ! -- local
465  integer(I4B) :: ip
466  integer(I4B) :: i
467  integer(I4B) :: n
468  integer(I4B) :: idiag
469  real(DP) :: rate
470  !
471  ! -- do for each flow package
472  do ip = 1, this%fmi%nflowpack
473  !
474  ! -- cycle if package is being managed as an advanced package
475  if (this%fmi%iatp(ip) /= 0) cycle
476  !
477  ! -- do for each boundary
478  do i = 1, this%fmi%gwfpackages(ip)%nbound
479  n = this%fmi%gwfpackages(ip)%nodelist(i)
480  if (n <= 0) cycle
481  call this%ssm_term(ip, i, rrate=rate)
482  idiag = this%dis%con%ia(n)
483  flowja(idiag) = flowja(idiag) + rate
484  !
485  end do
486  !
487  end do
488  !
489  ! -- Return
490  return

◆ ssm_cr()

subroutine, public tspssmmodule::ssm_cr ( type(tspssmtype), pointer  ssmobj,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout,
type(tspfmitype), intent(in), target  fmi,
real(dp), intent(in), pointer  eqnsclfac,
character(len=lenvarname), intent(in)  depvartype 
)

Create a new SSM package by defining names, allocating scalars and initializing the parser.

Parameters
ssmobjTspSsmType object
[in]name_modelname of the model
[in]inunitfortran unit for input
[in]ioutfortran unit for output
[in]fmiTransport FMI package
[in]eqnsclfacgoverning equation scale factor
[in]depvartypedependent variable type ('concentration' or 'temperature')

Definition at line 82 of file tsp-ssm.f90.

84  ! -- dummy
85  type(TspSsmType), pointer :: ssmobj !< TspSsmType object
86  character(len=*), intent(in) :: name_model !< name of the model
87  integer(I4B), intent(in) :: inunit !< fortran unit for input
88  integer(I4B), intent(in) :: iout !< fortran unit for output
89  type(TspFmiType), intent(in), target :: fmi !< Transport FMI package
90  real(DP), intent(in), pointer :: eqnsclfac !< governing equation scale factor
91  character(len=LENVARNAME), intent(in) :: depvartype !< dependent variable type ('concentration' or 'temperature')
92  !
93  ! -- Create the object
94  allocate (ssmobj)
95  !
96  ! -- create name and memory path
97  call ssmobj%set_names(1, name_model, 'SSM', 'SSM')
98  !
99  ! -- Allocate scalars
100  call ssmobj%allocate_scalars()
101  !
102  ! -- Set variables
103  ssmobj%inunit = inunit
104  ssmobj%iout = iout
105  ssmobj%fmi => fmi
106  ssmobj%eqnsclfac => eqnsclfac
107  !
108  ! -- Initialize block parser
109  call ssmobj%parser%Initialize(ssmobj%inunit, ssmobj%iout)
110  !
111  ! -- Store pointer to labels associated with the current model so that the
112  ! package has access to the corresponding dependent variable type
113  ssmobj%depvartype = depvartype
114  !
115  ! -- Return
116  return
Here is the caller graph for this function:

◆ ssm_da()

subroutine tspssmmodule::ssm_da ( class(tspssmtype this)

Deallocate the memory associated with this derived type

Parameters
thisTspSsmType object

Definition at line 684 of file tsp-ssm.f90.

685  ! -- modules
687  ! -- dummy
688  class(TspSsmType) :: this !< TspSsmType object
689  ! -- local
690  integer(I4B) :: ip
691  type(GwtSpcType), pointer :: ssmiptr
692  !
693  ! -- Deallocate the ssmi objects if package was active
694  if (this%inunit > 0) then
695  do ip = 1, size(this%ssmivec)
696  if (this%isrctype(ip) == 3 .or. this%isrctype(ip) == 4) then
697  ssmiptr => this%ssmivec(ip)
698  call ssmiptr%spc_da()
699  end if
700  end do
701  deallocate (this%ssmivec)
702  end if
703  !
704  ! -- Deallocate arrays if package was active
705  if (this%inunit > 0) then
706  call mem_deallocate(this%iauxpak)
707  call mem_deallocate(this%isrctype)
708  this%ibound => null()
709  this%fmi => null()
710  end if
711  !
712  ! -- output table object
713  if (associated(this%outputtab)) then
714  call this%outputtab%table_da()
715  deallocate (this%outputtab)
716  nullify (this%outputtab)
717  end if
718  !
719  ! -- Scalars
720  call mem_deallocate(this%nbound)
721  !
722  ! -- deallocate parent
723  call this%NumericalPackageType%da()
724  !
725  ! -- Return
726  return

◆ ssm_df()

subroutine tspssmmodule::ssm_df ( class(tspssmtype this)

This routine is called from gwt_df(), but does not do anything because df is typically used to set up dimensions. For the ssm package, the total number of ssm entries is defined by the flow model.

Parameters
thisTspSsmType object

Definition at line 125 of file tsp-ssm.f90.

126  ! -- modules
128  ! -- dummy
129  class(TspSsmType) :: this !< TspSsmType object
130  ! -- local
131  ! -- formats
132  !
133  ! -- Return
134  return

◆ ssm_fc()

subroutine tspssmmodule::ssm_fc ( class(tspssmtype this,
class(matrixbasetype), pointer  matrix_sln,
integer(i4b), dimension(:), intent(in)  idxglo,
real(dp), dimension(:), intent(inout)  rhs 
)

This routine adds the effects of the SSM to the matrix equations by updating the a matrix and right-hand side vector.

Definition at line 413 of file tsp-ssm.f90.

414  ! -- modules
415  ! -- dummy
416  class(TspSsmType) :: this
417  class(MatrixBaseType), pointer :: matrix_sln
418  integer(I4B), intent(in), dimension(:) :: idxglo
419  real(DP), intent(inout), dimension(:) :: rhs
420  ! -- local
421  integer(I4B) :: ip
422  integer(I4B) :: i
423  integer(I4B) :: n
424  integer(I4B) :: idiag
425  integer(I4B) :: nflowpack
426  integer(I4B) :: nbound
427  real(DP) :: hcofval
428  real(DP) :: rhsval
429  !
430  ! -- do for each flow package
431  nflowpack = this%fmi%nflowpack
432  do ip = 1, nflowpack
433  if (this%fmi%iatp(ip) /= 0) cycle
434  !
435  ! -- do for each entry in package (ip)
436  nbound = this%fmi%gwfpackages(ip)%nbound
437  do i = 1, nbound
438  n = this%fmi%gwfpackages(ip)%nodelist(i)
439  if (n <= 0) cycle
440  call this%ssm_term(ip, i, rhsval=rhsval, hcofval=hcofval)
441  idiag = idxglo(this%dis%con%ia(n))
442  call matrix_sln%add_value_pos(idiag, hcofval)
443  rhs(n) = rhs(n) + rhsval
444  !
445  end do
446  !
447  end do
448  !
449  ! -- Return
450  return

◆ ssm_ot_flow()

subroutine tspssmmodule::ssm_ot_flow ( class(tspssmtype this,
integer(i4b), intent(in)  icbcfl,
integer(i4b), intent(in)  ibudfl,
integer(i4b), intent(in)  icbcun 
)

Based on user-specified controls, print SSM mass flow rates to the GWT listing file and/or write the SSM mass flow rates to the GWT binary budget file.

Parameters
thisTspSsmType object
[in]icbcflflag for writing binary budget terms
[in]ibudflflag for printing budget terms to list file
[in]icbcunfortran unit number for binary budget file

Definition at line 555 of file tsp-ssm.f90.

556  ! -- modules
557  use tdismodule, only: kstp, kper
559  ! -- dummy
560  class(TspSsmType) :: this !< TspSsmType object
561  integer(I4B), intent(in) :: icbcfl !< flag for writing binary budget terms
562  integer(I4B), intent(in) :: ibudfl !< flag for printing budget terms to list file
563  integer(I4B), intent(in) :: icbcun !< fortran unit number for binary budget file
564  ! -- local
565  character(len=LINELENGTH) :: title
566  integer(I4B) :: node, nodeu
567  character(len=20) :: nodestr
568  integer(I4B) :: maxrows
569  integer(I4B) :: ip
570  integer(I4B) :: i, n2, ibinun
571  real(DP) :: rrate
572  real(DP) :: qssm
573  real(DP) :: cssm
574  integer(I4B) :: naux
575  real(DP), dimension(0, 0) :: auxvar
576  character(len=LENAUXNAME), dimension(0) :: auxname
577  ! -- for observations
578  character(len=LENBOUNDNAME) :: bname
579  ! -- formats
580  character(len=*), parameter :: fmttkk = &
581  &"(1X,/1X,A,' PERIOD ',I0,' STEP ',I0)"
582  !
583  ! -- set maxrows
584  maxrows = 0
585  if (ibudfl /= 0 .and. this%iprflow /= 0) then
586  call this%outputtab%set_kstpkper(kstp, kper)
587  do ip = 1, this%fmi%nflowpack
588  if (this%fmi%iatp(ip) /= 0) cycle
589  !
590  ! -- do for each boundary
591  do i = 1, this%fmi%gwfpackages(ip)%nbound
592  node = this%fmi%gwfpackages(ip)%nodelist(i)
593  if (node > 0) then
594  maxrows = maxrows + 1
595  end if
596  end do
597  end do
598  if (maxrows > 0) then
599  call this%outputtab%set_maxbound(maxrows)
600  end if
601  title = 'SSM PACKAGE ('//trim(this%packName)// &
602  ') FLOW RATES'
603  call this%outputtab%set_title(title)
604  end if
605  !
606  ! -- Set unit number for binary output
607  if (this%ipakcb < 0) then
608  ibinun = icbcun
609  else if (this%ipakcb == 0) then
610  ibinun = 0
611  else
612  ibinun = this%ipakcb
613  end if
614  if (icbcfl == 0) ibinun = 0
615  !
616  ! -- If cell-by-cell flows will be saved as a list, write header.
617  if (ibinun /= 0) then
618  naux = 0
619  call this%dis%record_srcdst_list_header(text, this%name_model, &
620  this%name_model, this%name_model, &
621  this%packName, naux, auxname, &
622  ibinun, this%nbound, this%iout)
623  end if
624  !
625  ! -- If no boundaries, skip flow calculations.
626  if (this%nbound > 0) then
627  !
628  ! -- Loop through each boundary calculating flow.
629  do ip = 1, this%fmi%nflowpack
630  if (this%fmi%iatp(ip) /= 0) cycle
631  !
632  ! -- do for each boundary
633  do i = 1, this%fmi%gwfpackages(ip)%nbound
634  !
635  ! -- Calculate rate for this entry
636  node = this%fmi%gwfpackages(ip)%nodelist(i)
637  if (node <= 0) cycle
638  call this%ssm_term(ip, i, rrate=rrate, qssm=qssm, cssm=cssm)
639  !
640  ! -- Print the individual rates if the budget is being printed
641  ! and PRINT_FLOWS was specified (this%iprflow<0)
642  if (ibudfl /= 0) then
643  if (this%iprflow /= 0) then
644  !
645  ! -- set nodestr and write outputtab table
646  nodeu = this%dis%get_nodeuser(node)
647  call this%dis%nodeu_to_string(nodeu, nodestr)
648  bname = this%fmi%gwfpackages(ip)%name
649  call this%outputtab%add_term(i)
650  call this%outputtab%add_term(trim(adjustl(nodestr)))
651  call this%outputtab%add_term(qssm)
652  call this%outputtab%add_term(cssm)
653  call this%outputtab%add_term(rrate)
654  call this%outputtab%add_term(bname)
655  end if
656  end if
657  !
658  ! -- If saving cell-by-cell flows in list, write flow
659  if (ibinun /= 0) then
660  n2 = i
661  call this%dis%record_mf6_list_entry(ibinun, node, n2, rrate, &
662  naux, auxvar(:, i), &
663  olconv2=.false.)
664  end if
665  !
666  end do
667  !
668  end do
669  end if
670  if (ibudfl /= 0) then
671  if (this%iprflow /= 0) then
672  write (this%iout, '(1x)')
673  end if
674  end if
675  !
676  ! -- Return
677  return
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter lenpackagename
maximum length of the package name
Definition: Constants.f90:22
integer(i4b), parameter lenauxname
maximum length of a aux variable
Definition: Constants.f90:34
integer(i4b), parameter lenboundname
maximum length of a bound name
Definition: Constants.f90:35
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:64
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

◆ ssm_rp()

subroutine tspssmmodule::ssm_rp ( class(tspssmtype this)

This routine is called from gwt_rp(). It is called at the beginning of each stress period. If any SPC input files are used to provide source and sink concentrations, then period blocks for the current stress period are read.

Parameters
thisTspSsmType object

Definition at line 199 of file tsp-ssm.f90.

200  ! -- modules
201  ! -- dummy
202  class(TspSsmType) :: this !< TspSsmType object
203  ! -- local
204  integer(I4B) :: ip
205  type(GwtSpcType), pointer :: ssmiptr
206  ! -- formats
207  !
208  ! -- Call the rp method on any ssm input files
209  do ip = 1, this%fmi%nflowpack
210  if (this%fmi%iatp(ip) /= 0) cycle
211  if (this%isrctype(ip) == 3 .or. this%isrctype(ip) == 4) then
212  ssmiptr => this%ssmivec(ip)
213  call ssmiptr%spc_rp()
214  end if
215  end do
216  !
217  ! -- Return
218  return

◆ ssm_term()

subroutine tspssmmodule::ssm_term ( class(tspssmtype this,
integer(i4b), intent(in)  ipackage,
integer(i4b), intent(in)  ientry,
real(dp), intent(out), optional  rrate,
real(dp), intent(out), optional  rhsval,
real(dp), intent(out), optional  hcofval,
real(dp), intent(out), optional  cssm,
real(dp), intent(out), optional  qssm 
)

This is the primary SSM routine that calculates the matrix coefficient and right-hand-side value for any package and package entry. It returns several different optional variables that are used throughout this package to update matrix terms, budget calculations, and output tables.

Parameters
thisTspSsmType
[in]ipackagepackage number
[in]ientrybound number
[out]rratecalculated mass flow rate
[out]rhsvalcalculated rhs value
[out]hcofvalcalculated hcof value
[out]cssmcalculated source concentration depending on flow direction
[out]qssmwater flow rate into model cell from boundary package

Definition at line 276 of file tsp-ssm.f90.

278  ! -- dummy
279  class(TspSsmType) :: this !< TspSsmType
280  integer(I4B), intent(in) :: ipackage !< package number
281  integer(I4B), intent(in) :: ientry !< bound number
282  real(DP), intent(out), optional :: rrate !< calculated mass flow rate
283  real(DP), intent(out), optional :: rhsval !< calculated rhs value
284  real(DP), intent(out), optional :: hcofval !< calculated hcof value
285  real(DP), intent(out), optional :: cssm !< calculated source concentration depending on flow direction
286  real(DP), intent(out), optional :: qssm !< water flow rate into model cell from boundary package
287  ! -- local
288  logical(LGP) :: lauxmixed
289  integer(I4B) :: n
290  integer(I4B) :: nbound_flow
291  real(DP) :: qbnd
292  real(DP) :: ctmp
293  real(DP) :: omega
294  real(DP) :: hcoftmp
295  real(DP) :: rhstmp
296  !
297  ! -- retrieve node number, qbnd and iauxpos
298  hcoftmp = dzero
299  rhstmp = dzero
300  ctmp = dzero
301  qbnd = dzero
302  nbound_flow = this%fmi%gwfpackages(ipackage)%nbound
303  n = this%fmi%gwfpackages(ipackage)%nodelist(ientry)
304  !
305  ! -- If cell is active (ibound > 0) then calculate values
306  if (this%ibound(n) > 0) then
307  !
308  ! -- retrieve qbnd and iauxpos
309  qbnd = this%fmi%gwfpackages(ipackage)%get_flow(ientry)
310  call this%get_ssm_conc(ipackage, ientry, nbound_flow, ctmp, lauxmixed)
311  !
312  ! -- assign values for hcoftmp, rhstmp, and ctmp for subsequent assignment
313  ! of hcof, rhs, and rate
314  if (.not. lauxmixed) then
315  !
316  ! -- If qbnd is positive, then concentration represents the inflow
317  ! concentration. If qbnd is negative, then the outflow concentration
318  ! is set equal to the simulated cell concentration
319  if (qbnd >= dzero) then
320  omega = dzero ! rhs
321  else
322  ctmp = this%cnew(n)
323  omega = done ! lhs
324  if (ctmp < dzero) then
325  omega = dzero ! concentration is negative, so set mass flux to zero
326  end if
327  end if
328  else
329  !
330  ! -- lauxmixed value indicates that this is a mixed sink type where
331  ! the concentration value represents the injected concentration if
332  ! qbnd is positive. If qbnd is negative, then the withdrawn water
333  ! is equal to the minimum of the aux concentration and the cell
334  ! concentration.
335  if (qbnd >= dzero) then
336  omega = dzero ! rhs (ctmp is aux value)
337  else
338  if (ctmp < this%cnew(n)) then
339  omega = dzero ! rhs (ctmp is aux value)
340  else
341  omega = done ! lhs (ctmp is cell concentration)
342  ctmp = this%cnew(n)
343  end if
344  end if
345  end if
346  !
347  ! -- Add terms based on qbnd sign
348  if (qbnd <= dzero) then
349  hcoftmp = qbnd * omega * this%eqnsclfac
350  else
351  rhstmp = -qbnd * ctmp * (done - omega) * this%eqnsclfac
352  end if
353  !
354  ! -- end of active ibound
355  end if
356  !
357  ! -- set requested values
358  if (present(hcofval)) hcofval = hcoftmp
359  if (present(rhsval)) rhsval = rhstmp
360  if (present(rrate)) rrate = hcoftmp * ctmp - rhstmp
361  if (present(cssm)) cssm = ctmp
362  if (present(qssm)) qssm = qbnd
363  !
364  ! -- Return
365  return

Variable Documentation

◆ ftype

character(len=lenftype) tspssmmodule::ftype = 'SSM'

Definition at line 28 of file tsp-ssm.f90.

28  character(len=LENFTYPE) :: ftype = 'SSM'

◆ text

character(len=lenpackagename) tspssmmodule::text = ' SOURCE-SINK MIX'

Definition at line 29 of file tsp-ssm.f90.

29  character(len=LENPACKAGENAME) :: text = ' SOURCE-SINK MIX'