0

Problem initializing water table at a given elevation in eos9

init init.out

Hello,

In the TOUGH2 manual, it says that an eos9 (Richard's Equation) model can be initialized with a specified water table elevation by setting CWET in REFCO to that elevation. However, this function seems to not be working in the version of iTOUGH2 that I am using. Would you mind helping me to figure out why this might be?

Attached are the input and output files from the initialization simulation that I'm trying to run. Because I set CWET to -6. in REFCO, the results should show a transition from saturated to unsaturated conditions at this elevation, but they do not. In fact, I tried changing CWET to -12., and the exact same result was achieved, meaning that this input value is effectively being ignored by the code.

Please let me know if you can help me to resolve this issue.

Thank you very much,

Mike

8 replies

null
    • Finsterle GeoConsulting
    • Stefan_Finsterle
    • 5 yrs ago
    • Reported - view

    Mike,

    Sorry, this option is simply not implemented for n iTOUGH2. I could certainly do it if there is a pressing need, i.e., you cannot use one of the other ways to accomplish the same result.

    Cheers,

    Stefan

      • Michelle_Plampin
      • 5 yrs ago
      • Reported - view

      Stefan,

      Thank you so much for your quick reply. It is a relief to know that it wasn't just me doing something wrong!

      If you could incorporate this functionality into iTOUGH2, it would be tremendously helpful to me. I am investigating recharge, for which I need a coupled model of the saturated and unsaturated zones that is initially in gravity-capillary equilibrium. I was able to create these initial conditions by hand for a 1-D model, but it will become significantly more laborious to do so in 2-D or 3-D, which will be a requirement of my project somewhat soon.

      Please let me know if you are willing to incorporate the water table elevation input capability into iTOUGH2.

      Thank you very much again,

      Mike

      • Mikey_Hannon
      • 5 yrs ago
      • Reported - view

      Mike,

      I could be completely off-base here, but could you set the pressure at the bottom of your mesh equal to the pressure at that elevation (p_atm + rho*g*h, where h is the depth of the water column)?  I believe that should do what you need.  If you're using a simple XYZ grid, this can be done with the AddBound executable that Stefan developed.  You would create and additional layer of elements underneath the bottom-most flow-domain layer, and AddBound can automatically search through those and convert them all to a single element, whose volume you can set large and whose distance to the bottom interface you can set small, so that the pressure is fixed throughout the course of the simulation.

      If you think that would do it, I (or Stefan) could connect you with that executable, if you don't have it already.

      Thanks!

      - Mikey

    • Michelle_Plampin
    • 5 yrs ago
    • Reported - view

    Mikey,

    Thank you for the advice, but adding a constant pressure boundary condition at the bottom of my model would not be appropriate in this case, because the goal of this project is to investigate changes in water table elevation (i.e., pressure changes throughout the domain) due to changes in flux at the top (infiltration). If I set a constant pressure at the bottom, the head will not change there and I won't get the output that I need.

    Mike

    • Michelle_Plampin
    • 5 yrs ago
    • Reported - view

    Mikey and Stefan,

    On second thought, I suppose I could use the constant pressure boundary solution that Mikey suggested, just to get the initial condition, then remove the boundary node when I go to start a transient solution. Sorry that I didn't think of this before.

    Mike

      • Mikey_Hannon
      • 5 yrs ago
      • Reported - view

      Makes sense to me!  Let me know if you'll need any help setting that up.

      • Michelle_Plampin
      • 5 yrs ago
      • Reported - view

      Mikey,

      I tried setting up the initialization with a single boundary node at the bottom of a 1-D model as you suggested, but it's not working for me. Since I am mainly interested in the unsaturated zone, I removed all  the nodes that would be under the water table and replaced them with a single boundary node at the bottom. I then set the pressure in this node to atmospheric and the final time of the simulation to a large number (1e12) so that it should run to steady state.

      No matter what I try, the simulation fails to converge right away. I have tried various "default/global" initial conditions in PARAM, including atmospheric pressure, or water saturation of 0.5, but I always seem to get the same result.

      I have attached my most recent input and output files. Any help would be greatly appreciated.

      Thanks,

      Mike

      output.txt input.txt
      • Michelle_Plampin
      • 5 yrs ago
      • Reported - view

      Never mind, sorry again. I forgot to change the - to a + when I changed the node volume. Looks like making that change fixed my problem

Login to reply

Content aside

  • 5 yrs agoLast active
  • 8Replies
  • 81Views
  • 3 Following