Loading ...
Sorry, an error occurred while loading the content.

Discussion archives on groundwater modelling

Expand Messages
  • C. P. Kumar
    Dear Friend, I have been compiling selected extracts on groundwater modelling discussion from various mailing lists for past few months, as given below. I
    Message 1 of 1 , Dec 14, 2004
      Dear Friend,

      I have been compiling selected extracts on groundwater
      modelling discussion from various mailing lists for
      past few months, as given below. I understand that in
      view of random and discontinuous placing, many parts
      may not be easily comprehensible. Still I think that
      it could be helpful to groundwater modellers to some

      C. P. KUMAR
      Scientist 'E1'
      National Institute of Hydrology
      Jal Vigyan Bhawan
      Roorkee - 247667 (Uttaranchal)

      Web Page : http://www.angelfire.com/nh/cpkumar/
      Kumar Links to Hydrology Resources:

      Google Groups in Hydrology:

      Software User Groups in Groundwater Modelling:


      * In my experience, PCG is faster than LMG for small
      models but LMG is faster for large models. Of course,
      for small models, speed generally is not so much of an

      * DAMP is used as a multiplier to scale the head
      change at each iteration, not RELAX. The DAMP
      parameter in the PCG solver package is analogous to
      the ACCL parameter in the SIP package. Although I
      cannot explain precisely what RELAX controls (since
      I'm no mathematician), the typical range for RELAX is
      0.97 to 1.0. In my experience, setting RELAX=0.97 can
      help to achieve convergence, and may speed up
      convergence on a "well-behaved" model, but not always.

      * I'm sure you've heard "All models are wrong but some
      are useful." If you replaced all the wells in a model
      with drains, you would expect to get different
      results. How is this any different? You have
      replaced one package with a different one that is
      based on different assumptions so of course, they get
      different results. If you want to know which one is
      better, examine the assumptions behind the different
      packages and see which assumpptions best match the
      conditions at the particular region to be modelled.
      If none of the available packages do what you need,
      write your own package to do what you need.

      * The newer LPF package is the way that most other
      models work, i.e., you enter Kz and the inter-cell
      conductances are computed internally. The way that
      MODFLOW88/96 did it with a pre-computed vertical
      conductance always struck me as a bit odd; however,
      that has been our happy MODFLOW world for about 20
      years. So now we have a new option that is
      mathematically more correct, but it is still
      disconcerting when the computed heads are different.
      I guess you'd have to ask a lawyer about who would
      win. My sense is that the differences could never be
      explained to a jury. It does present the MODFLOW
      modeling community with a problem, though. The way
      that you did it with MODFLOW-SURFACT is probably the
      way that the USGS should have done it as well - adding
      the input of Kz in the BCF Package.

      * Input of vertical hydraulic conductivity has its
      advantages, however, McDonald and Harbaugh got it
      correct in keeping a fixed thickness of a cell for the
      vertical direction conductance - for
      gravity-segregated vertical equilibrium (GSVE), you
      have a transmissivity that is conductivity times
      saturated thickness which is valid only for the
      horizontal direction and not for the vertical
      direction (in actuality, the vertical direction
      physics changes from GSVE to unsaturated flow which
      would lower the conductivity further and that would be
      closer to the "fixed thickness" answers that do not
      reduce separation distance between 2 cells vertically
      even for consideration of a homogeneous case, but that
      aside). Nor does the document state where or under
      what conditions the use of a "variable thickness for
      vertical conductances" option performs better. That's
      why I was asking to see if anyone knew why this was
      done in the first place.

      * I propose that the node location of a cell that
      contains the water table should be located at the
      water table instead of the center of the cell as is
      standard. Comparisons of this formulation to analytic
      calculations numerically implemented in the WATQ code
      indicate a much improved solution at shallow depths
      compared to the calculations of the current LPF
      formulation. This formulation does not, however, take
      into account the dynamics of drainage from the
      unsaturated region above the water table which can be
      a major influence. I want to repeat Richard Winston's
      point that "All models are wrong". A model is not
      intrinsically good or bad. That judgment depends on
      what the model is used for. What may be better in one
      situation may not be better in another. You have
      highlighted this concept in your discussion of GSVE.
      MODFLOW still exists for two contradictory reasons; 1)
      as the name implies it is modular and hence easily
      extensible to specialized situations. 2) It is widely
      known and accepted that the code has accurately solved
      groundwater flow problems. It is incorrect to extend
      the second point to expect MODFLOW to accurately solve
      all groundwater flow problems. Once again required
      accuracy like good and bad depends on the problem that
      needs to be solved. That MODFLOW produces two results
      for the same geometry and boundary conditions means
      that it approximates the physics of the flow situation
      differently. It proper care is taken, the ability to
      change the approximation of the physics of flow is a
      good thing because the code can be tailored to
      specific problems.

      * Yes, having flexibility in averaging schemes is a
      good thing. The LPF package provides one more
      averaging scheme in the vertical direction. It does a
      harmonic average of the conductances, as opposed to
      the BCF package's harmonic average of vertical
      conductivity times area divided by separation distance
      as provided in the MODFLOW documentation for the
      various fully three-dimensional, and
      quasi-three-dimensional alternatives. However, I am
      guessing that people are encountering drastically
      different and surprising or unexpected results is
      because of the vertical direction treatment.
      Specifically, the cell thickness varies as a function
      of its saturated thickness in the LPF package for
      vertical direction treatment. For a simple example,
      consider a homogeneous case with two cells where the
      cell below is getting drier. The vertical conductance
      quickly rises due to the cell below drying up, to such
      levels that drain and dry up the cell above (which
      also causes the conductance to increase in a cascading
      manner) which would be a very different result than
      with use of any of MODFLOW's BCF options. Like I
      stated earlier, infact, when things dry up, it should
      slow down the conductance due to physics change to
      unsaturated flow, and even a gravity flow assumption
      will not speed it up.

      Yes, there is a dilema in GSVE models for vertical
      direction treatment because GSVE by its definition,
      applies only to horizontal direction flow - I have
      seen this in saltwater-freshwater GSVE models with
      multiple aquifers as well. And there too, when you
      give it specialized treatment, it is due to certain
      physically motivated considerations. Otherwise, the
      way MODFLOW BCF did it (or a use of grid-block cell
      thickness in the LPF package) is the best you can do
      with the VE assumptions being used, since conductance
      may be smaller than saturated levels due to
      unsaturated zone physics - but you don't know by how
      much, unless you simulate unsaturated zone physics.
      Using saturated thickness instead of cell thickness in
      the vertical direction leads to results that are more
      removed from this physics.

      * This is one of the most challenging groundwater
      modelling projects you can find in real-world
      modelling... Ideally, you should use a
      saturated-unsaturated code, such as FEFLOW or Modflow
      surfact or similar. Such codes are able of simulating
      perched water conditions.

      If you are using traditional MODFLOW, and assuming
      saturated conditions, it is a bit more complicated to
      simulate such a scenario. You will probably have some
      dry cell issues (and perhaps some convergence
      problems), but it is possible to simulate the open pit
      crossing multiple layers using drain boundary

      For the upper layers, where seepage faces may occur,
      set the drain elevation slightly above the layer
      bottom elevation (say 0.1 m above cell bottom
      elevation). If you use Visual MODFLOW as the
      interface, that's very easy to do using formulas that
      will assign the correct elevation even in irregular
      layers. In that case, if the head in the pit wall is
      higher than the cell bottom elevation, Modflow's drain
      package will remove the outcroppng groundwater. The
      conductance value (another input to the drain package)
      is tipically a calibration parameter, obtained through
      matching drain removal rates with real values observed
      at the mine site.

      I suggest that you also use a solver that is not too
      aggressive, such as PCG2, and that you play a bit with
      re-wetting parameters and dumping factors until you
      obtain a stable solution. We've done several mine
      dewatering modelling projects using this approach and
      it definetely works for most cases, as long as you
      understand the code and it's limitations.

      * Dry cells in a multi-layer model are not necessarily
      a bad thing unless they do not accurately reflect the
      conditions in the field. Perhaps I am stating the
      obvious, but just for the record, dry cells occur in
      MODFLOW models when the calculated water table is
      below the bottom elevation of the grid cell. Since
      MODFLOW only simulates saturated groundwater flow, it
      does not represent a head value in the unsaturated
      cell, so the cell becomes 'dry' and MODFLOW assigns a
      'dummy' head value usually in the range of -1e-30 (or
      something else easily recognized as a dry cell).
      Although the presence of dry cells in a model is often
      a hassle because they create problems with the
      convergence of the numerical solution, dry cells do
      not necessarily indicate a 'problem' or error with the
      model. In your case, the dry cells don't appear to be
      causing a problem with the solution convergence (or
      atleast, you haven't mentioned this problem) and you
      are getting a good mass balance, which is a good sign
      but doesn't necessarily indicate 'correctness' of the
      solution. If the water levels in the second layer of
      your model match closely with the observed field
      conditions, then the occurance of dry cells in layer 1
      appear to be reasonable.

      However, if the dry cells are occurring where you
      would otherwise expect saturated conditions (based on
      water levels in the field) then you probably need to
      revisit your boundary conditions to make sure they are
      reasonable (conductance values are often the culprit),
      check your initial heads and make sure you are not
      starting the iterations with dry cells, and make sure
      you have enough recharge flux entering the system to
      accommodate flux leaving the system through wells,
      rivers, and other boundary conditions. If all this
      looks good, then Mark Wilsanac gave some very good
      suggestions on how to change the solver and rewetting
      package settings in order to deal with dry cells. -
      Patrick Delaney

      * If the river boundary cells are the only place in
      your model where water may leave the system, then one
      of the most likely causes of the mass balance problem
      is the conductance values you are using for the river
      cells. Check your conductance values to make sure
      they are within reason. The conductance values can be
      reasonably estimated using the following formula:

      C = (L x W x Kv)/B


      L = length of the river in a cell W = width of the
      river in a cell Kv = vertical hydraulic conductivity
      of the riverbed in a cell B = thickness of the
      riverbed in a cell

      Most of the popular graphical interfaces for MODFLOW
      allow you to enter these 'physical' parameters and
      will automatically calculate the conductance values
      for each river boundary cell.

      Another possible cause for mass balance problems could
      be the solver you are using (SIP often produces larger
      mass balance than other solvers like the PCG or LMG)
      or perhaps the convergence criteria for the solver is
      not 'tight enough'.

      Just as an afterthought, the fact that you have more
      water entering the system from the rivers than leaving
      the system through the rivers does not, in itself,
      indicate a problem with the water balance. This
      indicates you have 'losing rivers' because the water
      table adjacent to the river is lower than the river
      stage. The excess water entering the system from the
      rivers may be leaving through other avenues such as
      evapotranspiration, pumping wells, drains, specified
      head head boundaries, etc. - Patrick Delaney

      * Check out if you dont have any dry well cells, or
      any other dry boundary cell. If a cell goes dry during
      simulation, it becomes inactive for MODFLOW, and any
      boundary assigned to that dry cell will be ignored.
      That means that your wells may be inactive (even
      though you assigned them to be pumping all the time).
      That typically produces a noticeble change in the
      budget, but sometimes the heads may not change too
      much, unless you look at specific well cells very

      My suggestion is to first look into the budget for
      changes in well pumping/injection rates. If the total
      rate is different from those you entered into the
      model, then you have a dry well cell. Zoom in all well
      cells to identify which ones have become dry during
      that stress period (look in all model layers).

      * With a model set up with cell-to-cell size changes
      below 50%, there can be instabilities or failure to
      converge. This can be due to other reasons, or cell
      size change coupled to highly variable conditions in a
      transient solution. If a geologic unit pinches out,
      that should be represented primarily by changes to
      cell properties, and secondarily by changes in cell
      dimension (aspect-ratio rules apply vertically as well
      as horizontally). Remember to preserve aspect ratios
      of each cell below 10:1, and below 5:1 to be safe. It
      is acceptable to have a layer with variable properties
      from cell-to-cell. It is rarely acceptable to use
      MODFLOW cell dimensions to closely simulate physical
      variations in layers when they pinch out.

      * I have not used the WHS solver sold by waterloo
      hydrogeologic. while i therefore cannot comment on
      this solver and its claimed superiority to SIP, i feel
      quite sure that SIP has been the best performing
      solver for almost all models i have constructed over
      the past few years, and never had problems with mass
      balance which could not be successfully overcome by
      using SIP. i do recall a paper published in Ground
      Water years ago when WHS proponents at software
      developer Waterloo Hydrogeologic compared their WHS to
      SIP and PCG (among other public domain solvers). i
      felt that sip had been shortchanged and misrepresented
      in their analysis in that they had claimed that sip
      yielded a much greater mass balance error than that
      obtained when using their whs solver. what they had
      failed to do in their analysis, was to accordingly
      decrease the HCLOSE with their corresponding decrease
      in the acceleration parameter. this is the equivalent
      of not closing tightly enough on heads, and therefore,
      ending up with an excess mass balance error. this was
      an improper use of SIP.

      I have found that lowering of the acceleration
      parameter, and corresponding decrease in the HCLOSE
      will yield a very low (acceptable) to zero mass
      balance error, but one must be careful in that the
      acceleration parameter is not set so low that the
      model converges prematurely with insufficient head

      * Check all inputs in the flow model first. Make sure
      your flows are realistic. Make sure your K's and
      recharge are realistic. Refine your grid and use
      dispersivity values that will minimize numerical
      dispersion. Make sure your grid cell size complies
      with appropriate Peclet criteria (as a rule of thumb,
      for finite diferences, Pe should be ~ 2; that means
      your grid size should be roughly 2* your dispersivity
      value, in the area where, for MOC, this can be
      higher). I suggest your read some transport modelling
      books that explain this in detail. Be careful with
      using constant concentrations, they can act as a
      source, but also as sinks (if he concentration in
      neihbouring cells is higher than the ones you
      specified in the constant concentration node).

      * When she is referring to inflow, in plane, and
      outflow, I think what she is doing is looking at the
      color coding for the velocity vectors or pathlines.
      When she sees that there are very few arrows or
      pathlines with a green color (indicating 'inplane' or
      flat horizontal flow in the layer) she gets concerned
      for some reason. Most likely the lack of 'in plane'
      flow has very little to do with anything she is
      concerned about, it just reflects the fact that she
      has slight vertical gradients in the model.

      WRT to the vanishing concentrations, you could also
      suggest her to make sure she didn't input a default
      degradation rate of 0.5 or something very high and
      then forgot about it, or the other possibility is she
      has assigned a sorption coefficient of 1 (mistaking it
      for a retardation coefficient), or perhaps just a very
      large sorption coefficient of 0.01 L/mg which would
      cause the plume not to move at all, thus giving the
      impression that it is disappearing from the model,
      when in fact it hasn't moved at all. Another
      possibility is that she initially assigned
      non-conservative values of sorption and decay
      coefficients when she set up the transport model
      scenario, and she is now trying to modify these values
      in the Setup/Transport Engine dialog instead of in the
      Input mode (I've seen this problem a few times) or
      Visual MODFLOW.

      Finally, the fact that the solution solves in two
      iterations seems to indicate that the total run time
      of the MT3DMS solution may have been left as the
      default of 1 day (perhaps she assumed it automatically
      picks up the end time of the flow simulation). She
      can change this by selecting Run/MT3DMS/Output Times
      and entering the desired Simulation Time.

      * The simplest way of defining a starting
      concentration is to use the Recharge Coverage and
      specify a polygon shape that denotes the contaminant
      concentration to be applied to the model.

      * What to do when your contaminant plume does not
      migrate as expected.

      When using MT3D or RT3D for a contaminant transport
      simulation, if the plume fails to move as expected (or
      does not appear at all in your Visual MODFLOW output),
      here are a few common causes and their solutions.

      Problem #1: Overlay is not active, or masked by
      another overlay

      Solution #1: As with Head Equipotentials, Particles,
      and even Basemaps and other overlays, you must make
      the Concentration overlay active in order to see it.
      Click the F9-Overlay button at the bottom of your
      screen, and ensure the overlay has been checked off
      (to make it active). It is also possible that other
      overlays are masking your concentrations. To resolve
      this, change the Overlay Order setting to User
      Defined, then highlight the Concentration overlay, and
      use the arrow buttons to move it up the list.

      Problem #2: Simulation Output time

      Solution #2: In Visual MODFLOW, the Simulation Output
      time can be different for your flow and transport
      simulations. Check to see that you have defined the
      correct simulation time(s) in the Run Menu, by
      selecting MT3D (or RT3D) / Output - Time Steps.
      Additionally, in the Output Menu of Visual MODFLOW,
      ensure that you have selected the Concentration
      overlay as the active overlay, then select the Time
      button in the left toolbar, to change output times as

      Problem #3: No concentration input assigned, or not
      assigned properly

      Solution #3: Ensure that you have correctly defined at
      least one cell with an Initial Concentration, or a
      transport boundary condition of one of the following
      types: Constant Concentration, Recharge Concentration,
      Evapotranspiration Concentration, or Point Source.
      Please NOTE that Concentration values entered for
      Concentration Observation Wells are used for
      calibration purposes only; they do not contribute mass
      to a transport simulation.

      Problem #4: Inadequate flow gradient

      Solution #4: Before running your transport simulation,
      run just the flow simulation (MODFLOW). Using
      Pathlines or Velocity Vectors, ensure that there is a
      sufficient flow gradient to cause plume migration.

      Problem #5: Incorrect reaction parameters

      Solution #5: Check the reaction options in the Main
      Menu, by clicking on Setup / Numeric Engines /
      Transport. Check that the correct Sorption and
      Reaction options have been selected.

      First order decay rates (dissolved and sorbed phases)
      may have a tremendous impact on the plume mass, since
      during the simulation, some contaminant mass may be
      removed from the model. If the decay rate is too high,
      your plume will show very small concentrations, and/or
      will not be visible at all at later simulation times.

      Decay rates (due to biodegradation or other
      mass-removal processes) can be taken from literature
      values, but this should be done with caution, as this
      parameter is highly site-specific for most compounds.
      A good reference on decay rates used in natural
      attenuation studies can be found in:

      Wiedemeier et al., 1999: Technical Protocol for
      Implementing the Intrinsic Remediation with Long-Term
      Monitoring Option for Natural Attenuation of Dissolved
      Phase Fuel Contamination in Ground Water. Air Force
      Center for Environmental Excellence, Brooks, AFB.

      This document can be downloaded from:


      The decay rate (l, or sometimes called Kd, with units
      of 1/day), is typically obtained from half-life
      values, converted into appropriate units using the
      following relationship:

      l = ln(2)/t1/2

      Where t1/2 = half-life of the compound

      In addition, check that your Kd value is correctly
      defined in the Species Parameters tab. A very large Kd
      value will result in an extremely high retardation
      factor, which can result in lack of contaminant
      movement through your model.

      For typical organic compounds (such as TCE, DCE, PCE,
      BTEX, etc.) where linear, reversible sorption can be
      assumed, retardation will be calculated using the
      following formula:

      Retardation = 1+ (Bulk Density/porosity)*(Kd)

      Where Kd = Partition coefficient

      For hydrophobic organics, Kd can be determined using
      the relationship below:

      Kd = Koc*foc

      Koc - octanol-carbon coefficient
      foc - organic carbon fraction in the aquifer

      A more detailed overview on Kd's can be found in:

      US-EPA, 1999. Understanding variation in Partition
      Coefficient, Kd, values. EPA
      402-R-99-004A&B. US-EPA Office of Air and Radiation.

      This document can be downloaded from:


      * From the Modflow standpoint, water that exits the GW
      system through of a drain cell is accounted for in the
      budget, but is otherwise not simulated. The situation
      for river cells is similar--Modflow does not simulate
      any flow in the river, it only simulates exhange
      between the GW system and each individual river cell.
      So there is no way to route water from a drain cell to
      a river. The Streamflow Routing (STR) Package, on the
      other hand, does simulate flow in surface-water
      streams, accounting for flow rates in the stream and
      calculating stream stage, which is then used as the
      external head in the GW/SW exchange calculation. If
      you want to simulate exchange between drains and
      rivers, you could use the STR Package for both types
      of features.

      * How to import modflow model files to gwvistas 4?

      Select File/New and click the MODFLOW button. On the
      next dialog, click browse to go find the name file
      (*.nam). If this is an older MODFLOW88 dataset, then
      browse to find the *.bas file. Some datasets
      (especially those that were created without the use of
      a preprocessor) can be troublesome.

      * I want to input separately 2 sets of recharge into a
      modflow model: (1) Effective recharge from rainfall,
      and (2) urban leakage (same format as the (1)), which
      will overlap. Since both components are fairly complex
      and over a long time (1200 Stress Periods), I would
      like to keep them as independent as possible (avoiding
      the data compilation option, if you see what I mean).

      It is possible to do what you want using just the
      Recharge package. Here is how you do it.

      Define each recharge distribution for each stress
      period using parameters. You can use one set of
      parameters for the rainfall and another for the urban
      leakage. With each parameter, you can use multiplier
      arrays and zone arrays to define the spatial
      distribution of the recharge. Then, for each stress
      period, use the parameters for that stress period to
      apply recharge from both rainfall and urban leakage.
      MODFLOW-2000 will add the contributions from all the
      parameters you use in a stress period to determine the
      total recharge.

      You should be able to do this with any reasonably
      up-to-date GUI for MODFLOW-2000. If you are using a
      GUI and you can't figure out how to do this, contact
      the GUI developer for assistance. If you aren't using
      a GUI, you can check the input format in the Online
      Guide for MODFLOW-2000 at
      or in the original MODFLOW documentation as modified
      in "Time-varying-parameters.pdf" (distributed with

      * Since I know you are using Modflow VKD (based on
      MF96, not MF2K) and that the wells package is almost
      certainly being used for many abstraction wells, I
      suggest 2 routes:

      1. A FORTRAN utility to combine two recharge files - a
      days work? 2. To use a MF package that you aren't
      already using - but with some manipulation e.g. the
      rivers package. By setting the stage of the river as
      above any elevation in the model (1000m?) and the
      conductance equal to the quantity you want to recharge
      for that SP - with a 1 m difference between stage and
      bottom elevation you will get the appropriate amount
      leaking (recharge) to each cell.

      You can achieve the same result with the drains
      package by specifying 2 drains in each cell. The trick
      is to have one with a positive conductance and the
      other negative, but the same value. Make the
      difference in elevation 1m and the magnitude of the
      two conductance values the quantity you want to
      leak/abstract. If the elevation is outside the range
      of variation in the model, you end up with a constant
      recharge/abstraction from the model.

      For composite/trick boundary conditions of this type I
      suggest drawing figures for them like those in the
      MF88 book, it will help to understand what I'm trying
      to say.

      * Design the shoreline?

      Depending if you are using FEMWATER or MODFLOW, it is
      handled differently. With FEMWATER, you want the
      boundary of your 3D mesh to follow the shoreline. With
      MODFLOW, you are working with a 3D grid. So, you need
      to mark the cells that contain the shoreline with a
      fixed head boundary condition, and the cells outside
      this area will be disabled.

      If you are attempting to model the interaction with
      the sea bed, then you will need to set the top
      elevations of the mesh or grid to match that of the
      sea bed. However, generally you end your model at the

      * MT3D Performance

      What you need most depends on where your bottle-neck
      occurs. You need to monitor CPU and RAM usage on the
      computer during the MT3D run. If you run Windows 2000
      or XP, you can right click (for a right-handed
      pointing device) on the Taskbar (grey bar at bottom of
      screen), then select Task Manager in the menu pop-up.
      Once Windows Task Manager opens, select the tab
      labeled "Performance". "Performance" shows CPU usage
      graphically and memory usage in numeric form. If the
      "CPU Usage History" (top graph) stays maxed at 100%,
      then a CPU speed increase should help. In parallel,
      observe the "Physical Memory" values. If the
      "Available" value becomes a small fraction of the
      "Total", AND the "Commit Charge, Total" exceeds the
      "Physical Memory, Total" then you are likely to have a
      RAM-limited condition. In this case, more RAM would
      help. With 300,000 cells, consider 1 Gbyte or more. If
      both things are happening (unlikely), then upgrading
      CPU and RAM should help. Be aware that upgrading one
      component sometimes results in the other becoming the
      limiting condition. That is why it is "unlikely" both
      limits are happening simultaneously. Usually one tops
      out before the other. If you run Win98... consider
      upgrading the OS.

      * There are various things that you have to do to get
      Surfact running in Vista 3:

      1. You can only have a limited number of packages.
      For example, if you are doing contaminant transport
      modelling then you have to turn some other pacakges
      off, otherwise Surfact will not run e.g. you can turn
      the cell by cell flow packages off (just type 0 into
      the packages number in Model - Modflow - Packages).
      Another pace you can turn packages off is in the
      output control: Turn off the drawdown unit.

      2. You have to modify some things in the Model -
      Modflow options and some in the Modflow - Surfact
      options. These are:

      Model - Modflow - Packages: change modflow version to
      surfact; change solver to PCG4; and untick
      automatically reset package units

      Model - Modflow - BCF package: tick the box to use
      the BCF4 package

      Model - Surfact - Packages: add a unit number (one
      that you have not yet used) and tick the box for BCF4;
      add the unit number to PCG4 and tick the box (use same
      unit number as shown in Model - Modflow- Packages);
      add a unit number and tick box to RSF4 is you want the
      model to simulate surface seepages (and remeber to
      turn this unit on in Model - Surfact - RSF4 package).

      Model - Paths to models: set modflowwin32 option to
      DO NOT USE; put in correct path for surfact in the
      modflow box.

      I can't think of anything else right now - I always
      forget one thing and spend a couple of hours wondering
      why surfact won't run. If you have any difficulties
      it would be helpful to know if the surfact Dos window
      displays and if so, what it says in the window.

      * The most common problem when using SURFACT is
      setting up the program file. Select Model/Paths to
      Models. Change the MODFLOWwin32 Option at the top to
      "do not use". That tells Vistas that you are not
      using one of our model DLLs. Then change the MODFLOW
      program from MFWIN32.DLL to whatever your SURFACT
      program is (usually it is msft.exe).

      * I believe that if you are using the latest version
      of SURFACT (V2.2) that problem 1. below is no longer
      an issue. It used to be that SURFACT could only open
      a limited number of files but that is no longer the
      case. Also, if you upgrade to Vistas Version 4 and
      SURFACT V2.2 the link between the two programs is
      easier to establish. The same holds true for all of
      the other versions of MODFLOW that Vistas supports
      (MODFLOW88, MODFLOW96 in double precision,

      * VS2DI or VS2DT, being a Richards equation-based
      models have certain PROs and CONs when used for
      groundwater recharge estimation:

      Pros: the Richards equation models are the most
      theoretically based and allow representation of the
      flow processes in porous mediums more realistically
      than the water-balance models. They have been well
      tested against field and laboratory experiments and
      proved to provide good correlation with observed

      Cons: the major complicity of the Richards equation,
      comparing to the saturated flow models, is that its
      coefficients are non-linear functions of the pressure
      potential. To approximate these functions, three
      different algebraic equations are usually used (named
      after their authors): Brooks and Corey's, Gardner's,
      and van Genuchten's. Richards equation can be solved
      with a very limited number of boundary conditions. The
      condition for water at the boundary can be specified
      either as the flux of water across the boundary, the
      head at the boundary or as a combination of specified
      head and specified flux. The Richards equation with
      these boundary conditions is a nonlinear partial
      difference equation that has no close-form or
      analytical solution. To solve the Richards equation,
      a finite-difference method of numerical approximation
      is applied in VS2DI. Searching for a best numerical
      method for the Richards equation became a special
      field in hydrologic science, particularly in 1980's,
      however, almost all implemented approaches cannot
      assure that the model will unconditionally converge
      and simulation results will be obtained. It is an art
      to get the model to work correctly (produce results
      for the whole simulation period with a good water
      balance) when you have varying flow boundary
      conditions in your profile. Our experience also proved
      that applications of Richards equation-based models to
      highly heterogeneous soils with variable hydraulic
      properties can be extremely difficult to execute and
      time consuming.

      These days there more and more examples when less
      sophisticated in flow equation but much more advanced
      in terms of weather/plant effects water-balance model
      HELP is used to estimate groundwater recharge.

      * I will suggest you keep the north and south boundary
      as variable head boundary. This I am suggesting you
      treating there is no river, lake etc. If these
      features are present on north or south side of model
      domain, then you have look for another type of
      boundary. Regarding limited observation well data and
      that too is limited to layer1, then it will not be
      possible to do any realistic modelling. If there is no
      scope of generating additional data, i would like to
      suggest that you confine your modelling activity to
      first layer only.

      * What makes a good correlation depends on what your
      expectations are. If you are conducting exploratory
      research, maybe 0.2 isn�t bad. If you are testing a
      known process having some natural variability, 0.6
      might be acceptable. But if you are calibrating
      instruments, maybe 0.9 isn�t good enough.

      If you don�t have any clear expectations, how do you
      tell if a correlation is good? The answer comes in
      three parts.

      1. Value - Square the correlation coefficient (called
      the coefficient of determination or just R-square).
      R-square is the proportion of variation in the
      dependent variable (y) that can be accounted for by
      the independent variable (x). You might be able to
      decide how good your correlation is from a gut feel
      for how much of the variability you wanted your model
      to account for.

      2. Significance - Every calculated correlation
      coefficient is an estimate. The �real� value may be
      somewhat more or somewhat less. You can conduct a
      statistical test to determine if the correlation you
      calculated is different from zero (or some other
      number if it�s relevant). The larger the calculated
      correlation and the greater the number of samples, the
      more likely the correlation will be significantly
      different from zero. For example, a correlation of
      0.59 (R-square of 0.35) would be significantly greater
      than zero based on about 25 samples, but a correlation
      of 0.01 wouldn�t be significant with 250 samples.

      3. Plots - You should always plot the data used to
      calculate a correlation to ensure that the coefficient
      adequately represents the relationship. Specifically,
      the data should be linearly related and free of
      outliers. There are a few other things to look for
      (e.g., hidden trends, autocorrelation) that I won�t go
      into here.

      So, the reviewer who said that �the "R square" value
      of 0.35 has no significance and is just as good or as
      bad as say 0.01� would only be correct if she looked
      at a statistical test of whether the value was
      significantly different from zero. R square does not
      need to be more than any particular value to draw a
      comparison. Remember, too, that �no relationship� may
      also be an important finding.

      * If you are considering GUI's for your project, you
      may want to look at Visual MODFLOW Pro. In
      particular, if cell splitting, or grid refinement is
      an issue for you. Visual MODFLOW includes two
      distinct features that help users with this.

      1) Cell refinement tool
      2) Grid smoothing tool

      1) Cell refinement/splitting tool: You'll need no
      more than a few clicks to edit your grid in a very
      flexible manner. To split your cells, just click on
      'Refine by ...", and specify how many times you want
      the cell to be split into (2, 3, 4...). Then, select
      the area you need to refine with your mouse.You can
      apply this to one row or column, or all cells in the
      domain - whatever you wish. The same applies for
      splitting layers and coarsening the grid.

      2) Grid smoothing tool: After grid refinement, you may
      want to use the grid smoothing routine to produce a
      smooth transition between coarser areas to highly
      discretized ones. Again, just a few clicks. This is
      actually a perfect tool to help you better understand
      the effects grid size and refinement will have on your
      model results.

      Another feature that may be useful for you is Visual
      MODFLOW Pro's batch capabilities. It will allow you
      to prepare all data input files in different model
      projects and run them in 'batch mode'.

      * The interpretation method is different according to
      the geology. You have methods for porous media, others
      for fissured aquifer, others for double porosity
      media. Methods also vary according to the type of
      aquifer (confined or unconfined; semi-confined with a
      leaky aquifer...) and corrections can be made
      according to the type of well (partially penetrtaing
      well; pumping well with head loss; skin effect and so
      on). Finally, boundary effects (barrier or recharge)
      may affect significantly your interpretation.

      In order to remain pragmatic I would suggest you to
      start with the simplest method : the Theis formula,
      (or even the Jacob approximation if your pumping time
      is big enough) and check wether you get a good fit or
      not with realistic parameters. Even a fissured
      unconfined aquifer can be interpreted with Theis when
      your radius of influence is wide so that the fissure
      network can be assimilated to an homogeneous
      equivallent porous media and the depletion of the
      water level is moderate compared to the thickness of
      the aquifer. This should give you a rough estimate of
      your transmissivity, which is in many cases,

      If you need more accuracy you can use specific
      softwares. You will find on the net many softwares
      proposing many different methods. You don't need Pest
      as they have their own fitting engine. AQTESOLV for
      instance proposes methods with well storage effect.
      However, as soon as a software proposes its own
      adjustement, it's becoming very difficult to know
      what's you're doing and you end up playing a video
      game. So be carefull, especially if you lack
      experience in this domain, not to use complex methods
      with many parameters which will end in a perfect fit
      but with irrealistic results. In my company (French
      geological survey)we use ISAPE, a software which is no
      more commercialised, but whose philosophy is to let
      the user propose realistic values and where you can
      add or remove well effects or boundary influence. I
      wish this approach were more widespread.

      * I think your scenario will represent two different
      geologic materials with peculiar intrinsic properties.
      Hence, I suppose this will require two separate
      simulations. To the best of my knowledge, MODFLOW 2000
      only allow one K matrix values per layer per
      simulation and this is incorporated within the flow
      package. However, if the introduction of the permeable
      barrier is controlled along defined specific paths or
      unique relatively smaller area(s), then you can use
      one simulation with two stress periods and incorporate
      Horizontal Flow Barrier Package in the second period
      to account for the permeable barrier distribution.
      This is not suppose to be of much problem especially
      if you incorporate your modflow with GIS package like
      GRASS GIS.

      * Sophisticated modeling programs like Visual Modflow
      Pro, or any of a number of similar packages, will run
      properly with a variety of inputs that may not be
      physically realistic - and thus provide garbage for
      output. With that in mind, here are some thoughts on
      your questions:

      1) You may divide a single aquifer into more than one
      layer if there is some reason to look at vertical
      stratification within the aquifer. For example, a
      thick aquifer may have different hydraulic or
      transport properties in different vertical zones. You
      may also divide a single aquifer into several layers
      if you want to examine the vertical flow patterns due
      to pumping from different zones within the aquifer.

      2) The top of Layer 1 could be the ground surface,
      the water table (not as common because it can move
      vertically), or the top of the uppermost aquifer or
      aquitard of interest depending on the physical system
      being modeled.

      3) You'll need to assign boundary conditions no
      matter what the size of the physical system is. If
      you don't have natural physical boundaries, one
      strategy is to assign boundary conditions at a
      distance far enough away from the area of interest to
      minimize the effects of boundary assumptions on the

      * It is appropriate to use a single layer if you feel
      that it is reasonable to ignore vertical flow in your
      model. Assuming that you would be using the Layer
      Property Flow (LPF) Package of Modflow-2000, and if
      the aquifer represented by layer 1 is not overlain by
      a confining unit, then specifying the top of layer 1
      as ground surface is reasonable. The reason that you
      would be better off this way than specifying the water
      table elevation as the layer top is that the LPF
      Package only provides two options with respect to the
      confined/unconfined status of a layer: Confined or
      Convertible. If the head in a cell goes above the
      specified top of a convertible layer, the cell is
      treated as confined. Every model needs to have at
      least one fixed head, either as a constant-head cell
      or as a fixed external head (such as the stage at a
      river cell), or the solver will not be able to reach a
      solution. You may need to expand the domain of your
      problem to encompass a fixed head.

      * 1- Model design depends greatly on your project
      goals. For example, you could have one aquifer, but an
      overlying layer (and/or underlying layer) as well. Or,
      you could be modelling just 1 aquifer, but split it up
      into several layers in VMOD to provide more detailed
      flow information (i.e. if the aquifer is 100 feet
      thick, having just 1 layer will not provide very
      detailed output. You can subdivide the 1 "real-world"
      layer into 10 "model" layers, all with the same
      hydrological properties, in order to obtain more
      detailed output). Or, you can simply design a 1-layer
      model if that is all your project requires.

      2- In Visual MODFLOW, Layer 1 represents the uppermost
      layer of your profile, and the top of Layer 1 is
      therefore the ground surface. Layers in Visual MODFLOW
      represent soil layers (or conceptual soil
      layers)....the location of groundwater in these layers
      is determined by your model inputs, and the resulting
      calculations of the flow model.

      3- The assignment of boundary conditions is a question
      that often goes beyond the scope of Technical Support,
      and enters the realm of Extended Modeling Support.
      However, to provide you with some guidelines:

      -If you have no "obvious" water inputs (streams,
      rivers, lakes, injection wells, recharge from
      rainfall, etc.), you can try assigning an upgradient
      equipotential line, and downgradient
      equipotentialline, using Constant Head Cells at the
      upgradient and downgradient edges of your model, to
      represent the regional water table (other boundaries
      may be more appropriate, depending on your specific
      model domain and objectives). You can then add
      additional inputs (extraction wells, contaminant
      sources, etc.) to the model region, which will be
      influenced by your regional flow gradient.

      4- I would recommend taking some time to run through
      the Tutorials included with your Visual MODFLOW
      software (the PDF instruction files are located on
      your installation CD-Rom, and the files are
      automatically copied to the Tutorial subfolder in your
      Visual MODFLOW program folder). These tutorials are
      designed to teach you how to work with your software,
      and how to design a model, import data for model
      inputs, run the various packages, calibrate your
      model, examine your output in 2D and 3D, and export

      * You can use the field interpolator program in the
      PMWIN package. Using the field interpolator you can
      interpoltae the field data in txt file to the grid
      point of your model. On the other hand, you can import
      any dxf file which illustrate the changes in hydraulic
      properties to help you in the graphical interface.

      * Upwind discretization always has an artificial
      diffusion term no mater how refined your grid is (even
      thought it diminishes with grid refinement). Central
      discretization always produces an artificial
      dispersion term that causes "wigles" however you can't
      control it with grid refinement. For almost any Pellet
      number the upwind scheme will produce positive
      coefficients for your coefficient matrix, while that
      isn't true for central differences. However central
      differences has "second order accuracy" while upwind
      only has first order. You can visualize numerical
      diffusion and dispersion mathematically by deriving
      the "equivalent or modified equation" see
      for some cool notes on this (as good as any book I've

      Physically you can visualize numerical diffusion by
      imagining a 1D discretization of a domain | : | : | :
      | where | are faces and : is the center of the cells.

      If your initial property is 1 on a given cell say i
      and 0 on all others it would look like this:

      | 0 | 0 | 0 | 0 | 1 | 0 | 0 | 0 | 0

      So using upwind with a courant = v dt / dx of say 0.5
      and a velocity from left to right after the first time
      step the accurate solution would be something like
      (zooming on cell i):

      | 0 { 0.5 | 0.5 } 0 |

      Note that I have created a new cell delimited by the
      {} brackets this "new cell" that doesn't exist in our
      domain is exactly half of cell i and alf of cell i+1.
      The property's value will be 0.5 inside this new cell
      and 0 outside it (meaning that cell i and i+1 would
      have half with 0.5 and half with 0). This would be the
      accurate solution because a courant of 0.5 makes the
      property shift half a cell every time step. However we
      are bounded to our initial discretization so the
      solution we will obtain will be:

      | 0 | 0.5 | 0.5 |

      so you can see than instead of having a 0.5 of
      property in a volume the size of a cell we have 0.5 in
      a volume the size of 2 cell thus lowering the
      concentration. This happened because the definition of
      concentration is C = lim vol ->0 DM/DVOL and we
      implemented it as DM/DVOL. If we would diminish the
      size of our grid so that currant = 1 we could obtain
      the exact solution. So the smaller your grid is the
      more accurate your solution will be. This is the
      physically correct way of solving the diffusion
      problem another way is to use higher order schemes
      witch usually maintain the sharpness of the gradients
      but usually remove mass from all the wrong places
      (even thought they conserve it). Looking at the same
      example a central differences discretization would
      derive something like:

      | -0.25 | 0.75 | 0.25 |

      Cell i-1 will produce a negative property because the
      concentration al the face between i-1 and i is 0.5
      cell i +1 will obtain the accurate concentration of
      0.25 (0.5 in half the volume of a cell). But we are
      producing the famous "wigles" close to the properties
      gradient. There are smarter second order schemes like
      qwick or TVD schemes.

      * In visual MODFLOW there is provision for importing
      MODFLOW files (*.bas files). In PMWin there is
      provision for conversion of MODFLOW88/96 (*.nam). I
      have tried to import *.bas file in Visual MODFLOW and
      found that model data is not correctly coming to
      model. These options may be available in GMS but I am
      unable to locate the exact sub-menu. Number of options
      is available for opening different format files In the
      FILE menu. But It does not permit opening of *.bas or
      *.nam files. This option must be in different menu.


      Do you Yahoo!?
      Dress up your holiday email, Hollywood style. Learn more.
    Your message has been successfully submitted and will be delivered to recipients shortly.