Frequently Asked Questions (FAQ)

It should be noted that in a real flow system, a hydrostatic-pressure boundary condition will not be maintained under conditions of variable-density flow. Specification of a hydrostatic-pressure boundary condition using a uniform or other specified density distribution can cause disconcertingly large vertical flows to occur when the boundary density changes due to fluid compressibility or changes in the boundary temperature. Also, boundary-pressure values need to be specified to four or five significant digits to avoid small vertical flows caused by roundoff error.

Lateral boundaries are commonly specified as hydrostatic pressure and a geothermal temperature profile for the associated temperature. At an outflow boundary, a boundary condition of specified hydrostatic pressure and associated temperature distribution can become unrealistic. As the boundary temperature changes, the fluid density changes, so the pressure distribution no longer satisfies hydrostatic equilibrium. This can induce artificial vertical flows along the boundary. One possible remedy would be to extend the simulation domain so the problematic boundary is located away from the region of significant temperature changes. Then the assumption of a hydrostatic pressure boundary condition can be expected to hold during the simulation period.

Along a sloping land surface, with a boundary condition of a seepage surface in combination with precipitation recharge, the contact point between the two boundary conditions is transient and is determined as the simulation progresses. To avoid having to determine and iteratively relocate the contact point, the following approach can be used in HYDROTHERM. The region is first discretized giving a stair-step approximation to the sloping land surface. Then extra columns of thin facing cells are defined at the vertical boundary faces. These cells are specified as seepage boundary condition cells. The boundary cells along the stair-step boundary with land surface boundary faces that are only horizontal are specified as precipitation-recharge boundary cells. This approach removes the need to iteratively determine the contact point between the recharge part and the seepage-surface part of the boundary, and has the restriction that a cell not change type of boundary condition during the simulation.

In general, the direct band solver in the SSOR method should be used for all 1-dimensional and smaller (less than a few thousand nodes) 2-dimensional, vertical slice grids. The GMRES solver should be used for larger 2-dimensional, vertical slice grids, all horizontal slice grids, and all 3-dimensional grids.

The user should test that the convergence tolerance τGMRES is sufficiently small by rerunning the simulation with the tolerance reduced by a couple of orders of magnitude and verifying that the results do not change significantly. A short simulation period can be used to do this. Experience with test problems has shown that a tolerance of 10 ^-6 to 10 ^-12 is necessary to obtain solutions that have converged to 3 or 4 significant digits. The default value for τ_GMRES is 10 ^-10.

The GMRES solver does not seem to be very sensitive to the number of iterations between restarts. The default value of 5 seems to work satisfactorily for most problems. If solver convergence requires more than about 100 iterations, the restart interval should be increased to 10 or 15.

For ILUT preconditioning, the drop tolerance ( τ_ILU) and the maximum number of fill-in elements per row ( k) are required. An acceptable value for the drop tolerance must be determined by trial-and-error, because it is strongly problem dependent and a clearly optimal value may not exist. Limited experience has shown that τ_ILU should be in the range of 10 ^-2 to 10 ^-4 and that k should be in the range of 5 to 15. Limited testing shows that the value of k needs to be changed by about 5 to affect the convergence rate significantly. The default values for τ_ILU and k are 10 ^-3 and 5, respectively.

For ILUK preconditioning, the level-of-fill ( k ) is required. The workload on the preconditioner increases rapidly with k , and it rarely is worth using a high value. A suitable value is usually within the range of 0 to 10. The default value for level-of-fill for ILUK is 5.

If multiple time periods are being simulated with several different sets of output intervals to the plotfiles, the initial PRINT 6 record must have the selection of the visualization software and the format, even if no output to plotfiles is specified for the first time period. No changes to this selection may be made subsequently. Also note that the print interval for each type of plotfile written during the simulation must not be zero during the last simulation period or that plotfile will be deleted. If all printout to a given output file is turned off for the final simulation time period, that output file will be deleted at the end of the simulation.

Note that the units chosen for the distance measurements in the domain may affect the stability of the simulation. This is because HYDROTHERM is a three-dimensional simulator, so, for Cartesian coordinates, the third dimension of the domain, which is perpendicular to the two dimensions accessible with the GUI, is always one distance unit. Thus, a very thin simulation domain may result from using, for example, cm units to describe a field-scale domain. A very thin domain may necessitate substantial experimentation with the settings of the time-step, Newton-Raphson, and iterative-solver control parameters in order to obtain a stable simulation.

Check the value of your smallest print interval. The automatic time-step algorithm may have its maximum time-step length limited by the minimum print interval specified. Otherwise, check that the values of the maximum changes in pressure and enthalpy allowed per time step are reasonable.