MODFLOW 6  version 6.5.0.dev2
MODFLOW 6 Code Documentation
Data Types | Functions/Subroutines
gwfvscmodule Module Reference

Data Types

type  concentrationpointer
 
type  gwfvsctype
 

Functions/Subroutines

real(dp) function calc_visc (ivisc, viscref, dviscdc, cviscref, conc, a2, a3, a4)
 Generic function to calculate changes in fluid viscosity using a linear formulation. More...
 
subroutine, public vsc_cr (vscobj, name_model, inunit, iout)
 @ brief Create a new package object More...
 
subroutine vsc_df (this, dis, vsc_input)
 @ brief Define viscosity package options and dimensions More...
 
subroutine vsc_ar (this, ibound)
 @ brief Allocate and read method for viscosity package More...
 
subroutine vsc_ar_bnd (this, packobj)
 Activate viscosity in advanced packages. More...
 
subroutine set_npf_pointers (this)
 Set pointers to NPF variables. More...
 
subroutine vsc_rp (this)
 @ brief Read new period data in viscosity package More...
 
subroutine vsc_ad (this)
 @ brief Advance the viscosity package More...
 
subroutine vsc_ad_bnd (this, packobj, hnew)
 Advance the boundary packages when viscosity is active. More...
 
subroutine vsc_ad_standard_bnd (packobj, hnew, visc, viscref, locelev, locvisc, locconc, dviscdc, cviscref, ivisc, a2, a3, a4, ctemp)
 advance ghb while accounting for viscosity More...
 
subroutine vsc_ad_sfr (packobj, visc, viscref, elev, locvisc, locconc, dviscdc, cviscref, ivisc, a2, a3, a4, ctemp)
 Update sfr-related viscosity ratios. More...
 
subroutine vsc_ad_lak (packobj, visc, viscref, elev, locvisc, locconc, dviscdc, cviscref, ivisc, a2, a3, a4, ctemp)
 Update lak-related viscosity ratios. More...
 
subroutine vsc_ad_maw (packobj, visc, viscref, elev, locvisc, locconc, dviscdc, cviscref, ivisc, a2, a3, a4, ctemp)
 Update maw-related viscosity ratios. More...
 
real(dp) function update_bnd_cond (bndvisc, viscref, spcfdcond)
 Apply viscosity to the conductance term. More...
 
real(dp) function calc_vsc_ratio (viscref, bndvisc)
 calculate and return the viscosity ratio More...
 
real(dp) function calc_bnd_viscosity (n, locvisc, locconc, viscref, dviscdc, cviscref, ctemp, ivisc, a2, a3, a4, auxvar)
 @ brief Calculate the boundary viscosity More...
 
subroutine get_visc_ratio (this, n, m, gwhdn, gwhdm, viscratio)
 Calculate the viscosity ratio. More...
 
subroutine calc_q_visc (this, cellid, viscratio)
 Account for viscosity in the aquiferhorizontal flow barriers. More...
 
subroutine update_k_with_vsc (this)
 Appled the viscosity ratio (mu_o/mu) to the hydraulic conductivity. More...
 
subroutine vsc_set_changed_at (this, kper, kstp)
 Mark K changes as having occurred at (kper, kstp) More...
 
subroutine vsc_ot_dv (this, idvfl)
 @ brief Output viscosity package dependent-variable terms. More...
 
subroutine vsc_da (this)
 @ brief Deallocate viscosity package memory More...
 
subroutine read_dimensions (this)
 @ brief Read dimensions More...
 
subroutine read_packagedata (this)
 @ brief Read data for package More...
 
subroutine set_packagedata (this, input_data)
 Sets package data instead of reading from file. More...
 
subroutine vsc_calcvisc (this)
 Calculate fluid viscosity. More...
 
subroutine allocate_scalars (this)
 @ brief Allocate scalars More...
 
subroutine allocate_arrays (this, nodes)
 @ brief Allocate arrays More...
 
subroutine read_options (this)
 @ brief Read Options block More...
 
subroutine set_options (this, input_data)
 Sets options as opposed to reading them from a file. More...
 
subroutine set_concentration_pointer (this, modelname, conc, icbund, istmpr)
 @ brief Set pointers to concentration(s) More...
 

Function/Subroutine Documentation

◆ allocate_arrays()

subroutine gwfvscmodule::allocate_arrays ( class(gwfvsctype)  this,
integer(i4b), intent(in)  nodes 
)

Allocate and initialize arrays for the VSC package.

Definition at line 1271 of file gwf-vsc.f90.

1272  ! -- dummy
1273  class(GwfVscType) :: this
1274  integer(I4B), intent(in) :: nodes
1275  ! -- local
1276  integer(I4B) :: i
1277  !
1278  ! -- Allocate
1279  call mem_allocate(this%visc, nodes, 'VISC', this%memoryPath)
1280  call mem_allocate(this%ivisc, this%nviscspecies, 'IVISC', this%memoryPath)
1281  call mem_allocate(this%dviscdc, this%nviscspecies, 'DRHODC', &
1282  this%memoryPath)
1283  call mem_allocate(this%cviscref, this%nviscspecies, 'CRHOREF', &
1284  this%memoryPath)
1285  call mem_allocate(this%ctemp, this%nviscspecies, 'CTEMP', this%memoryPath)
1286  allocate (this%cmodelname(this%nviscspecies))
1287  allocate (this%cauxspeciesname(this%nviscspecies))
1288  allocate (this%modelconc(this%nviscspecies))
1289  !
1290  ! -- Initialize
1291  do i = 1, nodes
1292  this%visc(i) = this%viscref
1293  end do
1294  !
1295  ! -- Initialize nviscspecies arrays
1296  do i = 1, this%nviscspecies
1297  this%ivisc(i) = 1
1298  this%dviscdc(i) = dzero
1299  this%cviscref(i) = dzero
1300  this%ctemp(i) = dzero
1301  this%cmodelname(i) = ''
1302  this%cauxspeciesname(i) = ''
1303  end do
1304  !
1305  ! -- Return
1306  return

◆ allocate_scalars()

subroutine gwfvscmodule::allocate_scalars ( class(gwfvsctype)  this)
private

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

Definition at line 1228 of file gwf-vsc.f90.

1229  ! -- modules
1230  use constantsmodule, only: dzero, dten, dep3
1231  ! -- dummy
1232  class(GwfVscType) :: this
1233  !
1234  ! -- allocate scalars in NumericalPackageType
1235  call this%NumericalPackageType%allocate_scalars()
1236  !
1237  ! -- Allocate
1238  call mem_allocate(this%thermivisc, 'THERMIVISC', this%memoryPath)
1239  call mem_allocate(this%idxtmpr, 'IDXTMPR', this%memoryPath)
1240  call mem_allocate(this%ioutvisc, 'IOUTVISC', this%memoryPath)
1241  call mem_allocate(this%ireadelev, 'IREADELEV', this%memoryPath)
1242  call mem_allocate(this%iconcset, 'ICONCSET', this%memoryPath)
1243  call mem_allocate(this%viscref, 'VISCREF', this%memoryPath)
1244  call mem_allocate(this%a2, 'A2', this%memoryPath)
1245  call mem_allocate(this%a3, 'A3', this%memoryPath)
1246  call mem_allocate(this%a4, 'A4', this%memoryPath)
1247  !
1248  call mem_allocate(this%nviscspecies, 'NVISCSPECIES', this%memoryPath)
1249  !
1250  ! -- Initialize
1251  this%thermivisc = 0
1252  this%idxtmpr = 0
1253  this%ioutvisc = 0
1254  this%ireadelev = 0
1255  this%iconcset = 0
1256  this%viscref = dep3
1257  this%A2 = dten
1258  this%A3 = 248.37_dp
1259  this%A4 = 133.15_dp
1260  !
1261  this%nviscspecies = 0
1262  !
1263  ! -- Return
1264  return
constantsmodule
This module contains simulation constants.
Definition: Constants.f90:9
constantsmodule::dep3
real(dp), parameter dep3
real constant 1000
Definition: Constants.f90:87
constantsmodule::dzero
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:64
constantsmodule::dten
real(dp), parameter dten
real constant 10
Definition: Constants.f90:83

