R/wtConf.R

In KWB-R/kwb.wtaq: Interface to WTAQ Drawdown Model (http://water.usgs.gov/ogw/wtaq/)

# wtConfigureGeneral -----------------------------------------------------------
wtConfigureGeneral <- function # General parameter values
### Define general WTAQ parameters
(
  title = paste("Input file generated by kwb.wtaq-functions on", date()),
  ### Title of simulation; up to 70 characters in length. Leave this line blank 
  ### if no title is specified.
  format = "DIMENSIONAL"
  ### Analysis format. Enter DIMENSIONAL.
)
{
  ##seealso<< \code{\link{wtConfigure}}
  x <- list(title = title, format = format)
  x
  
  ### list with elements \emph{title} and \emph{format}, representing general 
  ### WTAQ parameters. See descriptions in Arguments section.
}

# wtConfigureAquifer -----------------------------------------------------------
wtConfigureAquifer <- function # Aquifer parameter values
### Define aquifer-related WTAQ parameters
(
  aqtype = "WATER TABLE",
  ### Type of aquifer being simulated. Two options are provided:
  ### AQTYPE = CONFINED or
  ### AQTYPE = WATER TABLE
  bb = 20,  
  ### Thickness or saturated thickness of aquifer at beginning of simulation, 
  ### in units of length.
  hkr = 4.10E-03, 
  ### Horizontal hydraulic conductivity of aquifer, in units of length per time.
  hkz = 1.74E-03, 
  ### Vertical hydraulic conductivity of aquifer, in units of length per time.
  ss = 3.76E-05, 
  ### Specific storage of aquifer, in units of inverse length.
  sy = 0.25
  ### Specific yield of aquifer, dimensionless. Enter 0 if 
  ### AQTYPE = CONFINED.
)
{
  ##seealso<< \code{\link{wtConfigure}}    
  x <- list(aqtype = aqtype, bb = bb, hkr = hkr, hkz = hkz, ss = ss, sy = sy)
  x
  
  ### list with elements \emph{aqtype}, \emph{bb}, \emph{hkr}, \emph{hkz},
  ### \emph{ss}, \emph{sy}, representing WTAQ parameters related to the aquifer.
  ### See descriptions in Arguments section.
}

# wtConfigureDrainage ----------------------------------------------------------
wtConfigureDrainage <- function # Drainage parameter values
### Define drainage-related WTAQ parameters
(
  idra = 2, 
  ### Type of drainage at water table. Enter 0 if AQTYPE = CONFINED. 
  ### Three options are provided:
  ### IDRA = 0: Instantaneous drainage.
  ### IDRA = 1: Gradual drainage.
  ### IDRA = 2: Drainage with unsaturated-zone characterization.
  alpha = c(), 
  ### Only used if idra = 1. Drainage constants, in units of inverse time. 
  ### Maximum of 5 values is allowed.
  acc = 5.0, 
  ### Only used if idra = 2. Soil-moisture retention exponent, in units of 
  ### inverse length.
  akk = 31.7, 
  ### Only used if idra = 2. Relative hydraulic-conductivity exponent, in units 
  ### of inverse length. The value specified must be greater than or equal to 
  ### that specified for ACC.
  amm = 100, 
  ### Only used if idra = 2. Initial unsaturated-zone thickness above the 
  ### capillary fringe, in units of length.
  axmm = 10
  ### Only used if idra = 2. The unsaturated-zone thickness above the capillary 
  ### fringe above which an assumption of an infinitely thick unsaturated-zone 
  ### thickness is assumed, in units of length.
)
{
  ##seealso<< \code{\link{wtConfigure}}
  x <- list(idra = idra, alpha = alpha, acc = acc, akk = akk, amm = amm, axmm = axmm)
  x  

  ### list with elements \emph{idra}, \emph{alpha}, \emph{acc}, \emph{akk},
  ### \emph{amm}, \emph{axmm}, representing WTAQ parameters related to drainage.
  ### See descriptions in Arguments section.
}

# wtConfigureTimes -------------------------------------------------------------
wtConfigureTimes <- function # Simulation time parameter values
### Define WTAQ parameters related to the times to be simulated
(
  its = 1, 
  ### Time specification:
  ### ITS = 0: Log-cycle time steps (use to generate theoretical curves).
  ### ITS = 1: User-specified times.
  tlast = 0.0, 
  ### Largest value of time. Enter 0 if ITS = 1.
  nlc = 0, 
  ### Number of logarithmic cycles on the time scale for which drawdown will 
  ### be calculated. Enter 0 if ITS = 1.
  nox = 0 
  ### Number of equally spaced times per logarithmic cycle for which drawdown 
  ### will be calculated. Enter 0 if ITS = 1.
)
{
  ##seealso<< \code{\link{wtConfigure}}   
  x <- list(its = its, tlast = tlast, nlc = nlc, nox = nox)
  x  

  ### list with elements \emph{its}, \emph{tlast}, \emph{nlc}, \emph{nox},
  ### representing WTAQ parameters defining the times to be simulated.
  ### See descriptions in Arguments section.
}

# wtConfigureSolution ----------------------------------------------------------
wtConfigureSolution <- function # Solver-related parameter values
### Define WTAQ parameters related to the solver algorithms
(
  isoln = 2, 
  ### Numerical-inversion solution type:
  ### ISOLN = 1: Solution by the Stehfest algorithm (must use this option for 
  ### confined aquifers).
  ### ISOLN = 2: Solution by the de Hoog algorithm (must use this option for 
  ### drainage with unsaturated-zone characterization [wtConfigureDrainage()$idra == 2]).
  rerrnr = 1.0E-10, 
  ### Relative error for Newton-Raphson iteration and finite summations of 
  ### drawdown for water-table aquifers. A value of 1.0E-10 is suggested. 
  ### Enter 0.0D0 for AQTYPE = CONFINED.
  rerrsum = 0.0, 
  ### Only used if isoln = 1. Relative error for finite summations of drawdown 
  ### for confined aquifers. Suggested value is 1.0E-07 to 1.0E-08. Enter 0
  ### if AQTYPE = WATER TABLE.
  nmax = 0, 
  ### Only used if isoln = 1. Maximum number of terms permitted in the finite 
  ### summations of drawdown for confined aquifers. Suggested value is 200. 
  ### Enter 0 if AQTYPE = WATER TABLE.
  ntms = 30, 
  ### Factor used to determine number of terms in the finite summations for 
  ### drawdown for water-table aquifers. Suggested values are 20 or 30. 
  ### Enter 0 if AQTYPE = CONFINED.
  ns = 8, 
  ### Only used if isoln = 1. Number of terms used in the Stehfest algorithm. 
  ### This must be an even integer, the value of which depends upon computer 
  ### precision. If the computer holds 16 significant figures in double 
  ### precision, let NS = 6 to 12. A value of 8 is recommended.
  error = 1.0E-04, 
  ### Only used if isoln = 2. Relative error sought for the accuracy of the 
  ### numerical inversion. A value of 1.0E-04 is suggested.
  nnn = 6, 
  ### Only used if isoln = 2. Number of terms used in the summation of the 
  ### Fourier series of the approximation to the inverse Laplace transform. 
  ### A value of 6 is suggested.
  method = 3
  ### Only used if isoln = 2. Indicates which method will be used to accelerate 
  ### convergence of the Fourier series. Options are 1, 2, or 3. Only METHOD = 3
  ### has been tested and was found to be satisfactory. Users can consult de 
  ### Hoog and others (1982) and John Knight's subroutine LAPADC for additional 
  ### details if needed.
) 
{
  ##seealso<< \code{\link{wtConfigure}}   
  x <- list(isoln = isoln, rerrnr = rerrnr, rerrsum = rerrsum, nmax = nmax, 
            ntms = ntms, ns = ns, error = error, nnn = nnn, method = method)
  x

  ### list with elements \emph{isoln}, \emph{rerrnr}, \emph{rerrsum}, \emph{nmax},
  ### \emph{ntms}, \emph{ns}, \emph{error}, \emph{nnn}, \emph{method},
  ### representing WTAQ parameters related to the solving algorithms used in WTAQ.
  ### See descriptions in Arguments section.
}

# wtConfigurePumpwell ----------------------------------------------------------
wtConfigurePumpwell <- function # Pumped well parameter values
### Define WTAQ parameters related to the pumped well
(
  ipws = 0, 
  ### Type of pumped well:
  ### IPWS = 0: Partially penetrating pumped well.
  ### IPWS = 1: Fully penetrating pumped well. Default: 0
  ipwd = 1, 
  ### Type of diameter of pumped well:
  ### IPWD = 0: Infinitesimal diameter (line-source theory).
  ### IPWD = 1: Finite diameter. Default: 1
  rw = 1.0, 
  ### Radius of pumped well screen, in units of length. Default: 1.0
  rc = ifelse(ipwd == 0, 0, rw), 
  ### Inside radius of pumped well in the interval where water levels are 
  ### changing during pumping, in units of length. Must be 0 if IPWD = 0.
  ### Default: value of rw, but 0 if ipwd = 0.
  zpd = 5.0, 
  ### Depth below top of aquifer or initial water table to the top of the 
  ### screened interval of the pumped well, in units of length. Default: 5.0
  zpl = 10.0,
  ### Depth below top of aquifer or initial water table to the bottom of the 
  ### screened interval of the pumped well, in units of length. Default: 10.0
  sw = 0.0, 
  ### Well-bore skin parameter, dimensionless. Default: 0.0
  qq = 0.05, 
  ### Pumping rate of well, in units of cubic length per time. Default: 0.05
  tspw = data.frame(t = 60 * 1:3, dd = 2 * c(0.4, 0.6, 0.7)),
  ### data.frame with column \emph{t} containing the user-specified times for 
  ### which drawdown at the pumped well will be calculated. If the data.frame
  ### has no rows, no drawdowns are calculated for the pumped well. The 
  ### data.frame can (optionally) contain the drawdowns that have been measured
  ### for the pumped well at the corresponding times in an additional column 
  ### \emph{dd}. 
  pwname = "PW",
  ### Name of pumping well. This parameter is not evaluated by WTAQ but will be 
  ### used in the result data frame returned by
  ### \code{\link{wtRunConfiguration}}. We introduce it in order to be able to
  ### identify the well within a wellfield. Default: "PW"
  ipump = 1,
  ### Option to suppress calculations of drawdown at pumped well:
  ### IPUMP = 0: Drawdown is not calculated at pumped well.
  ### IPUMP = 1: Drawdown is calculated at pumped well.  
  irun = 1
  ### Option to suppress drawdown calculations for the pumped well. Allows user 
  ### to specify time-drawdown data, but those data are ignored during the 
  ### simulation. Options are:
  ### IRUN = 0: Drawdowns not calculated.
  ### IRUN = 1: Drawdowns calculated.  
)
{
  ##seealso<< \code{\link{wtConfigure}}
  
  if (ipwd == 0 && rc != 0) {
    rc <- 0
    warning("rc was set to 0 since ipwd is 0 (Infinitesimal diameter)")
  }
  
  x <- list(
    irun = irun, ipws = ipws, ipwd = ipwd, ipump = ipump, qq = qq, rw = rw,
    rc = rc, zpd = zpd, zpl = zpl, sw = sw, tspw = tspw, pwname = pwname)
  x  

  ### list with elements \emph{irun}, \emph{ipws}, \emph{ipwd}, \emph{ipump}, 
  ### \emph{qq}, \emph{rw}, \emph{rc}, \emph{zpd}, \emph{zpl}, \emph{sw},
  ### \emph{tspw}, representing WTAQ parameters related to the pumped well. See
  ### descriptions in Arguments section.
}

# wtConfigureObservationWell ---------------------------------------------------
wtConfigureObservationWell <- function # Observation well parameter values
### Define WTAQ parameters related to an obervation well
(
  obname = "obs1", 
  ### Name of observation well or piezometer; up to 10 characters in length.
  r = 30.0, 
  ### Radial distance from axis of pumped well to observation well or 
  ### piezometer, in units of length.
  iows = 0, 
  ### Type of observation well or piezometer:
  ### IOWS = 0: Partially penetrating observation well.
  ### IOWS = 1: Fully penetrating observation well.
  ### IOWS = 2: Observation piezometer. Default: 0
  idpr = 0, 
  ### Options for delayed response of observation well.
  ### IDPR = 0: No delayed response. 
  ### IDPR = 1: Delayed response. Default: 0
  rp = ifelse(idpr == 0, 0, 0.1), 
  ### Inside radius of the observation well (or piezometer) standpipe in the 
  ### interval over which water levels are changing during pumping, in units of 
  ### length. Enter 0 if IDPR = 0. Default: 0.1, but 0 if idpr = 0.
  z1 = ifelse(iows == 2, 0, 5.0), 
  ### Depth below top of aquifer or initial water table to the top of screened 
  ### interval of observation well, in units of length. Use for IOWS = 0 or 1. 
  ### Enter 0 if IOWS = 2. Default: 5.0, but 0 if iows = 2.
  z2 = ifelse(iows == 2, 0, 10.0), 
  ### Depth below top of aquifer or initial water table to the bottom of 
  ### screened interval of observation well, in units of length. Use for 
  ### IOWS = 0 or 1. Enter 0 if IOWS = 2. Default: 10.0, but 0 if iows = 2.
  xll = ifelse(idpr == 0, 0, z2 - z1),
  ### Length of screened interval of observation well or piezometer, in units 
  ### of length. Must be 0 if IDPR = 0. Default: z2 - z1, but 0 if idpr = 0. 
  zp = ifelse(iows == 2, 3.0, 0.0),
  ### Depth below top of aquifer or initial water table to center of piezometer, 
  ### in units of length. Use for IOWS = 2. Enter 0.0 if IOWS = 0 or 1.
  ### Default: 3.0, but 0.0 if iows not equal to 2.
  tsobs = data.frame(t = 60 * 1:3, dd = c(0.4, 0.6, 0.7)),
  ### data.frame with column \emph{t} containing the user-specified times for 
  ### which drawdown at the observation well or piezometer will be calculated.
  ### If the data.frame has no rows, no drawdowns are calculated for the 
  ### observation well or piezometer. The data.frame can (optionally) contain 
  ### the measured drawdown for the observation well or piezometer at the 
  ### corresponding times in an additional column \emph{dd}.   
  irun = 1
  ### Option to suppress drawdown calculations for the observation well or 
  ### piezometer. Allows user to specify time-drawdown data, but those data are 
  ### ignored during the simulation. Options are:
  ### IRUN = 0: Drawdowns not calculated.
  ### IRUN = 1: Drawdowns calculated.  
)
{
  ##seealso<< \code{\link{wtConfigure}}    
  
  if (idpr == 0 && xll != 0) {
    xll <- 0
    warning("xll was set to 0 since idpr is 0 (No delayed response)")
  }
  
  if (idpr == 0 && rp != 0) {
    rp <- 0
    warning("rp was set to 0 since idpr is 0 (No delayed response)")
  }
  
  x <- list(irun = irun, obname = obname, iows = iows, idpr = idpr, r = r, z1 = z1, 
            z2 = z2, zp = zp, rp = rp, xll = xll, tsobs = tsobs)
  x  
  
  ### list with elements \emph{irun}, \emph{obname}, \emph{iows}, \emph{idpr},
  ### \emph{r}, \emph{z1}, \emph{z2}, \emph{zp}, \emph{rp}, \emph{xll}, \emph{tsobs},
  ### representing WTAQ parameters related to an observation well. 
  ### See description in Arguments section.
}

# wtConfiguredWellnames --------------------------------------------------------
wtConfiguredWellnames <- function
### Read all well names (pumping well \dQuote{PW} plus observation wells) from
### WTAQ configuration \emph{configuration}
(
  configuration
  ### WTAQ configuration as generated by \code{\link{wtConfigure}}
) 
{
  pumpingWellName <- configuration$pumpwell$pwname
  c(pumpingWellName, sapply(configuration$obswells, FUN = "[[", "obname"))  

  ### vector of character containing the name of the pumping well (as the first
  ### element) and the names of the observation wells as following elements
  ### (according to the order of their definition in \emph{configuration}
}

# wtConfiguredDistances --------------------------------------------------------
wtConfiguredDistances <- function # distances of observation wells
### Vector of distances between pumping well and observation wells (including
### 0 as first distance!)
(
  configuration
  ### WTAQ configuration as retrieved by \code{\link{wtConfigure}}.
) 
{
  distances <- c(0, sapply(configuration$obswells, FUN = "[[", "r"))
  names(distances) <- wtConfiguredWellnames(configuration)
  
  distances
  ### named vector of numeric containing zero as first element and the distances
  ### between pumping well and observation wells as following elements, in the
  ### order or their definition in the WTAQ configuration. The names of the
  ### elements correspond to the well names as retrieved also with 
  ### \code{\link{wtConfiguredWellnames}}
}

# wtCheckConfiguration ---------------------------------------------------------
wtCheckConfiguration <- structure(
  function # check WTAQ configuration
  ### Checks a WTAQ configuration for its structure
  (
    configuration, 
    ### WTAQ configuration (complete) or sub-configuration (general, aquifer, 
    ### rainage, times, solution, pumpwell, obswell)
    part = "complete",
    ### one of c("complete", "general", "aquifer", "drainage", "times",
    ### "solution", "pumpwell", "obswell")
    dbg = TRUE
    ### if TRUE, debug messages are shown, else not.
  ) 
  {
    msg <- ""
    
    # list of functions that can be used to obain the required list fields
    cfun <- list(complete = wtConfigure,
                 general  = wtConfigureGeneral,
                 aquifer  = wtConfigureAquifer,
                 drainage = wtConfigureDrainage,
                 times    = wtConfigureTimes,
                 solution = wtConfigureSolution,
                 pumpwell = wtConfigurePumpwell,
                 obswell  = wtConfigureObservationWell)
    
    parts <- names(cfun)
    if(!(part %in% parts)) {
      stop("part must be one of ", 
           paste('"', parts, '"', sep = "", collapse = ", "))
    }
    
    # required list elements
    req <- names(cfun[[part]]())
    
    if (length(intersect(names(configuration), req)) != length(req)) {
      msg <- paste(sprintf(
        paste("\n*** Configuration \"%s\" does not contain all of these",
              "required elements: %s.\n"), 
        part,
        paste('"', req, '"', sep = "", collapse = ", ")))
      cat(msg)
      return(msg)
    }
    
    # if configuration is a complete configuration (part == "complete") check all the
    # sub-configurations
    if (part == "complete") {
      
      # loop through part names (excluding "complete")
      for (p in parts[-1]) {
        
        kwb.utils::catIf(dbg, sprintf('Checking "%s" configuration... ', p))
        
        if (p == "obswell") {
          
          # loop through observation wells
          for (i in seq(1, by = 1, length.out=length(configuration$obswells))) {
            
            kwb.utils::catIf(dbg, sprintf('\n- Checking config of observation well %i...', i))
            
            msg <- wtCheckConfiguration(configuration$obswells[[i]], p)
            if (msg != "") {
              return(msg)
            } 
            
            kwb.utils::catIf(dbg, "ok.")
          }
        }
        else {
          msg <- wtCheckConfiguration(configuration[[p]], p)        
        }
        if (msg != "") {
          return(msg)
        }
        
        kwb.utils::catIf(dbg, "ok.\n")              
      }
    }
    
    msg
    
    ### Error message or "" if no error occurred
  }, ex = function() {
    
    # Generate default configuration
    cconf <- wtConfigure()  
    
    # Check (complete) configuration. No error -> ok
    wtCheckConfiguration(cconf) 
    
    # Generate default observation well configuration
    oconf <- wtConfigureObservationWell()
    
    # Check "obswell" configuration. No error -> ok
    wtCheckConfiguration(oconf, "obswell") 
    
    # Set obswells in complete configuration (forgetting that obswells must be
    # a list!):
    cconf$obswells <- oconf 
    
    # Check (complete) configuration again 
    # -> Error when "Checking config of observation well 1".
    wtCheckConfiguration(cconf) 
  }
)

# wtConfigure ------------------------------------------------------------------
wtConfigure <- function # Define parameter values for WTAQ model run
### Define parameter values for WTAQ model run
(
  general = wtConfigureGeneral(), 
  ### List of general parameters as retrieved by \code{\link{wtConfigureGeneral}}.
  aquifer = wtConfigureAquifer(), 
  ### List of aquifer parameters as retrieved by \code{\link{wtConfigureAquifer}}.  
  drainage = wtConfigureDrainage(),
  ### List of drainage parameters as retrieved by \code{\link{wtConfigureDrainage}}.  
  times = wtConfigureTimes(),
  ### List of simulation time parameters as retrieved by \code{\link{wtConfigureTimes}}.  
  solution = wtConfigureSolution(),
  ### List of solution parameters as retrieved by \code{\link{wtConfigureSolution}}.    
  pumpwell = wtConfigurePumpwell(),
  ### List of pumped well parameters as retrieved by \code{\link{wtConfigurePumpwell}}.    
  obswells = list(wtConfigureObservationWell(obname = "obs1"), wtConfigureObservationWell(obname = "obs2", r = 50))
  ### List of lists of observation well parameters each of which can be 
  ### retrieved by \code{\link{wtConfigureObservationWell}}.    
)
{
  ##seealso<< \code{\link{wtConfigurationExample2}, \link{wtConfigurationExample3}, \link{wtRunConfiguration}}  
  
  # check if obswells is a list
  if (!is.null(obswells) && mode(obswells) != "list") {
    msg <- paste("obswells must be a list but it is of mode: ", mode(obswells))
  }
  
  req <- names(wtConfigureObservationWell())
  for (i in 1:length(obswells)) {
    if(!is.null(obswells[[i]]) && mode(obswells[[i]]) != "list") {
      stop("Observation well ", i, " is not of mode list but of mode ",
           mode(obswells[[i]]))
    }
  }
  x <- list(general = general, aquifer = aquifer, drainage = drainage, 
            times = times, solution = solution, pumpwell = pumpwell, 
            obswells = obswells)

  class(x) <- "wtaqConfiguration"
  x
  
  ### list with elements \emph{general}, \emph{aquifer}, \emph{drainage},
  ### \emph{times}, \emph{solution}, \emph{pumpwell}, \emph{obswells},
  ### representing a full WTAQ configuration. For the different elements see the
  ### descriptions in the Arguments section.
}

# wtDefaultConfiguration -------------------------------------------------------
wtDefaultConfiguration <- function # Default WTAQ configuration
### Default WTAQ configuration. 
(
) 
{
  ##seealso<< \code{\link{wtConfigurationExample2}, \link{wtConfigurationExample3}}
  
  wtConfigurationExample2()
  ### Currently: configuration for WTAQ sample problem 2, as also returned by
  ### \code{\link{wtConfigurationExample2}}
}

# wtDefaultConfigurationSolution -----------------------------------------------
wtDefaultConfigurationSolution <- function # default WTAQ solver configuration
### default configuration for the WTAQ solver
(
  aqtype,
  ### Type of aquifer being simulated. Two options are provided:
  ### AQTYPE = CONFINED or
  ### AQTYPE = WATER TABLE
  idra = ifelse(aqtype == "CONFINED", 0, 1)
  ### Type of drainage at water table. Enter 0 if AQTYPE = CONFINED. 
  ### Three options are provided:
  ### IDRA = 0: Instantaneous drainage.
  ### IDRA = 1: Gradual drainage.
  ### IDRA = 2: Drainage with unsaturated-zone characterization.
)
{
  if (!aqtype  %in% c("CONFINED", "WATER TABLE")) {
    stop("aqtype must be one either 'CONFINED' or 'WATER TABLE'")
  }
  
  confined <- (aqtype == "CONFINED")
  waterTable <- (aqtype == "WATER TABLE")
  
  # Set ISOLN
  if (confined) {
    isoln <- 1
  }
  else if (idra == 2) {
    isoln <- 2
  }
  else {
    isoln <- 1
  }
  
  # start with default values
  rerrsum <- 0 
  nmax <- 0
  ntms <- 20
  ns <- 8
  error <- 1.0e-4
  nnn <- 6
  method <- 3

  if (isoln == 1) {
    rerrnr <- ifelse(confined, 0.0, 1.0e-10)
    rerrsum <- ifelse(waterTable, 0.0, 5.5e-8)
    nmax <- ifelse(waterTable, 0, 200)
    ntms <- ifelse(confined, 0, 20)
  }
  else {
    rerrnr <- 1.0e-10
  }
  
  wtConfigureSolution(
    isoln = isoln, 
    rerrnr = rerrnr, 
    rerrsum = rerrsum, 
    nmax = nmax, 
    ntms = ntms, 
    ns = ns, 
    error = error, 
    nnn = nnn, 
    method = method)
}

# wtRepairConfiguration --------------------------------------------------------
wtRepairConfiguration <- function # repair configuration
### repair configuration by resetting parameter values that depend on others
(
  configuration,
  ### WTAQ configuration, as retrieved by \code{\link{wtConfigure}}.  
  dbg = FALSE
  ### if TRUE, debug messages are shown, else not
)
{
  setUnlessEquals <- function(x, x.name, mustBeValue, reason) {
    if (x != mustBeValue) {
      kwb.utils::catIf(dbg, sprintf("%s was set to %s since %s.\n", x.name, mustBeValue, reason))
    }
    mustBeValue
  }

  dimensional <- "DIMENSIONAL"
  confined <- "CONFINED"
  waterTable <- "WATER TABLE"
  
  if (configuration$general$format != dimensional) {
    configuration$general$format <- dimensional
    kwb.utils::catIf("configuration$general$format was set to \"DIMENSIONAL\".\n")
  }
  
  if (configuration$aquifer$aqtype == confined) {
    reason <- paste("configuration$aquifer$aqtype ==", hsQuoteChr(confined))
    
    configuration$aquifer$sy <- setUnlessEquals(
      configuration$aquifer$sy, "configuration$aquifer$sy", 0.0, reason)
    
    configuration$drainage$idra <- setUnlessEquals(
      configuration$drainage$idra, "configuration$drainage$idra", 0, reason)
    
    configuration$solution$isoln <- setUnlessEquals(
      configuration$solution$isoln, "configuration$solution$isoln", 1, reason)
    
    configuration$solution$rerrnr <- setUnlessEquals(
      configuration$solution$rerrnr, "configuration$solution$rerrnr", 0.0, reason)

    configuration$solution$ntms <- setUnlessEquals(
      configuration$solution$ntms, "configuration$solution$ntms", 0, reason)
  }
  else if (configuration$aquifer$aqtype == waterTable) {
    reason <- paste("configuration$aquifer$aqtype ==", hsQuoteChr(waterTable))
    
    configuration$solution$rerrsum <- setUnlessEquals(
      configuration$solution$rerrsum, "configuration$solution$rerrsum", 0.0, reason)
    
    configuration$solution$nmax <- setUnlessEquals(
      configuration$solution$nmax, "configuration$solution$nmax", 0, reason)
  }  
  
  if (configuration$drainage$idra %in% c(0, 2)) {
    
    if (length(configuration$drainage$alpha) > 0) {
      
      configuration$drainage$alpha <- c()
      
      kwb.utils::catIf(dbg, sprintf("%s was set to an empty vector since %s.\n", 
                          "configuration$drainage$alpha", 
                          "configuration$drainage$idra is 0 or 2."))      
    }
    
    if (configuration$drainage$idra == 2) {
      configuration$solution$isoln <- setUnlessEquals(
        configuration$solution$isoln, "configuration$solution$isoln", 2, 
        reason = "configuration$drainage$idra == 2")
    }
  }
  else if (length(configuration$drainage$alpha) > 5) {
    
    configuration$drainage$alpha <- configuration$drainage$alpha[1:5]
    
    kwb.utils::catIf(dbg, sprintf("%s was reduced to its first five elements.\n", 
                        "configuration$drainage$alpha"))
  }
  
  if (configuration$times$its == 1) {
    reason <- "configuration$times$its == 1"
    
    configuration$times$tlast <- setUnlessEquals(
      configuration$times$tlast, "configuration$times$tlast", 0, reason)
    
    configuration$times$nlc <- setUnlessEquals(
      configuration$times$nlc, "configuration$times$nlc", 0, reason)
    
    configuration$times$nox <- setUnlessEquals(
      configuration$times$nox, "configuration$times$nox", 0, reason)
  }
  
  if (configuration$solution$isoln == 2 && configuration$solution$method != 3) {
    warning(paste("Only solution method 3 has been tested and was found",
                  "to be satisfactory. You should consider to set",
                  "configuration$solution$method to 3."))
  }
  
  numberOfObservationWells <- length(configuration$obswells)
    
  if (numberOfObservationWells > 0) {
    
    maxNumberOfObservationWells <- 24
    if (numberOfObservationWells > maxNumberOfObservationWells) {
      
      # choose 25 nearest wells to the pumping well
      selection <- 1:maxNumberOfObservationWells
      indices <- order(sapply(configuration$obswells, "[[", "r"))[selection]
      configuration$obswells <- configuration$obswells[indices]
      warning(paste("You specified more than 24 observation wells.",
                    "The application of WTAQ-2 is restricted to only 24",
                    "observation wells. The 24 observation wells nearest to",
                    "the pumping well have been chosen for the simulation."))
    }    
    
    for (i in 1:length(configuration$obswells)) {
      observationWell <- configuration$obswells[[i]]
      
      if (observationWell$iows == 2) {
        reason <- "observationWell$iows == 2"

        configuration$obswells[[i]]$z1 <- setUnlessEquals(
          observationWell$z1, 
          sprintf("configuration$obswells[[%d]]$z1", i), 0.0, reason)
      
        configuration$obswells[[i]]$z2 <- setUnlessEquals(
          observationWell$z2, 
          sprintf("configuration$obswells[[%d]]$z2", i), 0.0, reason)        
      }
      else if (observationWell$iows %in% c(0, 1)) {

        configuration$obswells[[i]]$zp <- setUnlessEquals(
          observationWell$zp, 
          sprintf("configuration$obswells[[%d]]$zp", i), 0.0, 
          reason = paste("observationWell$iows ==", observationWell$iows))
      }      
      
      if (observationWell$idpr == 0) {
        reason <- "observationWell$idpr == 0"
        
        configuration$obswells[[i]]$rp <- setUnlessEquals(
          observationWell$rp, 
          sprintf("configuration$obswells[[%d]]$rp", i), 0.0, reason)
        
        configuration$obswells[[i]]$xll <- setUnlessEquals(
          observationWell$xll, 
          sprintf("configuration$obswells[[%d]]$xll", i), 0.0, reason)
      }
    }
  }
    
  configuration
}