## References

Bayrock, L. A
1967. Catastrophic advance of the Steele Glacier, Yukon, Canada. Alberta University. Boreal Institute. Occasional Publication No. 3.

Carnahan, B, *and others*. 1969. Applied numerical methods, by B. Carnahan, H. A. Luther and J. O. Wilkes. New York, John Wiley and Sons, Inc.

Carslaw, H. S, *and*
Jaeger, J. C
1959. Conduction of heat in solids. Second edition. Oxford, Clarendon Press.

Clarke, G. K. C, *and*
Goodman, R In press. Radio echo soundings and ice temperature measurements in a surge-type glacier. Journal of Glaciology.

Classen, D. F, *and*
Clarke, G. K. C
1971. Basal hot spot on a surge type glacier. Nature, vol. 229, No. 5285, p. 481–83.

Classen, D. F, *and*
Clarke, G. K. C
1972. Thermal drilling and ice-temperature measurements in the Rusty Glacier. (*In*
Bushnell, V. C, and Ragle, R. H, ed. Icefield Ranges Research Project. Scientific results. Vol. 3. New York, American Geographical Society; Montreal, Arctic Institute of North America, p. 103– 16.)

Forsythe, G. E, *and*
Wasow, W. R
1960. Finite difference methods for partial differential equations. New York, John Wiley and Sons, Inc.

Jarvis, G. T Unpublished. Thermal studies related to surging glaciers. [M.Sc. thesis, University of British Columbia, 1973.]

Jarvis, G. T, *and*
Clarke, G. K. C Unpublished. The thermal regime of Trapridge Glacier and its relevance to glacier surging.

Nielsen, L. E
1969. The ice-dam, powder-flow theory of glacier surges. Canadian Journal of Earth Sciences, Vol. 6, No. 4, Pt. 2, p . 955–61.

Peaceman, D. W, *and*
Rachford, H. H Jr., 1955. The numerical solution of parabolic and elliptic differential equations. Journal of the Society for Industrial and Applied Mathematics, Vol. 3, No. 1, p. 28–41.

Sharp, R. P
1951. The glacial history of Wolf Creek, St. Elias Range, Canada. Journal of Geology, Vol. 59, No. 2, p. 97–117.

Stanley, A. D
1969. Observations of the surge of Steele Glacier, Yukon Territory, Canada. Canadian Journal of Earth Sciences, Vol. 6, No. 4, Pt. 2, p. 819–30.

Thomson, S
1972. Movement observations on the terminus areas of the Steele Glacier, Yukon, July 1967. (*In*
Bushnell, V . C
*and*
Ragle, R. H, ed. Icefield Ranges Research Project. Scientific results. Vol. 3. New York, American Geographical Society; Montreal, Arctic Institute of North America, p. 29–37.)

Weertman, J
1971. Theory of water-filled crevasses in glaciers applied to vertical magma transport beneath oceanic ridges. Journal of Geophysical Research, Vol. 76, No. 5, p. 1171–83.

Wood, W. A
1936. The Wood Yukon expedition of 1935: an experiment in photographic mapping. Geographical Review, Vol. 26, No. 2, p. 228–46.

Wood, W. A
1967[a]. Glaciology: chaos in nature. Explorers Journal, Vol. 45, No. 2, p. 79–87.

Wood, W. A
1967[b]. Steele Glacier surge. American Alpine Journal, Vol. 15, No. 2, p. 279–81.

Wood, W. A
1972. Steele Glacier, 1935–1968. (*In*
Bushnell, V. C, *and*
Ragle, R. H, ed. Icefield Ranges Research Project. Scientific results. Vol. 3. New York, American Geographical Society; Montreal. Arctic Institute of North America, p. 1–8.)

## Appendix A Freezing of a cylindrical water-filled hole in cold ice

To compute hole closure rates and cooling curves for a water-filled cylindrical hole in cold ice, diffusion equations of the form

must be solved in ice and water. At the ice–water interface the boundary conditions are

and

where *r*
_{c} is the radius of the water-filled cylinder and *T*
_{m} the melting temperature of ice. The remaining boundary conditions are

and