◆ calc_bnd_viscosity()

real(dp) function gwfvscmodule::calc_bnd_viscosity ( integer(i4b), intent(in)  n,
integer(i4b), intent(in)  locvisc,
integer(i4b), dimension(:), intent(in)  locconc,
real(dp), intent(in)  viscref,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
real(dp), dimension(:), intent(inout)  ctemp,
integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4,
real(dp), dimension(:, :), intent(in)  auxvar 
)
private

Return the viscosity of the boundary package using one of the options in the following order of priority:

  1. Assign as aux variable in column with name 'VISCOSITY'
  2. Calculate using viscosity equation and nviscspecies aux columns
  3. If neither of those, then assign as viscref !!

Definition at line 775 of file gwf-vsc.f90.

777  ! -- modules
778  ! -- dummy
779  integer(I4B), intent(in) :: n
780  integer(I4B), intent(in) :: locvisc
781  real(DP), intent(in) :: a2, a3, a4
782  integer(I4B), dimension(:), intent(in) :: ivisc
783  integer(I4B), dimension(:), intent(in) :: locconc
784  real(DP), intent(in) :: viscref
785  real(DP), dimension(:), intent(in) :: dviscdc
786  real(DP), dimension(:), intent(in) :: cviscref
787  real(DP), dimension(:), intent(inout) :: ctemp
788  real(DP), dimension(:, :), intent(in) :: auxvar
789  ! -- Return
790  real(DP) :: viscbnd
791  ! -- local
792  integer(I4B) :: i
793  !
794  ! -- assign boundary viscosity based on one of three options
795  if (locvisc > 0) then
796  ! -- assign viscosity to an aux column named 'VISCOSITY'
797  viscbnd = auxvar(locvisc, n)
798  else if (locconc(1) > 0) then
799  ! -- calculate viscosity using one or more concentration auxcolumns
800  do i = 1, size(locconc)
801  ctemp(i) = dzero
802  if (locconc(i) > 0) then
803  ctemp(i) = auxvar(locconc(i), n)
804  end if
805  end do
806  viscbnd = calc_visc(ivisc, viscref, dviscdc, cviscref, ctemp, a2, a3, a4)
807  else
808  ! -- neither of the above, so assign as viscref
809  viscbnd = viscref
810  end if
811  !
812  ! -- Return
813  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ calc_q_visc()

subroutine gwfvscmodule::calc_q_visc ( class(gwfvsctype)  this,
integer(i4b), intent(in)  cellid,
real(dp), intent(inout)  viscratio 
)

Will return the viscosity associated with the upgradient node (cell) to the HFB package for adjusting the hydraulic characteristic (hydchr) of the barrier

Definition at line 852 of file gwf-vsc.f90.

853  ! -- dummy variables
854  class(GwfVscType) :: this
855  integer(I4B), intent(in) :: cellid
856  ! -- Return
857  real(DP), intent(inout) :: viscratio
858  ! -- local
859  real(DP) :: visc
860  !
861  ! -- Retrieve viscosity for the passed node number
862  visc = this%visc(cellid)
863  !
864  ! -- Calculate the viscosity ratio for the
865  viscratio = calc_vsc_ratio(this%viscref, visc)
866  !
867  ! -- Return
868  return
Here is the call graph for this function:

◆ calc_visc()

real(dp) function gwfvscmodule::calc_visc ( integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  viscref,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
real(dp), dimension(:), intent(in)  conc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4 
)
private

Definition at line 97 of file gwf-vsc.f90.

99  ! -- dummy
100  integer(I4B), dimension(:), intent(in) :: ivisc
101  real(DP), intent(in) :: viscref
102  real(DP), dimension(:), intent(in) :: dviscdc
103  real(DP), dimension(:), intent(in) :: cviscref
104  real(DP), dimension(:), intent(in) :: conc
105  real(DP), intent(in) :: a2, a3, a4
106  ! -- return
107  real(DP) :: visc
108  ! -- local
109  integer(I4B) :: nviscspec
110  integer(I4B) :: i
111  real(DP) :: mu_t
112  real(DP) :: expon
113  !
114  nviscspec = size(dviscdc)
115  visc = viscref
116  !
117  do i = 1, nviscspec
118  if (ivisc(i) == 1) then
119  visc = visc + dviscdc(i) * (conc(i) - cviscref(i))
120  else
121  expon = -1 * a3 * ((conc(i) - cviscref(i)) / &
122  ((conc(i) + a4) * (cviscref(i) + a4)))
123  mu_t = viscref * a2**expon
124  ! Order matters!! (This assumes we apply the temperature correction after
125  ! accounting for solute concentrations)
126  ! If a nonlinear correction is applied, then b/c it takes into account
127  ! viscref, need to subtract it in this case
128  ! At most, there will only ever be 1 nonlinear correction
129  visc = (visc - viscref) + mu_t
130  end if
131  ! end if
132  end do
133  !
134  ! -- Return
135  return
Here is the caller graph for this function:

◆ calc_vsc_ratio()

real(dp) function gwfvscmodule::calc_vsc_ratio ( real(dp), intent(in)  viscref,
real(dp), intent(in)  bndvisc 
)
private

Definition at line 754 of file gwf-vsc.f90.

755  ! -- dummy
756  real(DP), intent(in) :: viscref
757  real(DP), intent(in) :: bndvisc
758  ! -- local
759  real(DP) :: viscratio
760  !
761  viscratio = viscref / bndvisc
762  !
763  ! -- Return
764  return
Here is the caller graph for this function:

◆ get_visc_ratio()

subroutine gwfvscmodule::get_visc_ratio ( class(gwfvsctype)  this,
integer(i4b), intent(in)  n,
integer(i4b), intent(in)  m,
real(dp), intent(in)  gwhdn,
real(dp), intent(in)  gwhdm,
real(dp), intent(inout)  viscratio 
)
private

Calculate the viscosity ratio applied to the hydraulic characteristic provided by the user. The viscosity ratio is applicable only when the hydraulic characteristic is specified as positive and will not be applied when the hydchr is negative

Definition at line 823 of file gwf-vsc.f90.

824  ! -- modules
825  use constantsmodule, only: done
826  ! -- dummy
827  class(GwfVscType) :: this
828  integer(I4B), intent(in) :: n, m
829  real(DP), intent(in) :: gwhdn, gwhdm
830  real(DP), intent(inout) :: viscratio
831  ! -- loca
832  integer(I4B) :: cellid
833  !
834  viscratio = done
835  if (gwhdm > gwhdn) then
836  cellid = m
837  else if (gwhdn >= gwhdm) then
838  cellid = n
839  end if
840  call this%calc_q_visc(cellid, viscratio)
841  !
842  ! -- Return
843  return
constantsmodule::done
real(dp), parameter done
real constant 1
Definition: Constants.f90:75

◆ read_dimensions()

subroutine gwfvscmodule::read_dimensions ( class(gwfvsctype), intent(inout)  this)
private

Read dimensions for the viscosity package

Definition at line 1017 of file gwf-vsc.f90.

