MODFLOW 6 version 6.9.0.dev0
USGS Modular Hydrologic Model
Loading...
Searching...
No Matches
Data Types | Functions/Subroutines | Variables
gwtistmodule Module Reference

– @ brief Immobile Storage and Transfer (IST) Module More...

Data Types

type  gwtisttype
 @ brief Immobile storage and transfer More...

Functions/Subroutines

subroutine, public ist_create (packobj, id, ibcnum, inunit, iout, namemodel, pakname, mempath, fmi, mst)
 @ brief Create a new package object
subroutine ist_ar (this)
 @ brief Allocate and read method for package
subroutine ist_rp (this)
 @ brief Read and prepare method for package
subroutine ist_ad (this)
 @ brief Advance the ist package
subroutine ist_fc (this, rhs, ia, idxglo, matrix_sln)
 @ brief Fill coefficient method for package
subroutine ist_cq (this, x, flowja, iadv)
 @ brief Calculate package flows.
subroutine ist_calc_csrb (this, cim)
 @ brief Calculate immobile sorbed concentration
subroutine ist_bd (this, model_budget)
 @ brief Add package flows to model budget.
subroutine ist_ot_model_flows (this, icbcfl, ibudfl, icbcun, imap)
 @ brief Output model flow terms.
subroutine ist_ot_dv (this, idvsave, idvprint)
 @ brief Output dependent variables.
subroutine output_immobile_concentration (this, idvsave, idvprint)
 @ brief Output immobile domain aqueous concentration.
subroutine output_immobile_sorbate_concentration (this, idvsave, idvprint)
 @ brief Output immobile domain sorbate concentration.
subroutine ist_ot_bdsummary (this, kstp, kper, iout, ibudfl)
 @ brief Output IST package budget summary.
subroutine ist_da (this)
 @ brief Deallocate package memory
subroutine allocate_scalars (this)
 @ brief Allocate package scalars
subroutine ist_allocate_arrays (this)
 @ brief Allocate package arrays
subroutine source_options (this)
 @ brief Source options for package
subroutine log_options (this, found, cim6_fname, budget_fname, budgetcsv_fname, sorbate_fname)
 Log user options to list file.
subroutine read_options (this)
 @ brief Read options for package
subroutine source_data (this)
 @ brief Source data for package
subroutine log_data (this, found)
 Log user data to list file.
subroutine ist_read_dimensions (this)
 @ brief Read dimensions for package
real(dp) function get_thetaim (this, node)
 @ brief Return thetaim
subroutine get_ddterm (thetaim, vcell, delt, swtpdt, volfracim, rhobim, kdnew, kdold, lambda1im, lambda2im, gamma1im, gamma2im, zetaim, ddterm, f)
 @ brief Calculate immobile domain equation terms
subroutine get_hcofrhs (ddterm, f, cimold, hcof, rhs)
 @ brief Calculate the hcof and rhs terms for immobile domain
real(dp) function get_ddconc (ddterm, f, cimold, cnew)
 @ brief Calculate the immobile domain concentration
subroutine accumulate_budterm (budterm, ddterm, cimnew, cimold, cnew, idcy)
 @ brief Calculate the immobile domain budget terms
real(dp) function get_freundlich_kd (conc, kf, a)
 @ brief Get effective Freundlich distribution coefficient
real(dp) function get_langmuir_kd (conc, kl, sbar)
 @ brief Get effective Langmuir distribution coefficient

Variables

character(len=lenftype) ftype = 'IST'
character(len=lenpackagename) text = ' IMMOBILE DOMAIN'
integer(i4b), parameter nbditems = 5
character(len=lenbudtxt), dimension(nbditems) budtxt

Detailed Description

The GwtIstModule is contains the GwtIstType, which is the derived type responsible for adding the effects of an immobile domain. In addition to representing transfer of mass between the mobile and immobile domain, the IST Package also represents the following processes within the immobile domain

  1. Changes in dissolved solute mass
  2. Decay of dissolved solute mass
  3. Sorption
  4. Decay of sorbed solute mass

Function/Subroutine Documentation

◆ accumulate_budterm()

subroutine gwtistmodule::accumulate_budterm ( real(dp), dimension(:, :), intent(inout) budterm,
real(dp), dimension(:), intent(in) ddterm,
real(dp), intent(in) cimnew,
real(dp), intent(in) cimold,
real(dp), intent(in) cnew,
integer(i4b), intent(in) idcy )
private

This subroutine calculates and accumulates the immobile domain budget terms into the budterm accumulator

Parameters
[in]ddtermterms comprising the balance equation of the immobile domain
[in]cimnewimmobile domain concenration at the end of this time step
[in]cimoldimmobile domain concentration at end of last time step
[in]cnewmobile domain concentration at the end of this time step
[in]idcyorder of decay rate (0:none, 1:first, 2:zero)

Definition at line 1501 of file gwt-ist.f90.

1502 ! -- modules
1503 ! -- dummy
1504 real(DP), dimension(:, :), intent(inout) :: budterm !<
1505 real(DP), dimension(:), intent(in) :: ddterm !< terms comprising the balance equation of the immobile domain
1506 real(DP), intent(in) :: cimnew !< immobile domain concenration at the end of this time step
1507 real(DP), intent(in) :: cimold !< immobile domain concentration at end of last time step
1508 real(DP), intent(in) :: cnew !< mobile domain concentration at the end of this time step
1509 integer(I4B), intent(in) :: idcy !< order of decay rate (0:none, 1:first, 2:zero)
1510 ! -- local
1511 real(DP) :: rate
1512 integer(I4B) :: i
1513 !
1514 ! -- calculate STORAGE-AQUEOUS
1515 i = 1
1516 rate = -ddterm(1) * cimnew + ddterm(2) * cimold
1517 if (rate > dzero) then
1518 budterm(1, i) = budterm(1, i) + rate
1519 else
1520 budterm(2, i) = budterm(2, i) - rate
1521 end if
1522 !
1523 ! -- calculate STORAGE-SORBED
1524 i = 2
1525 rate = -ddterm(3) * cimnew + ddterm(4) * cimold
1526 if (rate > dzero) then
1527 budterm(1, i) = budterm(1, i) + rate
1528 else
1529 budterm(2, i) = budterm(2, i) - rate
1530 end if
1531 !
1532 ! -- calculate DECAY-AQUEOUS
1533 i = 3
1534 rate = dzero
1535 if (idcy == 1) then
1536 rate = -ddterm(5) * cimnew
1537 else if (idcy == 2) then
1538 rate = -ddterm(7)
1539 else
1540 rate = dzero
1541 end if
1542 if (rate > dzero) then
1543 budterm(1, i) = budterm(1, i) + rate
1544 else
1545 budterm(2, i) = budterm(2, i) - rate
1546 end if
1547 !
1548 ! -- calculate DECAY-SORBED
1549 i = 4
1550 if (idcy == 1) then
1551 rate = -ddterm(6) * cimnew
1552 else if (idcy == 2) then
1553 rate = -ddterm(8)
1554 else
1555 rate = dzero
1556 end if
1557 if (rate > dzero) then
1558 budterm(1, i) = budterm(1, i) + rate
1559 else
1560 budterm(2, i) = budterm(2, i) - rate
1561 end if
1562 !
1563 ! -- calculate MOBILE-DOMAIN
1564 i = 5
1565 rate = ddterm(9) * cnew - ddterm(9) * cimnew
1566 if (rate > dzero) then
1567 budterm(1, i) = budterm(1, i) + rate
1568 else
1569 budterm(2, i) = budterm(2, i) - rate
1570 end if
1571 !
Here is the caller graph for this function:

◆ allocate_scalars()

subroutine gwtistmodule::allocate_scalars ( class(gwtisttype) this)
private

Allocate and initialize package scalars.

Parameters
thisGwtIstType object

Definition at line 832 of file gwt-ist.f90.

833 ! -- modules
834 use memorymanagermodule, only: mem_allocate, mem_setptr
835 use outputcontroldatamodule, only: ocd_cr
836 ! -- dummy
837 class(GwtIstType) :: this !< GwtIstType object
838 ! -- local
839 !
840 ! -- call standard BndType allocate scalars
841 call this%BndType%allocate_scalars()
842 !
843 ! -- Allocate
844 call mem_allocate(this%icimout, 'ICIMOUT', this%memoryPath)
845 call mem_allocate(this%ibudgetout, 'IBUDGETOUT', this%memoryPath)
846 call mem_allocate(this%ibudcsv, 'IBUDCSV', this%memoryPath)
847 call mem_allocate(this%ioutsorbate, 'IOUTSORBATE', this%memoryPath)
848 call mem_allocate(this%isrb, 'ISRB', this%memoryPath)
849 call mem_allocate(this%idcy, 'IDCY', this%memoryPath)
850 call mem_allocate(this%kiter, 'KITER', this%memoryPath)
851 !
852 ! -- Initialize
853 this%lstfmt = ''
854 this%icimout = 0
855 this%ibudgetout = 0
856 this%ibudcsv = 0
857 this%ioutsorbate = 0
858 this%isrb = 0
859 this%idcy = 0
860 this%kiter = 0
861 !
862 ! -- Create the ocd object, which is used to manage printing and saving
863 ! of the immobile domain concentrations
864 call ocd_cr(this%ocd)
memorymanagermodule::mem_allocate
Definition MemoryManager.f90:55
memorymanagermodule::mem_setptr
Definition MemoryManager.f90:93
outputcontroldatamodule
Output control data module.
Definition OutputControlData.f90:2
outputcontroldatamodule::ocd_cr
subroutine, public ocd_cr(ocdobj)
@ brief Create a new output control data type.
Definition OutputControlData.f90:47
Here is the call graph for this function:

◆ get_ddconc()

real(dp) function gwtistmodule::get_ddconc ( real(dp), dimension(:), intent(in) ddterm,
real(dp), intent(in) f,
real(dp), intent(in) cimold,
real(dp), intent(in) cnew )
private