where *T*
_{0} is the ambient ice temperature, *q*
_{1} (*t*) the strength of a line heating source, and *r*
_{0} its radius. Boundary condition (A5) allows the possibility of evaluating the effect of ohmic dissipation in the power cable leading to the thermal probe. In the calculations discussed above *q*
_{1} was negligible and the water phase was essentially isothermal at temperature *T*
_{m}. Solutions for large values of *q*
_{1} have also been computed to determine whether line heating can be used to inhibit hole closure during thermal drilling.

In passing to a finite-difference approximation of (A1) it is convenient to introduce a logarithmic grid by the transformation *R*= In *r*; thus (A1) becomes

with *T* = *T*(*R, t*), and (A2) gives

where *R*
_{c} = ln *r*
_{c}.

Following the Crank–Nicolson approach, the solutions of Equation (A6) for times *t* and *t* + *τ* are averaged to reduce the discretization error giving as the finite-difference equation in the *i*th medium

where *h* is the space increment of the logarithmic grid, *τ* the time increment, *λ*
_{i} = *κ*_{i}
*τ*/*h*^{2}
, and *θ* is an averaging parameter which is usually set to the value *θ* = 0.5. The variables *R* and *t* take discrete values *mh* and *nτ* respectively, where *m* and *n* are integers. In Equation (A8) the right-hand side terms arc known and the left-hand side terms are unknown. If similar equations are written at each grid point one obtains a tridiagonal set of linear equations which can be solved for the unknown temperatures *T*(*R*, *t* + *τ*).

Infinitely large grids are not feasible so that the boundary condition (A4) is replaced by *T*(*R*
_{max}, *t*) = *T*
_{
0
} where *R*
_{max} is some suitably large value of *R* =In *r*. At *R*
_{
0
}
*=* In *r*
_{0} we have the condition

and at *R*
_{c}
*=* In *r*
_{c}, the ice–water interface, *T* (*R*
_{c}, *t*) = *T*
_{m}. The finite-difference equations (A8) are solved subject to the above boundary conditions and the migration of the ice–water interface is evaluated at each time step by substituting finite-difference approximations of *•T*(*R*
_{c}, *t*)/*•t* and *•T*
_{w}(*R*
_{c}, *t*)/*•t* into Equation (A2). When the condition *R*
_{c}< *R*
_{0} is satisfied, the water phase is considered to vanish and a simple one-phase problem results.

## Appendix B Peaceman–Rachford numerical method

For a finite-difference grid with space intervals Δ*x*, Δ*y* and time step *τ*, the standard implicit finite-difference approximation to Equation (2) yields

where *λ*_{x}
=*kτ*/(Δ*x*)^{2}, *λ*_{y}
=*kτ*/(Δ*y*)^{2}, and the variables *x*, *y*, and *t* have the discrete values *x* = *i*Δ*x*, *y* = *j*Δ*y*, and *t = nτ* for integer values of *i*, *j*,and *n*. Equations (B1) are implicit in both *x*- and *y*-directions and have five unknowns per equation. Direct solution of this system of equations requires the inversion of a large five-band diagonal matrix and is computationally expensive.

In the Peaceman–Rachford implicit-alternating-dircction method, two systems of equations are used in turn over successive time steps of duration *τ/2*. The first equation is implicit in the *x*-direction only, the second in the *y*-direction. Using the notation *T**(*x,y*) to represent the intermediate values of *T* half-way through the time step *τ*, for implicit *x* we have

and for implicit *y*

The systems of Equations (B2) and (B3) have only three unknowns per equation and the implicit solution of each system merely involves the inversion of tridiagonal matrices for which simple and efficient algorithms are readily available.

Holding *y* constant, one equation of the form (B2) is written for each value of *x* and the resultant tridiagonal system of equations is solved simultaneously. Equations (B2) are solved in this manner once for each value of *y* to generate the complete solution *T**(*x,y*). Equations (B3) are now solved by substituting the solution *T**(*x,y*) obtained from (B2) into (B3). Holding *x* constant, one equation of the form (B3) can be written for each value of *y* and the new system of equations solved simultaneously. Equations (B3) are solved once for each value of *x* to generate *T*(*x*, *y*, *t* + *τ*), the temperature distribution advanced one full time step. This procedure is unconditionally stable for any value of τ and the discretization error is 0[*τ*
^{2}+ (Δ*x*)^{2}].