1018  ! -- modules
1019  ! -- dummy
1020  class(GwfVscType), intent(inout) :: this
1021  ! -- local
1022  character(len=LINELENGTH) :: errmsg, keyword
1023  integer(I4B) :: ierr
1024  logical :: isfound, endOfBlock
1025  !
1026  ! -- get dimensions block
1027  call this%parser%GetBlock('DIMENSIONS', isfound, ierr, &
1028  supportopenclose=.true.)
1029  !
1030  ! -- parse dimensions block if detected
1031  if (isfound) then
1032  write (this%iout, '(/1x,a)') 'Processing VSC DIMENSIONS block'
1033  do
1034  call this%parser%GetNextLine(endofblock)
1035  if (endofblock) exit
1036  call this%parser%GetStringCaps(keyword)
1037  select case (keyword)
1038  case ('NVISCSPECIES')
1039  this%nviscspecies = this%parser%GetInteger()
1040  write (this%iout, '(4x,a,i0)') 'NVISCSPECIES = ', this%nviscspecies
1041  case default
1042  write (errmsg, '(a,a)') &
1043  'Unknown VSC dimension: ', trim(keyword)
1044  call store_error(errmsg)
1045  call this%parser%StoreErrorUnit()
1046  end select
1047  end do
1048  write (this%iout, '(1x,a)') 'End of VSC DIMENSIONS block'
1049  else
1050  call store_error('Required VSC DIMENSIONS block not found.')
1051  call this%parser%StoreErrorUnit()
1052  end if
1053  !
1054  ! -- check dimension
1055  if (this%nviscspecies < 1) then
1056  call store_error('NVISCSPECIES must be greater than zero.')
1057  call this%parser%StoreErrorUnit()
1058  end if
1059  !
1060  ! -- Return
1061  return
Here is the call graph for this function:

◆ read_options()

subroutine gwfvscmodule::read_options ( class(gwfvsctype)  this)
private

Reads the options block inside the VSC package.

Definition at line 1313 of file gwf-vsc.f90.

1314  ! -- modules
1315  use openspecmodule, only: access, form
1316  use inputoutputmodule, only: urword, getunit, urdaux, openfile
1317  ! -- dummy
1318  class(GwfVscType) :: this
1319  ! -- local
1320  character(len=LINELENGTH) :: warnmsg, errmsg, keyword, keyword2
1321  character(len=MAXCHARLEN) :: fname
1322  integer(I4B) :: ierr
1323  logical :: isfound, endOfBlock
1324  ! -- formats
1325  character(len=*), parameter :: fmtfileout = &
1326  "(1x, 'VSC', 1x, a, 1x, 'Will be saved to file: ', &
1327  &a, /4x, 'opened on unit: ', I7)"
1328  character(len=*), parameter :: fmtlinear = &
1329  "(/,1x,'Viscosity will vary linearly with temperature &
1330  &change ')"
1331  character(len=*), parameter :: fmtnonlinear = &
1332  "(/,1x,'Viscosity will vary non-linearly with temperature &
1333  &change ')"
1334  !
1335  ! -- get options block
1336  call this%parser%GetBlock('OPTIONS', isfound, ierr, &
1337  supportopenclose=.true., blockrequired=.false.)
1338  !
1339  ! -- parse options block if detected
1340  if (isfound) then
1341  write (this%iout, '(1x,a)') 'Processing VSC OPTIONS block'
1342  do
1343  call this%parser%GetNextLine(endofblock)
1344  if (endofblock) exit
1345  call this%parser%GetStringCaps(keyword)
1346  select case (keyword)
1347  case ('VISCREF')
1348  this%viscref = this%parser%GetDouble()
1349  write (this%iout, '(4x,a,1pg15.6)') &
1350  'Reference viscosity has been set to: ', &
1351  this%viscref
1352  case ('VISCOSITY')
1353  call this%parser%GetStringCaps(keyword)
1354  if (keyword == 'FILEOUT') then
1355  call this%parser%GetString(fname)
1356  this%ioutvisc = getunit()
1357  call openfile(this%ioutvisc, this%iout, fname, 'DATA(BINARY)', &
1358  form, access, 'REPLACE')
1359  write (this%iout, fmtfileout) &
1360  'VISCOSITY', fname, this%ioutvisc
1361  else
1362  errmsg = 'Optional VISCOSITY keyword must be '// &
1363  'followed by FILEOUT'
1364  call store_error(errmsg)
1365  end if
1366  case ('TEMPERATURE_SPECIES_NAME')
1367  call this%parser%GetStringCaps(this%name_temp_spec)
1368  write (this%iout, '(4x, a)') 'Temperature species name set to: '// &
1369  trim(this%name_temp_spec)
1370  case ('THERMAL_FORMULATION')
1371  call this%parser%GetStringCaps(keyword2)
1372  if (trim(adjustl(keyword2)) == 'LINEAR') this%thermivisc = 1
1373  if (trim(adjustl(keyword2)) == 'NONLINEAR') this%thermivisc = 2
1374  select case (this%thermivisc)
1375  case (1)
1376  write (this%iout, fmtlinear)
1377  case (2)
1378  write (this%iout, fmtnonlinear)
1379  end select
1380  case ('THERMAL_A2')
1381  this%a2 = this%parser%GetDouble()
1382  if (this%thermivisc == 2) then
1383  write (this%iout, '(4x,a,1pg15.6)') &
1384  'A2 in nonlinear viscosity formulation has been set to: ', &
1385  this%a2
1386  else
1387  write (this%iout, '(4x,a,/,4x,a,/,4x,a)') 'THERMAL_A2 specified by user &
1388  &in VSC Package input file. LINEAR viscosity ', 'formulation also &
1389  &specified. THERMAL_A2 will not affect ', 'viscosity calculations.'
1390  end if
1391  case ('THERMAL_A3')
1392  this%a3 = this%parser%GetDouble()
1393  if (this%thermivisc == 2) then
1394  write (this%iout, '(4x,a,1pg15.6)') &
1395  'A3 in nonlinear viscosity formulation has been set to: ', &
1396  this%a3
1397  else
1398  write (this%iout, '(4x,a,/,4x,a,/,4x,a)') 'THERMAL_A3 specified by user &
1399  &in VSC Package input file. LINEAR viscosity ', 'formulation also &
1400  &specified. THERMAL_A3 will not affect ', 'viscosity calculations.'
1401  end if
1402  case ('THERMAL_A4')
1403  this%a4 = this%parser%GetDouble()
1404  if (this%thermivisc == 2) then
1405  write (this%iout, '(4x,a,1pg15.6)') &
1406  'A4 in nonlinear viscosity formulation has been set to: ', &
1407  this%a4
1408  else
1409  write (this%iout, '(4x,a,/,4x,a,/,4x,a)') 'THERMAL_A4 specified by user &
1410  &in VSC Package input file. LINEAR viscosity ', 'formulation also &
1411  &specified. THERMAL_A4 will not affect ', 'viscosity calculations.'
1412  end if
1413  case default
1414  write (errmsg, '(a,a)') 'Unknown VSC option: ', &
1415  trim(keyword)
1416  call store_error(errmsg)
1417  call this%parser%StoreErrorUnit()
1418  end select
1419  end do
1420  !
1421  if (this%thermivisc == 1) then
1422  if (this%a2 == 0.0) then
1423  write (errmsg, '(a)') 'LINEAR option selected for varying &
1424  &viscosity with temperature, but A1, a surrogate for &
1425  &dVISC/dT, set equal to 0.0'
1426  call store_error(errmsg)
1427  end if
1428  end if
1429  if (this%thermivisc > 1) then
1430  if (this%a2 == 0) then
1431  write (warnmsg, '(a)') 'NONLINEAR option selected for &
1432  &varying viscosity with temperature, but A2 set equal to &
1433  &zero which may lead to unintended values for viscosity'
1434  call store_warning(errmsg)
1435  end if
1436  if (this%a3 == 0) then
1437  write (warnmsg, '(a)') 'NONLINEAR option selected for &
1438  &varying viscosity with temperature,, but A3 set equal to &
1439  &zero which may lead to unintended values for viscosity'
1440  call store_warning(warnmsg)
1441  end if
1442  if (this%a4 == 0) then
1443  write (warnmsg, '(a)') 'NONLINEAR option selected for &
1444  &varying viscosity with temperature, BUT A4 SET EQUAL TO &
1445  &zero which may lead to unintended values for viscosity'
1446  call store_warning(warnmsg)
1447  end if
1448  end if
1449  end if
1450  !
1451  write (this%iout, '(/,1x,a)') 'end of VSC options block'
1452  !
1453  ! -- Return
1454  return
inputoutputmodule::urdaux
subroutine, public urdaux(naux, inunit, iout, lloc, istart, istop, auxname, line, text)
Read auxiliary variables from an input line.
Definition: InputOutput.f90:1395
inputoutputmodule::getunit
integer(i4b) function, public getunit()
Get a free unit number.
Definition: InputOutput.f90:176
inputoutputmodule::openfile
subroutine, public openfile(iu, iout, fname, ftype, fmtarg_opt, accarg_opt, filstat_opt, mode_opt)
Open a file.
Definition: InputOutput.f90:30
inputoutputmodule::urword
subroutine, public urword(line, icol, istart, istop, ncode, n, r, iout, in)
Extract a word from a string.
Definition: InputOutput.f90:427
openspecmodule::access
character(len=20) access
Definition: OpenSpec.f90:7
openspecmodule::form
character(len=20) form
Definition: OpenSpec.f90:7
Here is the call graph for this function:

◆ read_packagedata()

subroutine gwfvscmodule::read_packagedata ( class(gwfvsctype)  this)
private

Method to read data for the viscosity package.

Definition at line 1068 of file gwf-vsc.f90.

1069  ! -- dummy
1070  class(GwfVscType) :: this
1071  ! -- local
1072  character(len=LINELENGTH) :: errmsg
1073  character(len=LINELENGTH) :: line
1074  integer(I4B) :: ierr
1075  integer(I4B) :: iviscspec
1076  logical :: isfound, endOfBlock
1077  logical :: blockrequired
1078  integer(I4B), dimension(:), allocatable :: itemp
1079  character(len=10) :: c10
1080  character(len=16) :: c16
1081  ! -- format
1082  character(len=*), parameter :: fmterr = &
1083  "('Invalid value for IRHOSPEC (',i0,') detected in VSC Package. &
1084  &IRHOSPEC must be > 0 and <= NVISCSPECIES, and duplicate values &
1085  &are not allowed.')"
1086  !
1087  ! -- initialize
1088  allocate (itemp(this%nviscspecies))
1089  itemp(:) = 0
1090  !
1091  ! -- get packagedata block
1092  blockrequired = .true.
1093  call this%parser%GetBlock('PACKAGEDATA', isfound, ierr, &
1094  blockrequired=blockrequired, &
1095  supportopenclose=.true.)
1096  !
1097  ! -- parse packagedata block
1098  if (isfound) then
1099  write (this%iout, '(1x,a)') 'Procesing VSC PACKAGEDATA block'
1100  do
1101  call this%parser%GetNextLine(endofblock)
1102  if (endofblock) exit
1103  iviscspec = this%parser%GetInteger()
1104  if (iviscspec < 1 .or. iviscspec > this%nviscspecies) then
1105  write (errmsg, fmterr) iviscspec
1106  call store_error(errmsg)
1107  end if
1108  if (itemp(iviscspec) /= 0) then
1109  write (errmsg, fmterr) iviscspec
1110  call store_error(errmsg)
1111  end if
1112  itemp(iviscspec) = 1
1113  !
1114  this%dviscdc(iviscspec) = this%parser%GetDouble()
1115  this%cviscref(iviscspec) = this%parser%GetDouble()
1116  call this%parser%GetStringCaps(this%cmodelname(iviscspec))
1117  call this%parser%GetStringCaps(this%cauxspeciesname(iviscspec))
1118  !
1119  if (this%cauxspeciesname(iviscspec) == this%name_temp_spec) then
1120  if (this%idxtmpr > 0) then
1121  write (errmsg, '(a)') 'More than one species in VSC input identified &
1122  &as '//trim(this%name_temp_spec)//'. Only one species may be &
1123  &designated to represent temperature.'
1124  call store_error(errmsg)
1125  else
1126  this%idxtmpr = iviscspec
1127  if (this%thermivisc == 2) then
1128  this%ivisc(iviscspec) = 2
1129  end if
1130  end if
1131  end if
1132  end do
1133  else
1134  call store_error('Required VSC PACKAGEDATA block not found.')
1135  call this%parser%StoreErrorUnit()
1136  end if
1137  !
1138  ! -- Check for errors.
1139  if (count_errors() > 0) then
1140  call this%parser%StoreErrorUnit()
1141  end if
1142  !
1143  ! -- write packagedata information
1144  write (this%iout, '(/,1x,a)') 'Summary of species information in VSC Package'
1145  write (this%iout, '(1a11,5a17)') &
1146  'Species', 'DVISCDC', 'CVISCREF', 'Model', 'AUXSPECIESNAME'
1147  do iviscspec = 1, this%nviscspecies
1148  write (c10, '(i0)') iviscspec
1149  line = ' '//adjustr(c10)
1150 
1151  write (c16, '(g15.6)') this%dviscdc(iviscspec)
1152  line = trim(line)//' '//adjustr(c16)
1153  write (c16, '(g15.6)') this%cviscref(iviscspec)
1154  line = trim(line)//' '//adjustr(c16)
1155  write (c16, '(a)') this%cmodelname(iviscspec)
1156  line = trim(line)//' '//adjustr(c16)
1157  write (c16, '(a)') this%cauxspeciesname(iviscspec)
1158  line = trim(line)//' '//adjustr(c16)
1159  write (this%iout, '(a)') trim(line)
1160  end do
1161  !
1162  ! -- deallocate
1163  deallocate (itemp)
1164  !
1165  write (this%iout, '(/,1x,a)') 'End of VSC PACKAGEDATA block'
1166  !
1167  ! -- Return
1168  return
Here is the call graph for this function:

◆ set_concentration_pointer()

subroutine gwfvscmodule::set_concentration_pointer ( class(gwfvsctype)  this,
character(len=lenmodelname), intent(in)  modelname,
real(dp), dimension(:), pointer  conc,
integer(i4b), dimension(:), pointer  icbund,
integer(i4b), intent(in), optional  istmpr 
)
private

Pass in a gwt model name, concentration array, and ibound, and store a pointer to these in the VSC package so that viscosity can be calculated from them. This routine is called from the gwfgwt exchange in the exg_ar() method.

Definition at line 1477 of file gwf-vsc.f90.

1478  ! -- dummy
1479  class(GwfVscType) :: this
1480  character(len=LENMODELNAME), intent(in) :: modelname
1481  real(DP), dimension(:), pointer :: conc
1482  integer(I4B), dimension(:), pointer :: icbund
1483  integer(I4B), optional, intent(in) :: istmpr
1484  ! -- local
1485  integer(I4B) :: i
1486  logical :: found
1487  !
1488  this%iconcset = 1
1489  found = .false.
1490  do i = 1, this%nviscspecies
1491  if (this%cmodelname(i) == modelname) then
1492  this%modelconc(i)%conc => conc
1493  this%modelconc(i)%icbund => icbund
1494  found = .true.
1495  exit
1496  end if
1497  end do
1498  !
1499  ! -- Return
1500  return

◆ set_npf_pointers()

subroutine gwfvscmodule::set_npf_pointers ( class(gwfvsctype)  this)

Set array and variable pointers from the NPF package for access by VSC.

Definition at line 315 of file gwf-vsc.f90.

316  ! -- dummy
317  class(GwfVscType) :: this
318  ! -- local
319  character(len=LENMEMPATH) :: npfMemoryPath
320  !
321  ! -- Set pointers to other package variables
322  ! -- NPF
323  npfmemorypath = create_mem_path(this%name_model, 'NPF')
324  call mem_setptr(this%k11, 'K11', npfmemorypath)
325  call mem_setptr(this%k22, 'K22', npfmemorypath)
326  call mem_setptr(this%k33, 'K33', npfmemorypath)
327  call mem_setptr(this%k11input, 'K11INPUT', npfmemorypath)
328  call mem_setptr(this%k22input, 'K22INPUT', npfmemorypath)
329  call mem_setptr(this%k33input, 'K33INPUT', npfmemorypath)
330  call mem_setptr(this%kchangeper, 'KCHANGEPER', npfmemorypath)
331  call mem_setptr(this%kchangestp, 'KCHANGESTP', npfmemorypath)
332  call mem_setptr(this%nodekchange, 'NODEKCHANGE', npfmemorypath)
333  !
334  ! -- Return
335  return
Here is the call graph for this function:

◆ set_options()

subroutine gwfvscmodule::set_options ( class(gwfvsctype)  this,
type(gwfvscinputdatatype), intent(in)  input_data 
)
Parameters
[in]input_datathe input data to be set

Definition at line 1459 of file gwf-vsc.f90.

1460  ! -- dummy
1461  class(GwfVscType) :: this
1462  type(GwfVscInputDataType), intent(in) :: input_data !< the input data to be set
1463  !
1464  this%viscref = input_data%viscref
1465  !
1466  ! -- Return
1467  return

◆ set_packagedata()

subroutine gwfvscmodule::set_packagedata ( class(gwfvsctype)  this,
type(gwfvscinputdatatype), intent(in)  input_data 
)
private
Parameters
thisthis vscoancy pkg
[in]input_datathe input data to be set

Definition at line 1173 of file gwf-vsc.f90.

1174  ! -- dummy
1175  class(GwfVscType) :: this !< this vscoancy pkg
1176  type(GwfVscInputDataType), intent(in) :: input_data !< the input data to be set
1177  ! -- local
1178  integer(I4B) :: ispec
1179  !
1180  do ispec = 1, this%nviscspecies
1181  this%dviscdc(ispec) = input_data%dviscdc(ispec)
1182  this%cviscref(ispec) = input_data%cviscref(ispec)
1183  this%cmodelname(ispec) = input_data%cmodelname(ispec)
1184  this%cauxspeciesname(ispec) = input_data%cauxspeciesname(ispec)
1185  end do
1186  !
1187  ! -- Return
1188  return

◆ update_bnd_cond()

real(dp) function gwfvscmodule::update_bnd_cond ( real(dp), intent(in)  bndvisc,
real(dp), intent(in)  viscref,
real(dp), intent(in)  spcfdcond 
)

When the viscosity package is active apply the viscosity ratio to the active boundary package's conductance term.

Definition at line 734 of file gwf-vsc.f90.

735  ! -- dummy
736  real(DP), intent(in) :: viscref
737  real(DP), intent(in) :: bndvisc
738  real(DP), intent(in) :: spcfdcond
739  ! -- local
740  real(DP) :: vscratio
741  real(DP) :: updatedcond
742  !
743  vscratio = calc_vsc_ratio(viscref, bndvisc)
744  !
745  ! -- calculate new conductance here
746  updatedcond = vscratio * spcfdcond
747  !
748  ! -- Return
749  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ update_k_with_vsc()

subroutine gwfvscmodule::update_k_with_vsc ( class(gwfvsctype)  this)
private

This routine called after updating the viscosity values using the latest concentration and/or temperature values. The ratio mu_o/mu, reference viscosity divided by the updated viscosity value, is multiplied by K for each cell.

Definition at line 878 of file gwf-vsc.f90.

879  ! -- dummy
880  class(GwfVscType) :: this
881  ! -- local
882  integer(I4B) :: n
883  real(DP) :: viscratio
884  !
885  ! -- For viscosity-based K's, apply change of K to K11 by starting with
886  ! user-specified K values and not the K's leftover from the last viscosity
887  ! update.
888  do n = 1, this%dis%nodes
889  call this%calc_q_visc(n, viscratio)
890  this%k11(n) = this%k11input(n) * viscratio
891  this%k22(n) = this%k22input(n) * viscratio
892  this%k33(n) = this%k33input(n) * viscratio
893  this%nodekchange(n) = 1
894  end do
895  !
896  ! -- Flag kchange
897  call this%vsc_set_changed_at(kper, kstp)
898  !
899  ! -- Return
900  return

◆ vsc_ad()

subroutine gwfvscmodule::vsc_ad ( class(gwfvsctype)  this)

Advance data in the VSC package. The method sets or advances time series, time array series, and observation data.

Definition at line 379 of file gwf-vsc.f90.

380  ! -- dummy
381  class(GwfVscType) :: this
382  !
383  ! -- update viscosity using the latest concentration/temperature
384  call this%vsc_calcvisc()
385  !
386  ! -- Return
387  return

◆ vsc_ad_bnd()

subroutine gwfvscmodule::vsc_ad_bnd ( class(gwfvsctype)  this,
class(bndtype), pointer  packobj,
real(dp), dimension(:), intent(in)  hnew 
)
private

Update the conductance values associate with inflow from a boundary when VSC package is active.

Definition at line 395 of file gwf-vsc.f90.

396  ! -- modules
397  use bndmodule, only: bndtype
398  ! -- dummy
399  class(GwfVscType) :: this
400  class(BndType), pointer :: packobj
401  real(DP), intent(in), dimension(:) :: hnew
402  ! -- local
403  integer(I4B) :: i, j
404  integer(I4B) :: n, locvisc, locelev
405  integer(I4B), dimension(:), allocatable :: locconc
406  !
407  ! -- initialize
408  locvisc = 0
409  locelev = 0
410  allocate (locconc(this%nviscspecies))
411  locconc(:) = 0
412  !
413  ! -- Add viscosity terms for conductance-dependent boundaries
414  do n = 1, packobj%naux
415  if (packobj%auxname(n) == 'VISCOSITY') then
416  locvisc = n
417  else if (packobj%auxname(n) == 'ELEVATION') then
418  locelev = n
419  end if
420  end do
421  !
422  ! -- find aux columns for conc (or temp.) that affect viscosity
423  do i = 1, this%nviscspecies
424  locconc(i) = 0
425  do j = 1, packobj%naux
426  if (this%cauxspeciesname(i) == packobj%auxname(j)) then
427  locconc(i) = j
428  exit
429  end if
430  end do
431  if (locconc(i) == 0) then
432  ! -- one not found, so don't use and mark all as 0
433  locconc(:) = 0
434  exit
435  end if
436  end do
437  !
438  ! -- apply viscosity terms to inflow from boundary based on package type
439  select case (packobj%filtyp)
440  case ('GHB', 'DRN', 'RIV')
441  !
442  ! -- general head, drain, and river boundary
443  call vsc_ad_standard_bnd(packobj, hnew, this%visc, this%viscref, &
444  locelev, locvisc, locconc, this%dviscdc, &
445  this%cviscref, this%ivisc, this%a2, this%a3, &
446  this%a4, this%ctemp)
447  case ('LAK')
448  !
449  ! -- lake
450  ! Update 'viscratios' internal to lak such that they are
451  ! automatically applied in the LAK calc_cond() routine
452  call vsc_ad_lak(packobj, this%visc, this%viscref, this%elev, locvisc, &
453  locconc, this%dviscdc, this%cviscref, this%ivisc, &
454  this%a2, this%a3, this%a4, this%ctemp)
455  case ('SFR')
456  !
457  ! -- streamflow routing
458  ! Update 'viscratios' internal to sfr such that they are
459  ! automatically applied in the SFR calc_cond() routine
460  call vsc_ad_sfr(packobj, this%visc, this%viscref, this%elev, locvisc, &
461  locconc, this%dviscdc, this%cviscref, this%ivisc, &
462  this%a2, this%a3, this%a4, this%ctemp)
463  case ('MAW')
464  !
465  ! -- multi-aquifer well
466  call vsc_ad_maw(packobj, this%visc, this%viscref, this%elev, locvisc, &
467  locconc, this%dviscdc, this%cviscref, this%ivisc, &
468  this%a2, this%a3, this%a4, this%ctemp)
469  case ('UZF')
470  !
471  ! -- unsaturated-zone flow
472  case default
473  !
474  ! -- nothing
475  end select
476  !
477  ! -- deallocate
478  deallocate (locconc)
479  !
480  ! -- Return
481  return
bndmodule
This module contains the base boundary package.
Definition: BoundaryPackage.f90:8
bndmodule::bndtype
@ brief BndType
Definition: BoundaryPackage.f90:48
Here is the call graph for this function:

◆ vsc_ad_lak()

subroutine gwfvscmodule::vsc_ad_lak ( class(bndtype), pointer  packobj,
real(dp), dimension(:), intent(in)  visc,
real(dp), intent(in)  viscref,
real(dp), dimension(:), intent(in)  elev,
integer(i4b), intent(in)  locvisc,
integer(i4b), dimension(:), intent(in)  locconc,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4,
real(dp), dimension(:), intent(inout)  ctemp 
)

When the viscosity package is active, update the viscosity ratio that is applied to the lakebed conductance calculated in the LAK package

Definition at line 620 of file gwf-vsc.f90.

622  ! -- modules
623  use bndmodule, only: bndtype
624  use lakmodule, only: laktype
625  class(BndType), pointer :: packobj
626  ! -- dummy
627  real(DP), intent(in) :: viscref
628  real(DP), intent(in) :: a2, a3, a4
629  integer(I4B), intent(in) :: locvisc
630  integer(I4B), dimension(:), intent(in) :: locconc
631  integer(I4B), dimension(:), intent(in) :: ivisc
632  real(DP), dimension(:), intent(in) :: visc
633  real(DP), dimension(:), intent(in) :: elev
634  real(DP), dimension(:), intent(in) :: dviscdc
635  real(DP), dimension(:), intent(in) :: cviscref
636  real(DP), dimension(:), intent(inout) :: ctemp
637  ! -- local
638  integer(I4B) :: n
639  integer(I4B) :: node
640  real(DP) :: visclak
641  !
642  ! -- update viscosity ratios for updating hyd. cond (and conductance)
643  select type (packobj)
644  type is (laktype)
645  do n = 1, packobj%nbound
646  !
647  ! -- get gwf node number
648  node = packobj%nodelist(n)
649  !
650  ! -- Check if boundary cell is active, cycle if not
651  if (packobj%ibound(node) <= 0) cycle
652  !
653  ! --
654  !
655  ! -- calculate the viscosity associated with the boundary
656  visclak = calc_bnd_viscosity(n, locvisc, locconc, viscref, dviscdc, &
657  cviscref, ctemp, ivisc, a2, a3, a4, &
658  packobj%auxvar)
659  !
660  ! -- fill lak relative viscosity into column 1 of viscratios
661  packobj%viscratios(1, n) = calc_vsc_ratio(viscref, visclak)
662  !
663  ! -- fill gwf relative viscosity into column 2 of viscratios
664  packobj%viscratios(2, n) = calc_vsc_ratio(viscref, visc(node))
665  end do
666  end select
667  !
668  ! -- Return
669  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vsc_ad_maw()

subroutine gwfvscmodule::vsc_ad_maw ( class(bndtype), pointer  packobj,
real(dp), dimension(:), intent(in)  visc,
real(dp), intent(in)  viscref,
real(dp), dimension(:), intent(in)  elev,
integer(i4b), intent(in)  locvisc,
integer(i4b), dimension(:), intent(in)  locconc,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4,
real(dp), dimension(:), intent(inout)  ctemp 
)

When the viscosity package is active, update the viscosity ratio that is applied to the conductance calculated in the MAW package

Definition at line 677 of file gwf-vsc.f90.

679  ! -- modules
680  use bndmodule, only: bndtype
681  use mawmodule, only: mawtype
682  class(BndType), pointer :: packobj
683  ! -- dummy
684  real(DP), intent(in) :: viscref
685  real(DP), intent(in) :: a2, a3, a4
686  integer(I4B), intent(in) :: locvisc
687  integer(I4B), dimension(:), intent(in) :: locconc
688  integer(I4B), dimension(:), intent(in) :: ivisc
689  real(DP), dimension(:), intent(in) :: visc
690  real(DP), dimension(:), intent(in) :: elev
691  real(DP), dimension(:), intent(in) :: dviscdc
692  real(DP), dimension(:), intent(in) :: cviscref
693  real(DP), dimension(:), intent(inout) :: ctemp
694  ! -- local
695  integer(I4B) :: n
696  integer(I4B) :: node
697  real(DP) :: viscmaw
698  !
699  ! -- update viscosity ratios for updating hyd. cond (and conductance)
700  select type (packobj)
701  type is (mawtype)
702  do n = 1, packobj%nbound
703  !
704  ! -- get gwf node number
705  node = packobj%nodelist(n)
706  !
707  ! -- Check if boundary cell is active, cycle if not
708  if (packobj%ibound(node) <= 0) cycle
709  !
710  ! --
711  !
712  ! -- calculate the viscosity associated with the boundary
713  viscmaw = calc_bnd_viscosity(n, locvisc, locconc, viscref, dviscdc, &
714  cviscref, ctemp, ivisc, a2, a3, a4, &
715  packobj%auxvar)
716  !
717  ! -- fill lak relative viscosity into column 1 of viscratios
718  packobj%viscratios(1, n) = calc_vsc_ratio(viscref, viscmaw)
719  !
720  ! -- fill gwf relative viscosity into column 2 of viscratios
721  packobj%viscratios(2, n) = calc_vsc_ratio(viscref, visc(node))
722  end do
723  end select
724  !
725  ! -- Return
726  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vsc_ad_sfr()

subroutine gwfvscmodule::vsc_ad_sfr ( class(bndtype), pointer  packobj,
real(dp), dimension(:), intent(in)  visc,
real(dp), intent(in)  viscref,
real(dp), dimension(:), intent(in)  elev,
integer(i4b), intent(in)  locvisc,
integer(i4b), dimension(:), intent(in)  locconc,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4,
real(dp), dimension(:), intent(inout)  ctemp 
)

When the viscosity package is active, update the viscosity ratio that is applied to the hydraulic conductivity specified in the SFR package

Definition at line 563 of file gwf-vsc.f90.

565  ! -- modules
566  use bndmodule, only: bndtype
567  use sfrmodule, only: sfrtype
568  class(BndType), pointer :: packobj
569  ! -- dummy
570  real(DP), intent(in) :: viscref
571  real(DP), intent(in) :: a2, a3, a4
572  integer(I4B), intent(in) :: locvisc
573  integer(I4B), dimension(:), intent(in) :: locconc
574  integer(I4B), dimension(:), intent(in) :: ivisc
575  real(DP), dimension(:), intent(in) :: visc
576  real(DP), dimension(:), intent(in) :: elev
577  real(DP), dimension(:), intent(in) :: dviscdc
578  real(DP), dimension(:), intent(in) :: cviscref
579  real(DP), dimension(:), intent(inout) :: ctemp
580  ! -- local
581  integer(I4B) :: n
582  integer(I4B) :: node
583  real(DP) :: viscsfr
584  !
585  ! -- update viscosity ratios for updating hyd. cond (and conductance)
586  select type (packobj)
587  type is (sfrtype)
588  do n = 1, packobj%nbound
589  !
590  ! -- get gwf node number
591  node = packobj%nodelist(n)
592  !
593  ! -- Check if boundary cell is active, cycle if not
594  if (packobj%ibound(node) <= 0) cycle
595  !
596  ! --
597  !
598  ! -- calculate the viscosity associated with the boundary
599  viscsfr = calc_bnd_viscosity(n, locvisc, locconc, viscref, dviscdc, &
600  cviscref, ctemp, ivisc, a2, a3, a4, &
601  packobj%auxvar)
602  !
603  ! -- fill sfr relative viscosity into column 1 of viscratios
604  packobj%viscratios(1, n) = calc_vsc_ratio(viscref, viscsfr)
605  !
606  ! -- fill gwf relative viscosity into column 2 of viscratios
607  packobj%viscratios(2, n) = calc_vsc_ratio(viscref, visc(node))
608  end do
609  end select
610  !
611  ! -- Return
612  return
sfrmodule
This module contains the SFR package methods.
Definition: gwf-sfr.f90:7
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vsc_ad_standard_bnd()