This function calculates the concentration of the immobile domain.

Parameters
[in]ddtermterms comprising the balance equation of the immobile domain
[in]fthe f term used to calculate the immobile domain concentration
[in]cimoldimmobile domain concentration at end of last time step
[in]cnewconcentration of the mobile domain at the end of the time step
Returns
calculated concentration of the immobile domain

Definition at line 1479 of file gwt-ist.f90.

1480 ! -- dummy
1481 real(DP), dimension(:), intent(in) :: ddterm !< terms comprising the balance equation of the immobile domain
1482 real(DP), intent(in) :: f !< the f term used to calculate the immobile domain concentration
1483 real(DP), intent(in) :: cimold !< immobile domain concentration at end of last time step
1484 real(DP), intent(in) :: cnew !< concentration of the mobile domain at the end of the time step
1485 ! -- result
1486 real(DP) :: cimnew !< calculated concentration of the immobile domain
1487 ! -- local
1488 !
1489 ! -- calculate ddconc
1490 cimnew = (ddterm(2) + ddterm(4)) * cimold + ddterm(9) * cnew - ddterm(7) &
1491 - ddterm(8)
1492 cimnew = cimnew / f
Here is the caller graph for this function:

◆ get_ddterm()

subroutine gwtistmodule::get_ddterm ( real(dp), intent(in) thetaim,
real(dp), intent(in) vcell,
real(dp), intent(in) delt,
real(dp), intent(in) swtpdt,
real(dp), intent(in) volfracim,
real(dp), intent(in) rhobim,
real(dp), intent(in) kdnew,
real(dp), intent(in) kdold,
real(dp), intent(in) lambda1im,
real(dp), intent(in) lambda2im,
real(dp), intent(in) gamma1im,
real(dp), intent(in) gamma2im,
real(dp), intent(in) zetaim,
real(dp), dimension(:), intent(inout) ddterm,
real(dp), intent(inout) f )
private

This subroutine calculates the immobile domain (or dual domain) terms. The resulting ddterm and f terms are described in the GWT model report. The terms are concentration coefficients used in the balance equation for the immobile domain.

Parameters
[in]thetaimimmobile domain porosity
[in]vcellvolume of cell
[in]deltlength of time step
[in]swtpdtcell saturation at end of time step
[in]volfracimvolume fraction of immobile domain
[in]rhobimbulk density for the immobile domain (fim * rhob)
[in]kdneweffective distribution coefficient for new time
[in]kdoldeffective distribution coefficient for old time
[in]lambda1imfirst-order decay rate in aqueous phase
[in]lambda2imfirst-order decay rate in sorbed phase
[in]gamma1imzero-order decay rate in aqueous phase
[in]gamma2imzero-order decay rate in sorbed phase
[in]zetaimtransfer coefficient between mobile and immobile domains
[in,out]ddtermnine terms comprising the balance equation of the immobile domain
[in,out]fthe f term used to calculate the immobile domain concentration

Definition at line 1406 of file gwt-ist.f90.

1410 ! -- dummy
1411 real(DP), intent(in) :: thetaim !< immobile domain porosity
1412 real(DP), intent(in) :: vcell !< volume of cell
1413 real(DP), intent(in) :: delt !< length of time step
1414 real(DP), intent(in) :: swtpdt !< cell saturation at end of time step
1415 real(DP), intent(in) :: volfracim !< volume fraction of immobile domain
1416 real(DP), intent(in) :: rhobim !< bulk density for the immobile domain (fim * rhob)
1417 real(DP), intent(in) :: kdnew !< effective distribution coefficient for new time
1418 real(DP), intent(in) :: kdold !< effective distribution coefficient for old time
1419 real(DP), intent(in) :: lambda1im !< first-order decay rate in aqueous phase
1420 real(DP), intent(in) :: lambda2im !< first-order decay rate in sorbed phase
1421 real(DP), intent(in) :: gamma1im !< zero-order decay rate in aqueous phase
1422 real(DP), intent(in) :: gamma2im !< zero-order decay rate in sorbed phase
1423 real(DP), intent(in) :: zetaim !< transfer coefficient between mobile and immobile domains
1424 real(DP), dimension(:), intent(inout) :: ddterm !< nine terms comprising the balance equation of the immobile domain
1425 real(DP), intent(inout) :: f !< the f term used to calculate the immobile domain concentration
1426 ! -- local
1427 real(DP) :: tled
1428 !
1429 ! -- initialize
1430 tled = done / delt
1431 !
1432 ! -- Calculate terms. These terms correspond to the concentration
1433 ! coefficients in equation 7-4 of the GWT model report. However,
1434 ! an updated equation is presented as 9-9 in the supplemental technical
1435 ! information guide (mf6suptechinfo.pdf)
1436 ddterm(1) = thetaim * vcell * tled
1437 ddterm(2) = thetaim * vcell * tled
1438 ddterm(3) = volfracim * rhobim * vcell * kdnew * tled
1439 ddterm(4) = volfracim * rhobim * vcell * kdold * tled
1440 ddterm(5) = thetaim * lambda1im * vcell
1441 ddterm(6) = lambda2im * volfracim * rhobim * kdnew * vcell
1442 ddterm(7) = thetaim * gamma1im * vcell
1443 ddterm(8) = gamma2im * volfracim * rhobim * vcell
1444 ddterm(9) = vcell * swtpdt * zetaim
1445 !
1446 ! -- calculate denominator term, f
1447 f = ddterm(1) + ddterm(3) + ddterm(5) + ddterm(6) + ddterm(9)
Here is the caller graph for this function:

◆ get_freundlich_kd()

real(dp) function gwtistmodule::get_freundlich_kd ( real(dp), intent(in) conc,
real(dp), intent(in) kf,
real(dp), intent(in) a )
private
Parameters
[in]concsolute concentration
[in]kffreundlich constant
[in]afreundlich exponent
Returns
effective distribution coefficient

Definition at line 1576 of file gwt-ist.f90.

1577 ! -- dummy
1578 real(DP), intent(in) :: conc !< solute concentration
1579 real(DP), intent(in) :: kf !< freundlich constant
1580 real(DP), intent(in) :: a !< freundlich exponent
1581 real(DP), parameter :: K = 10.0_dp !< constant to limit derivative at c=0
1582 real(DP) :: eps !< small concentration offset
1583 ! -- return
1584 real(DP) :: kd !< effective distribution coefficient
1585 !
1586 eps = (k / (a * kf))**(1.0_dp / (a - 1.0_dp))
1587 if (conc > dzero) then
1588 kd = kf * (conc + eps)**(a - done) - kf * eps**(a - done)
1589 else
1590 kd = dzero
1591 end if
Here is the caller graph for this function:

◆ get_hcofrhs()

subroutine gwtistmodule::get_hcofrhs ( real(dp), dimension(:), intent(in) ddterm,
real(dp), intent(in) f,
real(dp), intent(in) cimold,
real(dp), intent(inout) hcof,
real(dp), intent(inout) rhs )
private

This subroutine calculates the hcof and rhs terms that must be added to the solution system of equations

Parameters
[in]ddtermterms comprising the balance equation of the immobile domain
[in]fthe f term used to calculate the immobile domain concentration
[in]cimoldimmobile domain concentration at end of last time step
[in,out]hcofcalculated contribution for the a-matrix diagonal position
[in,out]rhscalculated contribution for the solution right-hand side

Definition at line 1456 of file gwt-ist.f90.

1457 ! -- dummy
1458 real(DP), dimension(:), intent(in) :: ddterm !< terms comprising the balance equation of the immobile domain
1459 real(DP), intent(in) :: f !< the f term used to calculate the immobile domain concentration
1460 real(DP), intent(in) :: cimold !< immobile domain concentration at end of last time step
1461 real(DP), intent(inout) :: hcof !< calculated contribution for the a-matrix diagonal position
1462 real(DP), intent(inout) :: rhs !< calculated contribution for the solution right-hand side
1463 !
1464 ! -- calculate hcof
1465 hcof = ddterm(9)**2 / f - ddterm(9)
1466 !
1467 ! -- calculate rhs, and switch the sign because this term needs to
1468 ! be moved to the left hand side
1469 rhs = (ddterm(2) + ddterm(4)) * cimold - ddterm(7) - ddterm(8)
1470 rhs = rhs * ddterm(9) / f
1471 rhs = -rhs
Here is the caller graph for this function:

◆ get_langmuir_kd()

real(dp) function gwtistmodule::get_langmuir_kd ( real(dp), intent(in) conc,
real(dp), intent(in) kl,
real(dp), intent(in) sbar )
private
Parameters
[in]concsolute concentration
[in]kllangmuir constant
[in]sbarlangmuir sorption sites
Returns
effective distribution coefficient

Definition at line 1596 of file gwt-ist.f90.

1597 ! -- dummy
1598 real(DP), intent(in) :: conc !< solute concentration
1599 real(DP), intent(in) :: kl !< langmuir constant
1600 real(DP), intent(in) :: sbar !< langmuir sorption sites
1601 ! -- return
1602 real(DP) :: kd !< effective distribution coefficient
1603 !
1604 if (conc > dzero) then
1605 kd = (kl * sbar) / (done + kl * conc)
1606 else
1607 kd = dzero
1608 end if
Here is the caller graph for this function:

◆ get_thetaim()

real(dp) function gwtistmodule::get_thetaim ( class(gwtisttype) this,
integer(i4b), intent(in) node )
private

Calculate and return thetaim, volume of immobile voids per aquifer volume

Parameters
thisGwtIstType object
[in]nodenode number

Definition at line 1387 of file gwt-ist.f90.

1388 ! -- modules
1389 ! -- dummy
1390 class(GwtIstType) :: this !< GwtIstType object
1391 integer(I4B), intent(in) :: node !< node number
1392 ! -- return
1393 real(DP) :: thetaim
1394 !
1395 thetaim = this%volfrac(node) * this%porosity(node)

