ASIS:Solver temperature

Temperature solver overview
(Tony P initial description for distribution to interested parties)

The temperature (T) is a fairly standard diffusion-advection equation with the addition of a source term that incorporates heat generated by internal dissipation (

$$\Phi $$, the product of stress and strain-rate vectors). The temperature evolves according to (here in z-coordinates for clarity; -coordinates are used in the model)

$$\frac{\partial T}{\partial t}=\frac{k}{\rho c}\frac{\partial ^{2}T}{\partial z^{2}}-w\frac{\partial T}{\partial z}-u\nabla t+\frac{\Phi }{\rho c}$$

where w is vertical velocity, k the conductivity of ice and c its specific heat capacity. Scaling arguments justify the omission of horizontal diffusion. Vertical velocity is found by integrating horizontal velocity divergence from the ice sheet base to its upper surface (i.e, the mass continuity equation).

Boundary conditions are: prescribed mean-annual air temperature at the upper surface (Dirichlet) and, at the base, either the ocean-water temperature (in the case of floating ice) or the sum geothermal heat flux and friction generated by sliding over the bed (von Neumann). The ice sheet bed may reach melting point.

The code solves (7) using an implicit calculation in the vertical with the horizontal advection term coming in explicitly (using second-order upwinding). The implicit, vertical code is iterated over each horizontal grid point until temperatures converge. If, at the bed, temperatures exceed the pressure melting point, a melt rate is calculated and melt water is routed into a drainage system (see below). Conversely, if the ice wants to freeze but there is melt water present, a freeze-on rate is calculated and temperatures remain at melting point until the melt water reservoir is exhausted