subroutine gwfvscmodule::vsc_ad_standard_bnd ( class(bndtype), pointer  packobj,
real(dp), dimension(:), intent(in)  hnew,
real(dp), dimension(:), intent(in)  visc,
real(dp), intent(in)  viscref,
integer(i4b), intent(in)  locelev,
integer(i4b), intent(in)  locvisc,
integer(i4b), dimension(:), intent(in)  locconc,
real(dp), dimension(:), intent(in)  dviscdc,
real(dp), dimension(:), intent(in)  cviscref,
integer(i4b), dimension(:), intent(in)  ivisc,
real(dp), intent(in)  a2,
real(dp), intent(in)  a3,
real(dp), intent(in)  a4,
real(dp), dimension(:), intent(inout)  ctemp 
)

When flow enters from ghb boundary type, take into account the effects of viscosity on the user-specified conductance terms

Definition at line 489 of file gwf-vsc.f90.

492  ! -- modules
493  use bndmodule, only: bndtype
494  use drnmodule, only: drntype
495  use rivmodule, only: rivtype
496  use ghbmodule, only: ghbtype
497  class(BndType), pointer :: packobj
498  ! -- dummy
499  real(DP), intent(in), dimension(:) :: hnew
500  real(DP), intent(in), dimension(:) :: visc
501  real(DP), intent(in) :: a2, a3, a4
502  real(DP), intent(in) :: viscref
503  integer(I4B), intent(in) :: locelev
504  integer(I4B), intent(in) :: locvisc
505  integer(I4B), dimension(:), intent(in) :: locconc
506  integer(I4B), dimension(:), intent(in) :: ivisc
507  real(DP), dimension(:), intent(in) :: dviscdc
508  real(DP), dimension(:), intent(in) :: cviscref
509  real(DP), dimension(:), intent(inout) :: ctemp
510  ! -- local
511  integer(I4B) :: n
512  integer(I4B) :: node
513  real(DP) :: viscbnd
514  !
515  ! -- Process density terms for each GHB
516  do n = 1, packobj%nbound
517  node = packobj%nodelist(n)
518  !
519  ! -- Check if boundary cell is active, cycle if not
520  if (packobj%ibound(node) <= 0) cycle
521  !
522  ! -- calculate the viscosity associated with the boundary
523  viscbnd = calc_bnd_viscosity(n, locvisc, locconc, viscref, dviscdc, &
524  cviscref, ctemp, ivisc, a2, a3, a4, &
525  packobj%auxvar)
526  !
527  ! -- update boundary conductance based on viscosity effects
528  select case (packobj%filtyp)
529  case ('DRN')
530  select type (packobj)
531  type is (drntype)
532  packobj%cond(n) = update_bnd_cond(viscbnd, viscref, &
533  packobj%condinput(n))
534  end select
535  case ('GHB')
536  select type (packobj)
537  type is (ghbtype)
538  packobj%cond(n) = update_bnd_cond(viscbnd, viscref, &
539  packobj%condinput(n))
540  end select
541  case ('RIV')
542  select type (packobj)
543  type is (rivtype)
544  packobj%cond(n) = update_bnd_cond(viscbnd, viscref, &
545  packobj%condinput(n))
546  end select
547  case default
548  packobj%bound(2, n) = update_bnd_cond(viscbnd, viscref, &
549  packobj%condinput(n))
550  end select
551  !
552  end do
553  !
554  ! -- Return
555  return
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vsc_ar()

subroutine gwfvscmodule::vsc_ar ( class(gwfvsctype)  this,
integer(i4b), dimension(:), pointer  ibound 
)
private

Generic method to allocate and read static data for the viscosity package available within the GWF model type.

Definition at line 223 of file gwf-vsc.f90.

224  ! -- dummy
225  class(GwfVscType) :: this
226  integer(I4B), dimension(:), pointer :: ibound
227  !
228  ! -- store pointers to arguments that were passed in
229  this%ibound => ibound
230  !
231  ! -- Set pointers to npf variables
232  call this%set_npf_pointers()
233  !
234  ! -- Return
235  return

◆ vsc_ar_bnd()

subroutine gwfvscmodule::vsc_ar_bnd ( class(gwfvsctype)  this,
class(bndtype), pointer  packobj 
)
private

Viscosity ar_bnd routine to activate viscosity in the advanced packages. This routine is called from gwf_ar() as it moves through each package

Definition at line 244 of file gwf-vsc.f90.

245  ! -- modules
246  use bndmodule, only: bndtype
247  use drnmodule, only: drntype
248  use ghbmodule, only: ghbtype
249  use rivmodule, only: rivtype
250  use lakmodule, only: laktype
251  use sfrmodule, only: sfrtype
252  use mawmodule, only: mawtype
253  ! -- dummy
254  class(GwfVscType) :: this
255  class(BndType), pointer :: packobj
256  !
257  ! -- Add density terms based on boundary package type
258  select case (packobj%filtyp)
259  case ('DRN')
260  !
261  ! -- activate viscosity for the drain package
262  select type (packobj)
263  type is (drntype)
264  call packobj%bnd_activate_viscosity()
265  end select
266  case ('GHB')
267  !
268  ! -- activate viscosity for the drain package
269  select type (packobj)
270  type is (ghbtype)
271  call packobj%bnd_activate_viscosity()
272  end select
273  case ('RIV')
274  !
275  ! -- activate viscosity for the drain package
276  select type (packobj)
277  type is (rivtype)
278  call packobj%bnd_activate_viscosity()
279  end select
280  case ('LAK')
281  !
282  ! -- activate viscosity for lake package
283  select type (packobj)
284  type is (laktype)
285  call packobj%lak_activate_viscosity()
286  end select
287  case ('SFR')
288  !
289  ! -- activate viscosity for sfr package
290  select type (packobj)
291  type is (sfrtype)
292  call packobj%sfr_activate_viscosity()
293  end select
294  case ('MAW')
295  !
296  ! -- activate viscosity for maw package
297  select type (packobj)
298  type is (mawtype)
299  call packobj%maw_activate_viscosity()
300  end select
301  case default
302  !
303  ! -- nothing
304  end select
305  !
306  ! -- Return
307  return

◆ vsc_calcvisc()

subroutine gwfvscmodule::vsc_calcvisc ( class(gwfvsctype)  this)
private

Calculates fluid viscosity based on concentration or temperature

Definition at line 1196 of file gwf-vsc.f90.

1197  ! -- dummy
1198  class(GwfVscType) :: this
1199  ! -- local
1200  integer(I4B) :: n
1201  integer(I4B) :: i
1202  !
1203  ! -- Calculate the viscosity using the specified concentration and/or
1204  ! temperature arrays
1205  do n = 1, this%dis%nodes
1206  do i = 1, this%nviscspecies
1207  if (this%modelconc(i)%icbund(n) == 0) then
1208  this%ctemp = dzero
1209  else
1210  this%ctemp(i) = this%modelconc(i)%conc(n)
1211  end if
1212  end do
1213  !
1214  this%visc(n) = calc_visc(this%ivisc, this%viscref, this%dviscdc, &
1215  this%cviscref, this%ctemp, this%a2, &
1216  this%a3, this%a4)
1217  end do
1218  !
1219  ! -- Return
1220  return
Here is the call graph for this function:

◆ vsc_cr()

subroutine, public gwfvscmodule::vsc_cr ( type(gwfvsctype), pointer  vscobj,
character(len=*), intent(in)  name_model,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout 
)

Create a new VSC Package object.