◆ ist_ad()

subroutine gwtistmodule::ist_ad ( class(gwtisttype) this)
private

Advance the IST Package and handle the adaptive time stepping feature by copying from new to old or old to new accordingly

Parameters
thisGwtIstType object

Definition at line 246 of file gwt-ist.f90.

247 ! -- modules
248 use simvariablesmodule, only: ifailedstepretry
249 ! -- dummy variables
250 class(GwtIstType) :: this !< GwtIstType object
251 ! -- local variables
252 integer(I4B) :: n
253 !
254 ! -- Call parent advance
255 call this%BndType%bnd_ad()
256 !
257 ! -- set independent kiter counter to zero
258 this%kiter = 0
259 !
260 ! -- copy cimnew into cimold or vice versa if this is a repeat of
261 ! a failed time step
262 if (ifailedstepretry == 0) then
263 do n = 1, this%dis%nodes
264 this%cimold(n) = this%cimnew(n)
265 end do
266 else
267 do n = 1, this%dis%nodes
268 this%cimnew(n) = this%cimold(n)
269 end do
270 end if
simvariablesmodule
This module contains simulation variables.
Definition SimVariables.f90:9
simvariablesmodule::ifailedstepretry
integer(i4b) ifailedstepretry
current retry for this time step
Definition SimVariables.f90:43

◆ ist_allocate_arrays()

subroutine gwtistmodule::ist_allocate_arrays ( class(gwtisttype), intent(inout) this)
private

Allocate and initialize package arrays.

Parameters
[in,out]thisGwtIstType object

Definition at line 872 of file gwt-ist.f90.

873 ! -- modules
874 use memorymanagermodule, only: mem_allocate
875 ! -- dummy
876 class(GwtIstType), intent(inout) :: this !< GwtIstType object
877 ! -- local
878 integer(I4B) :: n
879 !
880 ! -- call standard BndType allocate scalars
881 ! nbound and maxbound are 0 in order to keep memory footprint low
882 call this%BndType%allocate_arrays()
883 !
884 ! -- allocate ist arrays of size nodes
885 call mem_allocate(this%strg, this%dis%nodes, 'STRG', this%memoryPath)
886 call mem_allocate(this%cim, this%dis%nodes, 'CIM', this%memoryPath)
887 call mem_allocate(this%cimnew, this%dis%nodes, 'CIMNEW', this%memoryPath)
888 call mem_allocate(this%cimold, this%dis%nodes, 'CIMOLD', this%memoryPath)
889 call mem_allocate(this%porosity, this%dis%nodes, 'POROSITY', this%memoryPath)
890 call mem_allocate(this%zetaim, this%dis%nodes, 'ZETAIM', this%memoryPath)
891 call mem_allocate(this%volfrac, this%dis%nodes, 'VOLFRAC', this%memoryPath)
892 if (this%isrb == 0) then
893 call mem_allocate(this%bulk_density, 1, 'BULK_DENSITY', this%memoryPath)
894 call mem_allocate(this%distcoef, 1, 'DISTCOEF', this%memoryPath)
895 call mem_allocate(this%sp2, 1, 'SP2', this%memoryPath)
896 call mem_allocate(this%cimsrb, 1, 'SORBATE', this%memoryPath)
897 else
898 call mem_allocate(this%bulk_density, this%dis%nodes, 'BULK_DENSITY', &
899 this%memoryPath)
900 call mem_allocate(this%distcoef, this%dis%nodes, 'DISTCOEF', &
901 this%memoryPath)
902 call mem_allocate(this%cimsrb, this%dis%nodes, 'SORBATE', &
903 this%memoryPath)
904 if (this%isrb == 1) then
905 call mem_allocate(this%sp2, 1, 'SP2', this%memoryPath)
906 else
907 call mem_allocate(this%sp2, this%dis%nodes, 'SP2', this%memoryPath)
908 end if
909 end if
910 if (this%idcy == 0) then
911 call mem_allocate(this%decay, 1, 'DECAY', this%memoryPath)
912 call mem_allocate(this%decaylast, 1, 'DECAYLAST', this%memoryPath)
913 else
914 call mem_allocate(this%decay, this%dis%nodes, 'DECAY', this%memoryPath)
915 call mem_allocate(this%decaylast, this%dis%nodes, 'DECAYLAST', &
916 this%memoryPath)
917 end if
918 if (this%isrb == 0 .and. this%idcy == 0) then
919 call mem_allocate(this%decayslast, 1, 'DECAYSLAST', this%memoryPath)
920 else
921 call mem_allocate(this%decayslast, this%dis%nodes, 'DECAYSLAST', &
922 this%memoryPath)
923 end if
924 call mem_allocate(this%decay_sorbed, 1, 'DECAY_SORBED', this%memoryPath)
925 !
926 ! -- initialize
927 do n = 1, this%dis%nodes
928 this%strg(n) = dzero
929 this%cim(n) = dzero
930 this%cimnew(n) = dzero
931 this%cimold(n) = dzero
932 this%porosity(n) = dzero
933 this%zetaim(n) = dzero
934 this%volfrac(n) = dzero
935 end do
936 do n = 1, size(this%bulk_density)
937 this%bulk_density(n) = dzero
938 this%distcoef(n) = dzero
939 this%cimsrb(n) = dzero
940 end do
941 do n = 1, size(this%sp2)
942 this%sp2(n) = dzero
943 end do
944 do n = 1, size(this%decay)
945 this%decay(n) = dzero
946 this%decaylast(n) = dzero
947 end do
948 do n = 1, size(this%decayslast)
949 this%decayslast(n) = dzero
950 end do
951 !
952 ! -- Set pointers
953 this%ocd%dis => this%dis

◆ ist_ar()

subroutine gwtistmodule::ist_ar ( class(gwtisttype), intent(inout) this)
private

Method to allocate and read static data for the package.

Parameters
[in,out]thisGwtIstType object

Definition at line 167 of file gwt-ist.f90.

