# Overview

Atlassian Sourcetree is a free Git and Mercurial client for Windows.

Atlassian Sourcetree is a free Git and Mercurial client for Mac.

This repository contains the unconfined (moving water table - aka

Neuman 1972) slug test solution given in

Malama, B., K.L. Kuhlman, W. Barrash, M. Cardiff & M. Thoma,

2011. "Modeling slug tests in unconfined aquifers taking into account

water table kinematics, wellbore skin and inertial effects", *Journal
of Hydrology*, 408(1–2), 113–126.

http://dx.doi.org/10.1016/j.jhydrol.2011.07.028

which is updated and applied to cross-hole slug tests in

Malama, B., K.L. Kuhlman, R. Brauchler & P. Bayer, 2016. "Modeling

cross-hole slug tests in an unconfined aquifer", *Journal of
Hydrology*, 540:784–796.

http://dx.doi.org/10.1016/j.jhydrol.2016.06.060

# 2015 Update

The model has been modified since publication in 2011, to

compute the response at an observation location outside the source

well (i.e., cross-hole slug tests), to use the finite Hankel transform

for inversion, include effective source/observation well skin, and to

include observation well storage and inertial effects. The current

version of the code still works for its original purpose as a

traditional slug-test analyzer.

Cross-hole slug test data from Widen, Switzerland and the associated

parameter estimation results (PEST and DREAM) are located in the

"swiss-data" directory of the repository.

# Code License

The code is available under the MIT open-source license, which allows

free use and redistribution; please acknowledge the original authors

by citing our papers.

# Compilation Directions

The code requires gfortran (>= 4.6) compiler to properly compile (one

that includes libquadmath), but I have gotten it to work with some

minor tweaks using the latest version of the Intel Fortran compiler

(which supports quad precision, but not for the complex hyperbolic

trig functions -- so they must be implemented in Fortran, rather than

in the compiler library).

The code can be compiled, given you have the correct compiler set up

already, by calling the compiler-specific makefile:

make -f Makefile.gfortran debug_driver

or

make -f Makefile.gfortran driver

these two makefile targets build the debug target (slower, checks for

bounds overflows, etc) and the standard target (builds using OpenMP

for a parallel build which is optimized and much faster.

# Input Instructions

The input.dat file has all the input parameters read by the 2011

program. The finite-Hankel cross-hole slug test program (2015) uses

input-finite.dat. Parameters include physical properties of the

formation, well geometry parameters, numerical convergence parameters

and the locations/times to compute the solution at. Anything to the

right is a comment and is ignored by the program when reading in the

input. Units are arbitrary but consistent.

For the cross-hole version of the program, consult

"cross-hole-slugtest-input-explanation.txt" for discussion about each

variable appearing in the input file (input-finite.dat).

# Numerical Approach

The numerical approach used to invert the Laplace- and **infinite Hankel
transform** numerically is a slightly modified and updated version of

the approach documented in Appendix B (repeated below citation) of:

Malama, B., K.L. Kuhlman, and A. Revil, 2009. "Theory of transient

streaming potentials associated with axial-symmetric flow in unconfined

aquifers", Geophysical Journal International, 179(2), 990-1003.

http://dx.doi.org/10.1111/j.1365-246X.2009.04336.x

The numerical approach used to invert the double Laplace-Hankel

inversion is an improvement on the method in Malama et al. (2009). In

this study the Laplace transform is now inverted using the Stehfest

algorithm (Stehfest 1970) and a portion of the Hankel integral is

integrated using an accelerated form of tanh-sinh quadrature (Takahasi

& Mori 1974). The infinite integral in eq. (A2) for inverting the

Hankel transform is evaluated in two steps (Wieder 1999). One finite

(0 <= a <= j_{0,1}), the other approaching infinite length (j_{0,1} <= a <=

j_{0,m}), where j_{0,n} is the nth zero of J_0(ar_D) and m is typically

10-12. The finite portion is evaluated using tanh-sinh quadrature

Takahasi & Mori (1974) where each level k of abscissa-weight pairs uses

a step size h = 4/2k and a number of abscissa N = 2k , similar to the

approach of Bailey et al. (2000). Evaluating the function at Carissa

related to level k, the function evaluations at lower levels k-j

are simply the 2^{k-j} th terms in the series (only requiring the

weights to be computed for each level, reusing the function

evaluations from the largest k). Richardson’s deferred approach to the

limit (Press et al. 2007) is used to extrapolate the approximate

solutions to the limit as h -> 0 using a polynomial extrapolation. The

infinite portion of the Hankel integral is integrated between

successive j_{0,n} using Gauss-Lobatto quadrature and these alternating

sums are accelerated using Wynn’s epsilon algorithm.

The **finite Hankel transform** used in the 2016 paper is a simple sum, and

is therefore much simpler than the above approach.

# Contact Code Author

If you have any questions or require assistance using the code, please

contact me (Kris Kuhlman) at klkuhlm at sandia dot gov. I will

try to help you. I mostly work in Linux and on a Mac, but I have

gotten codes to work in Windows before.

July 2016