Definition at line 142 of file gwf-vsc.f90.

143  ! -- dummy
144  type(GwfVscType), pointer :: vscobj
145  character(len=*), intent(in) :: name_model
146  integer(I4B), intent(in) :: inunit
147  integer(I4B), intent(in) :: iout
148  !
149  ! -- Create the object
150  allocate (vscobj)
151  !
152  ! -- create name and memory path
153  call vscobj%set_names(1, name_model, 'VSC', 'VSC')
154  !
155  ! -- Allocate scalars
156  call vscobj%allocate_scalars()
157  !
158  ! -- Set variables
159  vscobj%inunit = inunit
160  vscobj%iout = iout
161  !
162  ! -- Initialize block parser
163  call vscobj%parser%Initialize(vscobj%inunit, vscobj%iout)
164  !
165  ! -- Return
166  return
Here is the caller graph for this function:

◆ vsc_da()

subroutine gwfvscmodule::vsc_da ( class(gwfvsctype)  this)
private

Deallocate viscosity package scalars and arrays.

Definition at line 967 of file gwf-vsc.f90.

968  ! -- dummy
969  class(GwfVscType) :: this
970  !
971  ! -- Deallocate arrays if package was active
972  if (this%inunit > 0) then
973  call mem_deallocate(this%visc)
974  call mem_deallocate(this%ivisc)
975  call mem_deallocate(this%dviscdc)
976  call mem_deallocate(this%cviscref)
977  call mem_deallocate(this%ctemp)
978  deallocate (this%cmodelname)
979  deallocate (this%cauxspeciesname)
980  deallocate (this%modelconc)
981  end if
982  !
983  ! -- Scalars
984  call mem_deallocate(this%thermivisc)
985  call mem_deallocate(this%idxtmpr)
986  call mem_deallocate(this%ioutvisc)
987  call mem_deallocate(this%ireadelev)
988  call mem_deallocate(this%iconcset)
989  call mem_deallocate(this%viscref)
990  call mem_deallocate(this%nviscspecies)
991  call mem_deallocate(this%a2)
992  call mem_deallocate(this%a3)
993  call mem_deallocate(this%a4)
994  !
995  ! -- Nullify pointers to other package variables
996  nullify (this%k11)
997  nullify (this%k22)
998  nullify (this%k33)
999  nullify (this%k11input)
1000  nullify (this%k22input)
1001  nullify (this%k33input)
1002  nullify (this%kchangeper)
1003  nullify (this%kchangestp)
1004  nullify (this%nodekchange)
1005  !
1006  ! -- deallocate parent
1007  call this%NumericalPackageType%da()
1008  !
1009  ! -- Return
1010  return

◆ vsc_df()

subroutine gwfvscmodule::vsc_df ( class(gwfvsctype)  this,
class(disbasetype), intent(in), pointer  dis,
type(gwfvscinputdatatype), intent(in), optional  vsc_input 
)
private

Define viscosity package options and dimensions

Parameters
thisthis viscosity package
[in]dispointer to discretization
[in]vsc_inputoptional vsc input data, otherwise read from file

Definition at line 173 of file gwf-vsc.f90.

174  ! -- dummy
175  class(GwfVscType) :: this !< this viscosity package
176  class(DisBaseType), pointer, intent(in) :: dis !< pointer to discretization
177  type(GwfVscInputDataType), optional, intent(in) :: vsc_input !< optional vsc input data, otherwise read from file
178  ! -- formats
179  character(len=*), parameter :: fmtvsc = &
180  "(1x,/1x,'VSC -- Viscosity Package, version 1, 11/15/2022', &
181  &' input read from unit ', i0, //)"
182  !
183  ! --print a message identifying the viscosity package
184  write (this%iout, fmtvsc) this%inunit
185  !
186  ! -- store pointers to arguments that were passed in
187  this%dis => dis
188  !
189  if (.not. present(vsc_input)) then
190  !
191  ! -- Read viscosity options
192  call this%read_options()
193  !
194  ! -- Read viscosity dimensions
195  call this%read_dimensions()
196  else
197  ! set from input data instead
198  call this%set_options(vsc_input)
199  this%nviscspecies = vsc_input%nviscspecies
200  end if
201  !
202  ! -- Allocate arrays
203  call this%allocate_arrays(dis%nodes)
204  !
205  if (.not. present(vsc_input)) then
206  !
207  ! -- Read viscosity packagedata
208  call this%read_packagedata()
209  else
210  ! set from input data instead
211  call this%set_packagedata(vsc_input)
212  end if
213  !
214  ! -- Return
215  return

◆ vsc_ot_dv()

subroutine gwfvscmodule::vsc_ot_dv ( class(gwfvsctype)  this,
integer(i4b), intent(in)  idvfl 
)
private

Save calculated viscosity array to binary file

Definition at line 925 of file gwf-vsc.f90.

926  ! -- dummy
927  class(GwfVscType) :: this
928  integer(I4B), intent(in) :: idvfl
929  ! -- local
930  character(len=1) :: cdatafmp = ' ', editdesc = ' '
931  integer(I4B) :: ibinun
932  integer(I4B) :: iprint
933  integer(I4B) :: nvaluesp
934  integer(I4B) :: nwidthp
935  real(DP) :: dinact
936  !
937  ! -- Set unit number for viscosity output
938  if (this%ioutvisc /= 0) then
939  ibinun = 1
940  else
941  ibinun = 0
942  end if
943  if (idvfl == 0) ibinun = 0
944  !
945  ! -- save viscosity array
946  if (ibinun /= 0) then
947  iprint = 0
948  dinact = dhnoflo
949  !
950  ! -- write viscosity to binary file
951  if (this%ioutvisc /= 0) then
952  ibinun = this%ioutvisc
953  call this%dis%record_array(this%visc, this%iout, iprint, ibinun, &
954  ' VISCOSITY', cdatafmp, nvaluesp, &
955  nwidthp, editdesc, dinact)
956  end if
957  end if
958  !
959  ! -- Return
960  return

◆ vsc_rp()

subroutine gwfvscmodule::vsc_rp ( class(gwfvsctype)  this)
private

Method to read and prepare period data for the VSC package.

Definition at line 342 of file gwf-vsc.f90.

343  ! -- modules
344  use tdismodule, only: kstp, kper
345  ! -- dummy
346  class(GwfVscType) :: this
347  ! -- local
348  character(len=LINELENGTH) :: errmsg
349  integer(I4B) :: i
350  ! -- formats
351  character(len=*), parameter :: fmtc = &
352  "('Viscosity Package does not have a concentration set &
353  &for species ',i0,'. One or more model names may be specified &
354  &incorrectly in the PACKAGEDATA block or a GWF-GWT exchange may need &
355  &to be activated.')"
356  !
357  ! -- Check to make sure all concentration pointers have been set
358  if (kstp * kper == 1) then
359  do i = 1, this%nviscspecies
360  if (.not. associated(this%modelconc(i)%conc)) then
361  write (errmsg, fmtc) i
362  call store_error(errmsg)
363  end if
364  end do
365  if (count_errors() > 0) then
366  call this%parser%StoreErrorUnit()
367  end if
368  end if
369  !
370  ! -- Return
371  return
tdismodule::kstp
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
tdismodule::kper
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23
Here is the call graph for this function:

◆ vsc_set_changed_at()

subroutine gwfvscmodule::vsc_set_changed_at ( class(gwfvsctype)  this,
integer(i4b), intent(in)  kper,
integer(i4b), intent(in)  kstp 
)
private

Procedure called by VSC code when K updated due to viscosity changes. K values changed at (kper, kstp).

Definition at line 908 of file gwf-vsc.f90.

909  ! -- dummy variables
910  class(GwfVscType) :: this
911  integer(I4B), intent(in) :: kper
912  integer(I4B), intent(in) :: kstp
913  !
914  this%kchangeper = kper
915  this%kchangestp = kstp
916  !
917  ! -- Return
918  return