168 ! -- modules
169 use simmodule, only: store_error, count_errors, store_error_filename
170 use budgetmodule, only: budget_cr
171 ! -- dummy
172 class(GwtIstType), intent(inout) :: this !< GwtIstType object
173 ! -- local
174 integer(I4B) :: n
175 !
176 ! -- Allocate arrays
177 call this%ist_allocate_arrays()
178 !
179 ! -- Now that arrays are allocated, check in the cimnew array to
180 ! the output control manager for subsequent printing/saving
181 call this%ocd%init_dbl('CIM', this%cimnew, this%dis, 'PRINT LAST ', &
182 'COLUMNS 10 WIDTH 11 DIGITS 4 GENERAL ', &
183 this%iout, dhnoflo)
184
185 ! -- apply user override if provided
186 if (this%lstfmt /= '') then
187 call this%ocd%set_prnfmt(trim(this%lstfmt)//" ", 0)
188 end if
189 !
190 ! -- source the data block
191 call this%source_data()
192 !
193 ! -- set cimnew to the cim start values read from input
194 do n = 1, this%dis%nodes
195 this%cimnew(n) = this%cim(n)
196 end do
197 !
198 ! -- add volfrac to the volfracim accumulator in mst package
199 call this%mst%addto_volfracim(this%volfrac)
200 !
201 ! -- setup the immobile domain budget
202 call budget_cr(this%budget, this%memoryPath)
203 call this%budget%budget_df(nbditems, 'MASS', 'M', bdzone=this%packName)
204 call this%budget%set_ibudcsv(this%ibudcsv)
205 !
206 ! -- Create isotherm object if sorption is active
207 this%isotherm => create_isotherm(this%isrb, this%distcoef, this%sp2)
208 !
209 ! -- Perform a check to ensure that sorption and decay are set
210 ! consistently between the MST and IST packages.
211 if (this%idcy /= this%mst%idcy) then
212 call store_error('DECAY must be activated consistently between the &
213 &MST and IST Packages. Activate or deactivate DECAY for both &
214 &Packages.')
215 end if
216 if (this%isrb /= this%mst%isrb) then
217 call store_error('SORPTION must be activated consistently between the &
218 &MST and IST Packages. Activate or deactivate SORPTION for both &
219 &Packages. If activated, the same type of sorption (LINEAR, &
220 &FREUNDLICH, or LANGMUIR) must be specified in the options block of &
221 &both the MST and IST Packages.')
222 end if
223 if (count_errors() > 0) then
224 call store_error_filename(this%input_fname)
225 end if
budgetmodule
This module contains the BudgetModule.
Definition Budget.f90:20
budgetmodule::budget_cr
subroutine, public budget_cr(this, name_model)
@ brief Create a new budget object
Definition Budget.f90:85
simmodule
This module contains simulation methods.
Definition Sim.f90:10
simmodule::store_error
subroutine, public store_error(msg, terminate)
Store an error message.
Definition Sim.f90:92
simmodule::count_errors
integer(i4b) function, public count_errors()
Return number of errors.
Definition Sim.f90:59
simmodule::store_error_filename
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
Definition Sim.f90:204
Here is the call graph for this function:

◆ ist_bd()

subroutine gwtistmodule::ist_bd ( class(gwtisttype) this,
type(budgettype), intent(inout) model_budget )
private

Add the flow between IST package and the model (ratin and ratout) to the model budget.

Parameters
thisGwtIstType object
[in,out]model_budgetmodel budget object

Definition at line 573 of file gwt-ist.f90.

574 ! -- modules
575 use tdismodule, only: delt
576 use budgetmodule, only: budgettype, rate_accumulator
577 ! -- dummy
578 class(GwtIstType) :: this !< GwtIstType object
579 type(BudgetType), intent(inout) :: model_budget !< model budget object
580 ! -- local
581 real(DP) :: ratin
582 real(DP) :: ratout
583 integer(I4B) :: isuppress_output
584 isuppress_output = 0
585 call rate_accumulator(this%strg(:), ratin, ratout)
586 call model_budget%addentry(ratin, ratout, delt, this%text, &
587 isuppress_output, this%packName)
budgetmodule::rate_accumulator
subroutine, public rate_accumulator(flow, rin, rout)
@ brief Rate accumulator subroutine
Definition Budget.f90:634
tdismodule::delt
real(dp), pointer, public delt
length of the current time step
Definition tdis.f90:32
budgetmodule::budgettype
Derived type for the Budget object.
Definition Budget.f90:40
Here is the call graph for this function:

◆ ist_calc_csrb()

subroutine gwtistmodule::ist_calc_csrb ( class(gwtisttype) this,
real(dp), dimension(:), intent(in) cim )
private
Parameters
thisGwtMstType object
[in]cimimmobile domain aqueous concentration at end of this time step

Definition at line 548 of file gwt-ist.f90.

549 ! -- dummy
550 class(GwtIstType) :: this !< GwtMstType object
551 real(DP), intent(in), dimension(:) :: cim !< immobile domain aqueous concentration at end of this time step
552 ! -- local
553 integer(I4B) :: n
554 real(DP) :: csrb
555
556 ! Calculate sorbed concentration
557 do n = 1, size(cim)
558 csrb = dzero
559 if (this%ibound(n) > 0 .and. this%isrb /= sorption_off) then
560 csrb = this%isotherm%value(cim, n)
561 end if
562 this%cimsrb(n) = csrb
563 end do
564

◆ ist_cq()

subroutine gwtistmodule::ist_cq ( class(gwtisttype), intent(inout) this,
real(dp), dimension(:), intent(in) x,
real(dp), dimension(:), intent(inout), contiguous flowja,
integer(i4b), intent(in), optional iadv )
private
Parameters
[in,out]thisGwtIstType object
[in]xcurrent dependent-variable value
[in,out]flowjaflow between two connected control volumes
[in]iadvflag that indicates if this is an advance package

Definition at line 402 of file gwt-ist.f90.

403 ! modules
404 use tdismodule, only: delt
405 use constantsmodule, only: dzero
406 ! dummy
407 class(GwtIstType), intent(inout) :: this !< GwtIstType object
408 real(DP), dimension(:), intent(in) :: x !< current dependent-variable value
409 real(DP), dimension(:), contiguous, intent(inout) :: flowja !< flow between two connected control volumes
410 integer(I4B), optional, intent(in) :: iadv !< flag that indicates if this is an advance package
411 ! local
412 integer(I4B) :: idiag
413 integer(I4B) :: n
414 real(DP) :: rate
415 real(DP) :: swt, swtpdt
416 real(DP) :: hhcof, rrhs
417 real(DP) :: vcell
418 real(DP) :: thetaim
419 real(DP) :: zetaim
420 real(DP) :: kdnew
421 real(DP) :: kdold
422 real(DP) :: volfracim
423 real(DP) :: rhobim
424 real(DP) :: lambda1im
425 real(DP) :: lambda2im
426 real(DP) :: gamma1im
427 real(DP) :: gamma2im
428 real(DP) :: cimnew
429 real(DP) :: cimold
430 real(DP) :: f
431 real(DP) :: cimsrbold
432 real(DP) :: cimsrbnew
433 real(DP), dimension(9) :: ddterm
434 ! formats
435
436 ! initialize
437 this%budterm(:, :) = dzero
438
439 ! Calculate immobile domain transfer rate
440 do n = 1, this%dis%nodes
441
442 ! skip if transport inactive
443 rate = dzero
444 cimnew = dzero
445 if (this%ibound(n) > 0) then
446
447 ! calculate new and old water volumes
448 vcell = this%dis%area(n) * (this%dis%top(n) - this%dis%bot(n))
449 swtpdt = this%fmi%gwfsat(n)
450 swt = this%fmi%gwfsatold(n, delt)
451 thetaim = this%get_thetaim(n)
452
453 ! set exchange coefficient
454 zetaim = this%zetaim(n)
455
456 ! Calculate exchange with immobile domain
457 rate = dzero
458 hhcof = dzero
459 rrhs = dzero
460 kdnew = dzero
461 kdold = dzero
462 volfracim = dzero
463 rhobim = dzero
464 lambda1im = dzero
465 lambda2im = dzero
466 gamma1im = dzero
467 gamma2im = dzero
468 if (this%idcy == 1) lambda1im = this%decay(n)
469 if (this%idcy == 2) then
470 gamma1im = get_zero_order_decay(this%decay(n), this%decaylast(n), 0, &
471 this%cimold(n), this%cimnew(n), delt)
472 end if
473
474 ! setup sorption variables
475 if (this%isrb > 0) then
476
477 ! initialize sorption variables
478 volfracim = this%volfrac(n)
479 rhobim = this%bulk_density(n)
480
481 ! set isotherm dependent sorption variables
482 cimsrbnew = this%isotherm%value(this%cimnew, n)
483 cimsrbold = this%isotherm%value(this%cimold, n)
484 select case (this%isrb)
485 case (sorption_linear)
486 kdnew = this%distcoef(n)
487 kdold = this%distcoef(n)
488 case (sorption_freund)
489 kdnew = get_freundlich_kd(this%cimnew(n), this%distcoef(n), &
490 this%sp2(n))
491 kdold = get_freundlich_kd(this%cimold(n), this%distcoef(n), &
492 this%sp2(n))
493 case (sorption_lang)
494 kdnew = get_langmuir_kd(this%cimnew(n), this%distcoef(n), &
495 this%sp2(n))
496 kdold = get_langmuir_kd(this%cimold(n), this%distcoef(n), &
497 this%sp2(n))
498 end select
499
500 ! set decay of sorbed solute parameters
501 if (this%idcy == 1) then
502 lambda2im = this%decay_sorbed(n)
503 else if (this%idcy == 2) then
504 gamma2im = get_zero_order_decay(this%decay_sorbed(n), &
505 this%decayslast(n), &
506 0, cimsrbold, &
507 cimsrbnew, delt)
508 end if
509 end if
510
511 ! calculate the terms and then get the hcof and rhs contributions
512 call get_ddterm(thetaim, vcell, delt, swtpdt, &
513 volfracim, rhobim, kdnew, kdold, lambda1im, lambda2im, &
514 gamma1im, gamma2im, zetaim, ddterm, f)
515 cimold = this%cimold(n)
516 call get_hcofrhs(ddterm, f, cimold, hhcof, rrhs)
517
518 ! calculate rate from hcof and rhs
519 rate = hhcof * x(n) - rrhs
520
521 ! calculate immobile domain concentration
522 cimnew = get_ddconc(ddterm, f, cimold, x(n))
523
524 ! accumulate the budget terms
525 call accumulate_budterm(this%budterm, ddterm, cimnew, cimold, x(n), &
526 this%idcy)
527 end if
528
529 ! store rate and add to flowja
530 this%strg(n) = rate
531 idiag = this%dis%con%ia(n)
532 flowja(idiag) = flowja(idiag) + rate
533
534 ! store immobile domain concentration
535 this%cimnew(n) = cimnew
536
537 end do
538
539 ! calculate csrb
540 if (this%isrb /= 0) then
541 call this%ist_calc_csrb(this%cimnew)
542 end if
543
constantsmodule
This module contains simulation constants.
Definition Constants.f90:9
constantsmodule::dzero
real(dp), parameter dzero
real constant zero
Definition Constants.f90:65
Here is the call graph for this function:

◆ ist_create()

subroutine, public gwtistmodule::ist_create ( class(bndtype), pointer packobj,
integer(i4b), intent(in) id,
integer(i4b), intent(in) ibcnum,
integer(i4b), intent(in) inunit,
integer(i4b), intent(in) iout,
character(len=*), intent(in) namemodel,
character(len=*), intent(in) pakname,
character(len=*), intent(in) mempath,
type(tspfmitype), intent(in), pointer fmi,
type(gwtmsttype), intent(in), pointer mst )

Create a new IST object

Parameters
packobjBndType pointer that will point to new IST Package
[in]idname of the model
[in]ibcnumconsecutive package number
[in]inunitunit number of package input file
[in]ioutunit number of model listing file
[in]namemodelname of the model
[in]paknamename of the package

Definition at line 119 of file gwt-ist.f90.

121 ! -- dummy
122 class(BndType), pointer :: packobj !< BndType pointer that will point to new IST Package
123 integer(I4B), intent(in) :: id !< name of the model
124 integer(I4B), intent(in) :: ibcnum !< consecutive package number
125 integer(I4B), intent(in) :: inunit !< unit number of package input file
126 integer(I4B), intent(in) :: iout !< unit number of model listing file
127 character(len=*), intent(in) :: namemodel !< name of the model
128 character(len=*), intent(in) :: pakname !< name of the package
129 character(len=*), intent(in) :: mempath
130 type(TspFmiType), pointer, intent(in) :: fmi
131 type(GwtMstType), pointer, intent(in) :: mst
132 ! -- local
133 type(GwtIstType), pointer :: istobj
134 !
135 ! -- allocate the object and assign values to object variables
136 allocate (istobj)
137 packobj => istobj
138 !
139 ! -- create name and memory path
140 call packobj%set_names(ibcnum, namemodel, pakname, ftype, mempath)
141 packobj%text = text
142 !
143 ! -- allocate scalars
144 call packobj%allocate_scalars()
145 !
146 ! -- initialize package
147 call packobj%pack_initialize()
148 !
149 ! -- store values
150 packobj%inunit = inunit
151 packobj%iout = iout
152 packobj%id = id
153 packobj%ibcnum = ibcnum
154 packobj%ncolbnd = 1
155 packobj%iscloc = 1
156 !
157 ! -- Point IST specific variables
158 istobj%fmi => fmi
159 istobj%mst => mst
Here is the caller graph for this function:

◆ ist_da()

subroutine gwtistmodule::ist_da ( class(gwtisttype) this)
private

Deallocate package scalars and arrays.

Parameters
thisGwtIstType object

Definition at line 777 of file gwt-ist.f90.

778 ! -- modules
779 use memorymanagermodule, only: mem_deallocate
780 ! -- dummy
781 class(GwtIstType) :: this !< GwtIstType object
782 !
783 ! -- Deallocate arrays if package was active
784 if (this%inunit > 0) then
785 call mem_deallocate(this%icimout)
786 call mem_deallocate(this%ibudgetout)
787 call mem_deallocate(this%ibudcsv)
788 call mem_deallocate(this%ioutsorbate)
789 call mem_deallocate(this%idcy)
790 call mem_deallocate(this%isrb)
791 call mem_deallocate(this%kiter)
792 call mem_deallocate(this%cim)
793 call mem_deallocate(this%cimnew)
794 call mem_deallocate(this%cimold)
795 call mem_deallocate(this%cimsrb)
796 call mem_deallocate(this%zetaim)
797 call mem_deallocate(this%porosity)
798 call mem_deallocate(this%volfrac)
799 call mem_deallocate(this%bulk_density)
800 call mem_deallocate(this%distcoef)
801 call mem_deallocate(this%sp2)
802 call mem_deallocate(this%decay)
803 call mem_deallocate(this%decaylast)
804 call mem_deallocate(this%decayslast)
805 call mem_deallocate(this%decay_sorbed)
806 call mem_deallocate(this%strg)
807 this%fmi => null()
808 this%mst => null()
809 end if
810 !
811 ! -- Scalars
812 !
813 ! -- Objects
814 call this%budget%budget_da()
815 deallocate (this%budget)
816 call this%ocd%ocd_da()
817 deallocate (this%ocd)
818 if (associated(this%isotherm)) then
819 deallocate (this%isotherm)
820 nullify (this%isotherm)
821 end if
822 !
823 ! -- deallocate parent
824 call this%BndType%bnd_da()
memorymanagermodule::mem_deallocate
Definition MemoryManager.f90:129

◆ ist_fc()

subroutine gwtistmodule::ist_fc ( class(gwtisttype) this,
real(dp), dimension(:), intent(inout) rhs,
integer(i4b), dimension(:), intent(in) ia,
integer(i4b), dimension(:), intent(in) idxglo,
class(matrixbasetype), pointer matrix_sln )
private
Parameters
thisGwtIstType object
[in,out]rhsright-hand side vector for model
[in]iasolution CRS row pointers
[in]idxglomapping vector for model (local) to solution (global)
matrix_slnsolution coefficient matrix

Definition at line 275 of file gwt-ist.f90.

276 ! modules
277 use tdismodule, only: delt
278 ! dummy
279 class(GwtIstType) :: this !< GwtIstType object
280 real(DP), dimension(:), intent(inout) :: rhs !< right-hand side vector for model
281 integer(I4B), dimension(:), intent(in) :: ia !< solution CRS row pointers
282 integer(I4B), dimension(:), intent(in) :: idxglo !< mapping vector for model (local) to solution (global)
283 class(MatrixBaseType), pointer :: matrix_sln !< solution coefficient matrix
284 ! local
285 integer(I4B) :: n, idiag
286 real(DP) :: tled
287 real(DP) :: hhcof, rrhs
288 real(DP) :: swt, swtpdt
289 real(DP) :: vcell
290 real(DP) :: thetaim
291 real(DP) :: zetaim
292 real(DP) :: kdnew
293 real(DP) :: kdold
294 real(DP) :: volfracim
295 real(DP) :: rhobim
296 real(DP) :: lambda1im
297 real(DP) :: lambda2im
298 real(DP) :: gamma1im
299 real(DP) :: gamma2im
300 real(DP) :: cimold
301 real(DP) :: f
302 real(DP) :: cimsrbold
303 real(DP) :: cimsrbnew
304 real(DP), dimension(9) :: ddterm
305
306 ! set variables
307 tled = done / delt
308 this%kiter = this%kiter + 1
309
310 ! loop through each node and calculate immobile domain contribution
311 ! to hcof and rhs
312 do n = 1, this%dis%nodes
313
314 ! skip if transport inactive
315 if (this%ibound(n) <= 0) cycle
316
317 ! calculate new and old water volumes
318 vcell = this%dis%area(n) * (this%dis%top(n) - this%dis%bot(n))
319 swtpdt = this%fmi%gwfsat(n)
320 swt = this%fmi%gwfsatold(n, delt)
321 thetaim = this%get_thetaim(n)
322 idiag = ia(n)
323
324 ! set exchange coefficient
325 zetaim = this%zetaim(n)
326
327 ! Add dual domain mass transfer contributions to rhs and hcof
328 ! dcimsrbdc = DZERO
329 kdnew = dzero
330 kdold = dzero
331 volfracim = dzero
332 rhobim = dzero
333 lambda1im = dzero
334 lambda2im = dzero
335 gamma1im = dzero
336 gamma2im = dzero
337
338 ! set variables for decay of aqueous solute
339 if (this%idcy == 1) lambda1im = this%decay(n)
340 if (this%idcy == 2) then
341 gamma1im = get_zero_order_decay(this%decay(n), this%decaylast(n), &
342 this%kiter, this%cimold(n), &
343 this%cimnew(n), delt)
344 this%decaylast(n) = gamma1im
345 end if
346
347 ! setup sorption variables
348 if (this%isrb > 0) then
349
350 ! initialize sorption variables
351 volfracim = this%volfrac(n)
352 rhobim = this%bulk_density(n)
353
354 ! set isotherm dependent sorption variables
355 cimsrbnew = this%isotherm%value(this%cimnew, n)
356 cimsrbold = this%isotherm%value(this%cimold, n)
357 select case (this%isrb)
358 case (sorption_linear)
359 kdnew = this%distcoef(n)
360 kdold = this%distcoef(n)
361 case (sorption_freund)
362 kdnew = get_freundlich_kd(this%cimnew(n), this%distcoef(n), &
363 this%sp2(n))
364 kdold = get_freundlich_kd(this%cimold(n), this%distcoef(n), &
365 this%sp2(n))
366 case (sorption_lang)
367 kdnew = get_langmuir_kd(this%cimnew(n), this%distcoef(n), &
368 this%sp2(n))
369 kdold = get_langmuir_kd(this%cimold(n), this%distcoef(n), &
370 this%sp2(n))
371 end select
372
373 ! set decay of sorbed solute parameters
374 if (this%idcy == 1) then
375 lambda2im = this%decay_sorbed(n)
376 else if (this%idcy == 2) then
377 gamma2im = get_zero_order_decay(this%decay_sorbed(n), &
378 this%decayslast(n), &
379 this%kiter, cimsrbold, &
380 cimsrbnew, delt)
381 this%decayslast(n) = gamma2im
382 end if
383 end if
384
385 ! calculate dual domain terms and get the hcof and rhs contributions
386 call get_ddterm(thetaim, vcell, delt, swtpdt, &
387 volfracim, rhobim, kdnew, kdold, lambda1im, lambda2im, &
388 gamma1im, gamma2im, zetaim, ddterm, f)
389 cimold = this%cimold(n)
390 call get_hcofrhs(ddterm, f, cimold, hhcof, rrhs)
391
392 ! update solution accumulators
393 call matrix_sln%add_value_pos(idxglo(idiag), hhcof)
394 rhs(n) = rhs(n) + rrhs
395
396 end do
397
Here is the call graph for this function:

◆ ist_ot_bdsummary()

subroutine gwtistmodule::ist_ot_bdsummary ( class(gwtisttype) this,
integer(i4b), intent(in) kstp,
integer(i4b), intent(in) kper,
integer(i4b), intent(in) iout,
integer(i4b), intent(in) ibudfl )
private

Output advanced boundary package budget summary. This method only needs to be overridden for advanced packages that save budget summaries to the model listing file.

Parameters
thisGwtIstType object
[in]kstptime step number
[in]kperperiod number
[in]ioutflag and unit number for the model listing file
[in]ibudflflag indicating budget should be written

Definition at line 746 of file gwt-ist.f90.

747 ! -- modules
748 use tdismodule, only: delt, totim
749 ! -- dummy variables
750 class(GwtIstType) :: this !< GwtIstType object
751 integer(I4B), intent(in) :: kstp !< time step number
752 integer(I4B), intent(in) :: kper !< period number
753 integer(I4B), intent(in) :: iout !< flag and unit number for the model listing file
754 integer(I4B), intent(in) :: ibudfl !< flag indicating budget should be written
755 ! -- local
756 integer(I4B) :: isuppress_output = 0
757 !
758 ! -- Fill budget terms
759 call this%budget%reset()
760 call this%budget%addentry(this%budterm, delt, budtxt, isuppress_output)
761 !
762 ! -- Write budget to list file
763 call this%budget%finalize_step(delt)
764 if (ibudfl /= 0) then
765 call this%budget%budget_ot(kstp, kper, iout)
766 end if
767 !
768 ! -- Write budget csv
769 call this%budget%writecsv(totim)
tdismodule::totim
real(dp), pointer, public totim
time relative to start of simulation
Definition tdis.f90:35
tdismodule::kstp
integer(i4b), pointer, public kstp
current time step number
Definition tdis.f90:27
tdismodule::kper
integer(i4b), pointer, public kper
current stress period number
Definition tdis.f90:26

◆ ist_ot_dv()

subroutine gwtistmodule::ist_ot_dv ( class(gwtisttype) this,
integer(i4b), intent(in) idvsave,
integer(i4b), intent(in) idvprint )
private
Parameters
thisBndType object
[in]idvsaveflag and unit number for dependent-variable output
[in]idvprintflag indicating if dependent-variable should be written to the model listing file

Definition at line 659 of file gwt-ist.f90.

660 ! dummy variables
661 class(GwtIstType) :: this !< BndType object
662 integer(I4B), intent(in) :: idvsave !< flag and unit number for dependent-variable output
663 integer(I4B), intent(in) :: idvprint !< flag indicating if dependent-variable should be written to the model listing file
664
665 call this%output_immobile_concentration(idvsave, idvprint)
666 call this%output_immobile_sorbate_concentration(idvsave, idvprint)
667

◆ ist_ot_model_flows()

subroutine gwtistmodule::ist_ot_model_flows ( class(gwtisttype) this,
integer(i4b), intent(in) icbcfl,
integer(i4b), intent(in) ibudfl,
integer(i4b), intent(in) icbcun,
integer(i4b), dimension(:), intent(in), optional imap )
private

Output flow terms between the IST package and model to a binary file and/or print flows to the model listing file.

Parameters
thisGwtIstType object
[in]icbcflflag for cell-by-cell output
[in]ibudflflag indication if cell-by-cell data should be saved
[in]icbcununit number for cell-by-cell output
[in]imapmapping vector

Definition at line 596 of file gwt-ist.f90.

597 ! -- modules
598 use constantsmodule, only: dzero
599 ! -- dummy
600 class(GwtIstType) :: this !< GwtIstType object
601 integer(I4B), intent(in) :: icbcfl !< flag for cell-by-cell output
602 integer(I4B), intent(in) :: ibudfl !< flag indication if cell-by-cell data should be saved
603 integer(I4B), intent(in) :: icbcun !< unit number for cell-by-cell output
604 integer(I4B), dimension(:), optional, intent(in) :: imap !< mapping vector
605 ! -- local
606 integer(I4B) :: n
607 integer(I4B) :: ibinun
608 integer(I4B) :: nbound
609 integer(I4B) :: naux
610 real(DP) :: rate
611 real(DP), dimension(0) :: auxrow
612 !
613 ! -- Set unit number for binary output
614 if (this%ipakcb < 0) then
615 ibinun = icbcun
616 elseif (this%ipakcb == 0) then
617 ibinun = 0
618 else
619 ibinun = this%ipakcb
620 end if
621 if (icbcfl == 0) ibinun = 0
622 !
623 ! -- Record the storage rate if requested
624 !
625 ! -- If cell-by-cell flows will be saved as a list, write header.
626 if (ibinun /= 0) then
627 nbound = this%dis%nodes
628 naux = 0
629 call this%dis%record_srcdst_list_header(this%text, this%name_model, &
630 this%name_model, this%name_model, &
631 this%packName, naux, this%auxname, &
632 ibinun, nbound, this%iout)
633 end if
634 !
635 ! -- Calculate immobile domain rhs and hcof
636 do n = 1, this%dis%nodes
637 !
638 ! -- skip if transport inactive
639 rate = dzero
640 if (this%ibound(n) > 0) then
641 !
642 ! -- set rate from this%strg
643 rate = this%strg(n)
644 end if
645 !
646 ! -- If saving cell-by-cell flows in list, write flow
647 if (ibinun /= 0) then
648 call this%dis%record_mf6_list_entry(ibinun, n, n, rate, &
649 naux, auxrow, &
650 olconv=.true., &
651 olconv2=.true.)
652 end if
653 !
654 end do

◆ ist_read_dimensions()

subroutine gwtistmodule::ist_read_dimensions ( class(gwtisttype), intent(inout) this)
private

Read dimensions for package.

Parameters
[in,out]thisGwtIstType object

Definition at line 1375 of file gwt-ist.f90.

1376 ! -- dummy
1377 class(GwtIstType), intent(inout) :: this !< GwtIstType object
1378 ! -- local
1379 ! -- format

◆ ist_rp()

subroutine gwtistmodule::ist_rp ( class(gwtisttype), intent(inout) this)
private

Method to read and prepare package data

Parameters
[in,out]thisGwtIstType object

Definition at line 233 of file gwt-ist.f90.

234 ! -- dummy
235 class(GwtIstType), intent(inout) :: this !< GwtIstType object
236 ! -- local
237 ! -- format

◆ log_data()

subroutine gwtistmodule::log_data ( class(gwtisttype), intent(inout) this,
type(gwtistparamfoundtype), intent(in) found )
private

Definition at line 1334 of file gwt-ist.f90.

1335 use gwtistinputmodule, only: gwtistparamfoundtype
1336 class(GwTIstType), intent(inout) :: this
1337 type(GwtIstParamFoundType), intent(in) :: found
1338
1339 write (this%iout, '(1x,a)') 'PROCESSING GRIDDATA'
1340 if (found%porosity) then
1341 write (this%iout, '(4x,a)') 'MOBILE DOMAIN POROSITY set from input file'
1342 end if
1343 if (found%volfrac) then
1344 write (this%iout, '(4x,a)') 'VOLFRAC set from input file'
1345 end if
1346 if (found%zetaim) then
1347 write (this%iout, '(4x,a)') 'ZETAIM set from input file'
1348 end if
1349 if (found%cim) then
1350 write (this%iout, '(4x,a)') 'CIM set from input file'
1351 end if
1352 if (found%decay) then
1353 write (this%iout, '(4x,a)') 'DECAY RATE set from input file'
1354 end if
1355 if (found%decay_sorbed) then
1356 write (this%iout, '(4x,a)') 'DECAY SORBED RATE set from input file'
1357 end if
1358 if (found%bulk_density) then
1359 write (this%iout, '(4x,a)') 'BULK DENSITY set from input file'
1360 end if
1361 if (found%distcoef) then
1362 write (this%iout, '(4x,a)') 'DISTRIBUTION COEFFICIENT set from input file'
1363 end if
1364 if (found%sp2) then
1365 write (this%iout, '(4x,a)') 'SECOND SORPTION PARAM set from input file'
1366 end if
1367 write (this%iout, '(1x,a)') 'END PROCESSING GRIDDATA'

◆ log_options()

subroutine gwtistmodule::log_options ( class(gwtisttype), intent(inout) this,
type(gwtistparamfoundtype), intent(in) found,
character(len=*), intent(in) cim6_fname,
character(len=*), intent(in) budget_fname,
character(len=*), intent(in) budgetcsv_fname,
character(len=*), intent(in) sorbate_fname )
private

Definition at line 1064 of file gwt-ist.f90.

1066 use gwtistinputmodule, only: gwtistparamfoundtype
1067 class(GwTIstType), intent(inout) :: this
1068 type(GwtIstParamFoundType), intent(in) :: found
1069 character(len=*), intent(in) :: cim6_fname
1070 character(len=*), intent(in) :: budget_fname
1071 character(len=*), intent(in) :: budgetcsv_fname
1072 character(len=*), intent(in) :: sorbate_fname
1073 ! -- formats
1074 character(len=*), parameter :: fmtisvflow = &
1075 "(4x,'CELL-BY-CELL FLOW INFORMATION WILL BE SAVED TO BINARY FILE &
1076 &WHENEVER ICBCFL IS NOT ZERO.')"
1077 character(len=*), parameter :: fmtlinear = &
1078 &"(4x,'LINEAR SORPTION IS SELECTED. ')"
1079 character(len=*), parameter :: fmtfreundlich = &
1080 &"(4x,'FREUNDLICH SORPTION IS ACTIVE. ')"
1081 character(len=*), parameter :: fmtlangmuir = &
1082 &"(4x,'LANGMUIR SORPTION IS ACTIVE. ')"
1083 character(len=*), parameter :: fmtidcy1 = &
1084 &"(4x,'FIRST-ORDER DECAY IS ACTIVE. ')"
1085 character(len=*), parameter :: fmtidcy2 = &
1086 &"(4x,'ZERO-ORDER DECAY IS ACTIVE. ')"
1087 character(len=*), parameter :: fmtistbin = &
1088 "(4x, 'IST ', 1x, a, 1x, ' WILL BE SAVED TO FILE: ', a, &
1089 &/4x, 'OPENED ON UNIT: ', I0)"
1090
1091 write (this%iout, '(1x,a)') 'PROCESSING IMMOBILE STORAGE AND TRANSFER &
1092 &OPTIONS'
1093 if (found%save_flows) then
1094 write (this%iout, fmtisvflow)
1095 end if
1096 if (found%budgetfile) then
1097 write (this%iout, fmtistbin) 'BUDGET', trim(adjustl(budget_fname)), &
1098 this%ibudgetout
1099 end if
1100 if (found%budgetcsvfile) then
1101 write (this%iout, fmtistbin) 'BUDGET CSV', trim(adjustl(budgetcsv_fname)), &
1102 this%ibudcsv
1103 end if
1104 if (found%sorption) then
1105 select case (this%isrb)
1106 case (sorption_linear)
1107 write (this%iout, fmtlinear)
1108 case (sorption_freund)
1109 write (this%iout, fmtfreundlich)
1110 case (sorption_lang)
1111 write (this%iout, fmtlangmuir)
1112 end select
1113 end if
1114 if (found%order1_decay) then
1115 write (this%iout, fmtidcy1)
1116 end if
1117 if (found%order0_decay) then
1118 write (this%iout, fmtidcy2)
1119 end if
1120 if (found%sorbatefile) then
1121 write (this%iout, fmtistbin) &
1122 'SORBATE', sorbate_fname, this%ioutsorbate
1123 end if
1124 write (this%iout, '(1x,a)') 'END OF IMMOBILE STORAGE AND TRANSFER &
1125 &OPTIONS'

◆ output_immobile_concentration()

subroutine gwtistmodule::output_immobile_concentration ( class(gwtisttype) this,
integer(i4b), intent(in) idvsave,
integer(i4b), intent(in) idvprint )
private
Parameters
thisBndType object
[in]idvsaveflag and unit number for dependent-variable output
[in]idvprintflag indicating if dependent-variable should be written to the model listing file

Definition at line 672 of file gwt-ist.f90.

673 ! modules
674 use tdismodule, only: kstp, endofperiod
675 ! dummy variables
676 class(GwtIstType) :: this !< BndType object
677 integer(I4B), intent(in) :: idvsave !< flag and unit number for dependent-variable output
678 integer(I4B), intent(in) :: idvprint !< flag indicating if dependent-variable should be written to the model listing file
679 ! local
680 integer(I4B) :: ipflg
681 integer(I4B) :: ibinun
682 !
683 ! Save cim to a binary file. ibinun is a flag where 1 indicates that
684 ! cim should be written to a binary file if a binary file is open for it.
685 ipflg = 0
686 ibinun = 1
687 if (idvsave == 0) ibinun = 0
688 if (ibinun /= 0) then
689 call this%ocd%ocd_ot(ipflg, kstp, endofperiod, this%iout, &
690 iprint_opt=0, isav_opt=ibinun)
691 end if
692 !
693 ! Print immobile domain concentrations to listing file
694 if (idvprint /= 0) then
695 call this%ocd%ocd_ot(ipflg, kstp, endofperiod, this%iout, &
696 iprint_opt=idvprint, isav_opt=0)
697 end if
698
tdismodule::endofperiod
logical(lgp), pointer, public endofperiod
flag indicating end of stress period
Definition tdis.f90:30

◆ output_immobile_sorbate_concentration()

subroutine gwtistmodule::output_immobile_sorbate_concentration ( class(gwtisttype) this,
integer(i4b), intent(in) idvsave,
integer(i4b), intent(in) idvprint )
private
Parameters
thisBndType object
[in]idvsaveflag and unit number for dependent-variable output
[in]idvprintflag indicating if dependent-variable should be written to the model listing file

Definition at line 703 of file gwt-ist.f90.

704 ! modules
705 ! dummy
706 class(GwtIstType) :: this !< BndType object
707 integer(I4B), intent(in) :: idvsave !< flag and unit number for dependent-variable output
708 integer(I4B), intent(in) :: idvprint !< flag indicating if dependent-variable should be written to the model listing file
709 ! local
710 character(len=1) :: cdatafmp = ' ', editdesc = ' '
711 integer(I4B) :: ibinun
712 integer(I4B) :: iprint, nvaluesp, nwidthp
713 real(DP) :: dinact
714
715 ! Save cimsrb to a binary file. ibinun is a flag where 1 indicates that
716 ! cim should be written to a binary file if a binary file is open for it.
717 ! Set unit number for sorbate output
718 if (this%ioutsorbate /= 0) then
719 ibinun = 1
720 else
721 ibinun = 0
722 end if
723 if (idvsave == 0) ibinun = 0
724
725 ! save sorbate concentration array
726 if (ibinun /= 0) then
727 iprint = 0
728 dinact = dhnoflo
729 if (this%ioutsorbate /= 0) then
730 ibinun = this%ioutsorbate
731 call this%dis%record_array(this%cimsrb, this%iout, iprint, ibinun, &
732 ' SORBATE', cdatafmp, nvaluesp, &
733 nwidthp, editdesc, dinact)
734 end if
735 end if
736

◆ read_options()

subroutine gwtistmodule::read_options ( class(gwtisttype), intent(inout) this)
private

Read options for boundary packages.

Parameters
[in,out]thisGwtIstType object

Definition at line 1133 of file gwt-ist.f90.

1134 ! -- modules
1135 ! -- dummy
1136 class(GwtIstType), intent(inout) :: this !< GwtIstType object
1137
1138 ! -- source options
1139 call this%source_options()

◆ source_data()

subroutine gwtistmodule::source_data ( class(gwtisttype) this)
private

Method to source data for the package.

Definition at line 1146 of file gwt-ist.f90.

1147 ! -- modules
1148 use constantsmodule, only: linelength
1149 use simvariablesmodule, only: errmsg, warnmsg
1150 use simmodule, only: count_errors, store_error, store_warning, &
1151 store_error_filename
1152 use memorymanagermodule, only: get_isize, mem_reallocate
1153 use memorymanagerextmodule, only: mem_set_value
1154 use gwtistinputmodule, only: gwtistparamfoundtype
1155 ! -- dummy
1156 class(GwtIsttype) :: this
1157 ! -- locals
1158 !character(len=LINELENGTH) :: errmsg
1159 type(GwtIstParamFoundType) :: found
1160 integer(I4B) :: asize
1161 integer(I4B), dimension(:), pointer, contiguous :: map
1162 !
1163 ! -- set map to convert user input data into reduced data
1164 map => null()
1165 if (this%dis%nodes < this%dis%nodesuser) map => this%dis%nodeuser
1166 !
1167 ! -- reallocate
1168 if (this%isrb == 0) then
1169 call get_isize('BULK_DENSITY', this%input_mempath, asize)
1170 if (asize > 0) &
1171 call mem_reallocate(this%bulk_density, this%dis%nodes, &
1172 'BULK_DENSITY', this%memoryPath)
1173 call get_isize('DISTCOEF', this%input_mempath, asize)
1174 if (asize > 0) &
1175 call mem_reallocate(this%distcoef, this%dis%nodes, 'DISTCOEF', &
1176 this%memoryPath)
1177 end if
1178 if (this%idcy == 0) then
1179 call get_isize('DECAY', this%input_mempath, asize)
1180 if (asize > 0) &
1181 call mem_reallocate(this%decay, this%dis%nodes, 'DECAY', this%memoryPath)
1182 end if
1183 call get_isize('DECAY_SORBED', this%input_mempath, asize)
1184 if (asize > 0) then
1185 call mem_reallocate(this%decay_sorbed, this%dis%nodes, &
1186 'DECAY_SORBED', this%memoryPath)
1187 end if
1188 if (this%isrb < 2) then
1189 call get_isize('SP2', this%input_mempath, asize)
1190 if (asize > 0) &
1191 call mem_reallocate(this%sp2, this%dis%nodes, 'SP2', this%memoryPath)
1192 end if
1193 !
1194 ! -- update defaults with memory sourced values
1195 call mem_set_value(this%porosity, 'POROSITY', this%input_mempath, map, &
1196 found%porosity)
1197 call mem_set_value(this%volfrac, 'VOLFRAC', this%input_mempath, map, &
1198 found%volfrac)
1199 call mem_set_value(this%zetaim, 'ZETAIM', this%input_mempath, map, &
1200 found%zetaim)
1201 call mem_set_value(this%cim, 'CIM', this%input_mempath, map, &
1202 found%cim)
1203 call mem_set_value(this%decay, 'DECAY', this%input_mempath, map, &
1204 found%decay)
1205 call mem_set_value(this%decay_sorbed, 'DECAY_SORBED', this%input_mempath, &
1206 map, found%decay_sorbed)
1207 call mem_set_value(this%bulk_density, 'BULK_DENSITY', this%input_mempath, &
1208 map, found%bulk_density)
1209 call mem_set_value(this%distcoef, 'DISTCOEF', this%input_mempath, map, &
1210 found%distcoef)
1211 call mem_set_value(this%sp2, 'SP2', this%input_mempath, map, &
1212 found%sp2)
1213
1214 ! -- log options
1215 if (this%iout > 0) then
1216 call this%log_data(found)
1217 end if
1218
1219 ! -- Check for required sorption variables
1220 if (this%isrb > 0) then
1221 if (.not. found%bulk_density) then
1222 write (errmsg, '(a)') 'Sorption is active but BULK_DENSITY &
1223 &not specified. BULK_DENSITY must be specified in GRIDDATA block.'
1224 call store_error(errmsg)
1225 end if
1226 if (.not. found%distcoef) then
1227 write (errmsg, '(a)') 'Sorption is active but distribution &
1228 &coefficient not specified. DISTCOEF must be specified in &
1229 &GRIDDATA block.'
1230 call store_error(errmsg)
1231 end if
1232 if (this%isrb > 1) then
1233 if (.not. found%sp2) then
1234 write (errmsg, '(a)') 'Freundlich or langmuir sorption is active &
1235 &but SP2 not specified. SP2 must be specified in &
1236 &GRIDDATA block.'
1237 call store_error(errmsg)
1238 end if
1239 end if
1240 else
1241 if (found%bulk_density) then
1242 write (warnmsg, '(a)') 'Sorption is not active but &
1243 &BULK_DENSITY was specified. BULK_DENSITY will have no affect on &
1244 &simulation results.'
1245 call store_warning(warnmsg)
1246 write (this%iout, '(1x,a)') 'WARNING. '//warnmsg
1247 end if
1248 if (found%distcoef) then
1249 write (warnmsg, '(a)') 'Sorption is not active but &
1250 &distribution coefficient was specified. DISTCOEF will have &
1251 &no affect on simulation results.'
1252 call store_warning(warnmsg)
1253 write (this%iout, '(1x,a)') 'WARNING. '//warnmsg
1254 end if
1255 if (found%sp2) then
1256 write (warnmsg, '(a)') 'Sorption is not active but &
1257 &SP2 was specified. SP2 will have &
1258 &no affect on simulation results.'
1259 call store_warning(warnmsg)
1260 write (this%iout, '(1x,a)') 'WARNING. '//warnmsg
1261 end if
1262 end if
1263
1264 ! -- Check for required decay/production rate coefficients
1265 if (this%idcy > 0) then
1266 if (.not. found%decay) then
1267 write (errmsg, '(a)') 'First or zero order decay is &
1268 &active but the first rate coefficient is not specified. DECAY &
1269 &must be specified in GRIDDATA block.'
1270 call store_error(errmsg)
1271 end if
1272 if (.not. found%decay_sorbed) then
1273 !
1274 ! -- If DECAY_SORBED not specified and sorption is active, then
1275 ! terminate with an error
1276 write (errmsg, '(a)') 'DECAY_SORBED not provided in GRIDDATA &
1277 &block but decay and sorption are active. Specify DECAY_SORBED &
1278 &in GRIDDATA block.'
1279 call store_error(errmsg)
1280 end if
1281 else
1282 if (found%decay) then
1283 write (warnmsg, '(a)') 'First- or zero-order decay &
1284 &is not active but decay was specified. DECAY will &
1285 &have no affect on simulation results.'
1286 call store_warning(warnmsg)
1287 write (this%iout, '(1x,a)') 'WARNING. '//warnmsg
1288 end if
1289 if (found%decay_sorbed) then
1290 write (warnmsg, '(a)') 'First- or zero-order decay &
1291 &is not active but DECAY_SORBED was specified. &
1292 &DECAY_SORBED will have no affect on simulation results.'
1293 call store_warning(warnmsg)
1294 write (this%iout, '(1x,a)') 'WARNING. '//warnmsg
1295 end if
1296 end if
1297
1298 ! -- Check for required dual domain arrays or warn if they are specified
1299 ! but won't be used.
1300 if (.not. found%cim) then
1301 write (this%iout, '(1x,a)') 'Warning. Dual domain is active but &
1302 &initial immobile domain concentration was not specified. &
1303 &Setting CIM to zero.'
1304 end if
1305 if (.not. found%zetaim) then
1306 write (errmsg, '(a)') 'Dual domain is active but dual &
1307 &domain mass transfer rate (ZETAIM) was not specified. ZETAIM &
1308 &must be specified in GRIDDATA block.'
1309 call store_error(errmsg)
1310 end if
1311 if (.not. found%porosity) then
1312 write (errmsg, '(a)') 'Dual domain is active but &
1313 &immobile domain POROSITY was not specified. POROSITY &
1314 &must be specified in GRIDDATA block.'
1315 call store_error(errmsg)
1316 end if
1317 if (.not. found%volfrac) then
1318 write (errmsg, '(a)') 'Dual domain is active but &
1319 &immobile domain VOLFRAC was not specified. VOLFRAC &
1320 &must be specified in GRIDDATA block. This is a new &
1321 &requirement for MODFLOW versions later than version &
1322 &6.4.1.'
1323 call store_error(errmsg)
1324 end if
1325
1326 ! -- terminate if errors
1327 if (count_errors() > 0) then
1328 call store_error_filename(this%input_fname)
1329 end if
memorymanagerextmodule::mem_set_value
Definition MemoryManagerExt.f90:17
memorymanagermodule::mem_reallocate
Definition MemoryManager.f90:82
constantsmodule::linelength
integer(i4b), parameter linelength
maximum length of a standard line
Definition Constants.f90:45
memorymanagermodule::get_isize
subroutine, public get_isize(name, mem_path, isize)
@ brief Get the number of elements for this variable
Definition MemoryManager.f90:276
simmodule::store_warning
subroutine, public store_warning(msg, substring)
Store warning message.
Definition Sim.f90:237
simvariablesmodule::errmsg
character(len=maxcharlen) errmsg
error message string
Definition SimVariables.f90:28
simvariablesmodule::warnmsg
character(len=maxcharlen) warnmsg
warning message string
Definition SimVariables.f90:29
Here is the call graph for this function:

◆ source_options()

subroutine gwtistmodule::source_options ( class(gwtisttype), intent(inout) this)
private

Method to source options for the package.

Definition at line 960 of file gwt-ist.f90.

961 ! -- modules
962 use constantsmodule, only: lenvarname, linelength, mnormal, lenbigline
963 use simmodule, only: store_error, store_error_filename
964 use openspecmodule, only: access, form
965 use inputoutputmodule, only: getunit, assign_iounit, openfile
966 use memorymanagerextmodule, only: mem_set_value
967 use gwtistinputmodule, only: gwtistparamfoundtype
968 ! -- dummy
969 class(GwtIstType), intent(inout) :: this
970 ! -- locals
971 character(len=LINELENGTH) :: prnfmt
972 integer(I4B), pointer :: columns, width, digits
973 type(GwtIstParamFoundType) :: found
974 character(len=LENVARNAME), dimension(3) :: sorption_method = &
975 &[character(len=LENVARNAME) :: 'LINEAR', 'FREUNDLICH', 'LANGMUIR']
976 character(len=linelength) :: sorbate_fname, cim6_fname, budget_fname, &
977 budgetcsv_fname
978 allocate (columns)
979 allocate (width)
980 allocate (digits)
981 !
982 ! -- update defaults with memory sourced values
983 call mem_set_value(this%ipakcb, 'SAVE_FLOWS', this%input_mempath, &
984 found%save_flows)
985 call mem_set_value(budget_fname, 'BUDGETFILE', this%input_mempath, &
986 found%budgetfile)
987 call mem_set_value(budgetcsv_fname, 'BUDGETCSVFILE', this%input_mempath, &
988 found%budgetcsvfile)
989 call mem_set_value(this%isrb, 'SORPTION', this%input_mempath, &
990 sorption_method, found%sorption)
991 call mem_set_value(this%idcy, 'ORDER1_DECAY', this%input_mempath, &
992 found%order1_decay)
993 call mem_set_value(this%idcy, 'ORDER0_DECAY', this%input_mempath, &
994 found%order0_decay)
995 call mem_set_value(cim6_fname, 'CIMFILE', this%input_mempath, &
996 found%cimfile)
997 call mem_set_value(sorbate_fname, 'SORBATEFILE', this%input_mempath, &
998 found%sorbatefile)
999 call mem_set_value(columns, 'COLUMNS', this%input_mempath, &
1000 found%columns)
1001 call mem_set_value(width, 'WIDTH', this%input_mempath, &
1002 found%width)
1003 call mem_set_value(digits, 'DIGITS', this%input_mempath, &
1004 found%digits)
1005 call mem_set_value(prnfmt, 'FORMAT', this%input_mempath, &
1006 found%format)
1007
1008 ! -- found side effects
1009 if (found%save_flows) this%ipakcb = -1
1010 if (found%budgetfile) then
1011 call assign_iounit(this%ibudgetout, this%inunit, "BUDGET fileout")
1012 call openfile(this%ibudgetout, this%iout, budget_fname, 'DATA(BINARY)', &
1013 form, access, 'REPLACE', mode_opt=mnormal)
1014 end if
1015 if (found%budgetcsvfile) then
1016 call assign_iounit(this%ibudcsv, this%inunit, "BUDGETCSV fileout")
1017 call openfile(this%ibudcsv, this%iout, budgetcsv_fname, 'CSV', &
1018 filstat_opt='REPLACE')
1019 end if
1020 if (found%sorption) then
1021 if (this%isrb == 0) then
1022 call store_error('Unknown sorption type was specified. &
1023 &Sorption must be specified as LINEAR, &
1024 &FREUNDLICH, or LANGMUIR.')
1025 call store_error_filename(this%input_fname)
1026 end if
1027 end if
1028 if (found%order1_decay) this%idcy = 1
1029 if (found%order0_decay) this%idcy = 2
1030 if (found%cimfile) then
1031 call this%ocd%set_ocfile(cim6_fname, this%iout)
1032 end if
1033 if ((found%columns .or. found%width .or. found%digits) .and. &
1034 .not. (found%columns .and. found%width .and. found%digits)) then
1035 call store_error('PRINT_FORMAT must include COLUMNS, WIDTH and DIGITS. &
1036 &Syntax is: COLUMNS <columns> WIDTH <width> &
1037 &DIGITS <digits> [<format>].')
1038 call store_error_filename(this%input_fname)
1039 end if
1040 if (found%columns .and. found%width .and. &
1041 found%digits .and. found%format) then
1042 write (this%lstfmt, '(a,i0,a,i0,a,i0,a)') 'COLUMNS ', columns, &
1043 ' WIDTH ', width, ' DIGITS ', digits, ' '//trim(prnfmt)
1044 end if
1045 if (found%sorbatefile) then
1046 this%ioutsorbate = getunit()
1047 call openfile(this%ioutsorbate, this%iout, sorbate_fname, &
1048 'DATA(BINARY)', form, access, 'REPLACE', mode_opt=mnormal)
1049 end if
1050 !
1051 ! -- log options
1052 if (this%iout > 0) then
1053 call this%log_options(found, cim6_fname, budget_fname, &
1054 budgetcsv_fname, sorbate_fname)
1055 end if
1056
1057 deallocate (columns)
1058 deallocate (width)
1059 deallocate (digits)
constantsmodule::mnormal
@ mnormal
normal output mode
Definition Constants.f90:207
constantsmodule::lenbigline
integer(i4b), parameter lenbigline
maximum length of a big line
Definition Constants.f90:15
constantsmodule::lenvarname
integer(i4b), parameter lenvarname
maximum length of a variable name
Definition Constants.f90:17
inputoutputmodule::assign_iounit
subroutine, public assign_iounit(iounit, errunit, description)
@ brief assign io unit number
Definition InputOutput.f90:185
inputoutputmodule::getunit
integer(i4b) function, public getunit()
Get a free unit number.
Definition InputOutput.f90:166
inputoutputmodule::openfile
subroutine, public openfile(iu, iout, fname, ftype, fmtarg_opt, accarg_opt, filstat_opt, mode_opt)
Open a file.
Definition InputOutput.f90:30
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:

Variable Documentation

◆ budtxt

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

Definition at line 38 of file gwt-ist.f90.

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

◆ ftype

character(len=lenftype) gwtistmodule::ftype = 'IST'
private

Definition at line 35 of file gwt-ist.f90.

◆ nbditems

integer(i4b), parameter gwtistmodule::nbditems = 5
private

Definition at line 37 of file gwt-ist.f90.

◆ text

character(len=lenpackagename) gwtistmodule::text = ' IMMOBILE DOMAIN'
private

Definition at line 36 of file gwt-ist.f90.