  In this paper we introduce adaptive time step control for simulation of
evolution of ice sheets. The discretization error in the approximations is
estimated using "Milne's device" by comparing the result from two different
methods in a predictor-corrector pair. Using a predictor-corrector pair the
expensive part of the procedure, the solution of the velocity and pressure
equations, is performed only once per time step and an estimate of the local
error is easily obtained. The stability of the numerical solution is maintained
and the accuracy is controlled by keeping the local error below a given
threshold using PI-control. Depending on the threshold, the time step $\Delta
t$ is bound by stability requirements or accuracy requirements. Our method
takes a shorter $\Delta t$ than an implicit method but with less work in each
time step and the solver is simpler. The method is analyzed theoretically with
respect to stability and applied to the simulation of a 2D ice slab and a 3D
circular ice sheet. %The automatically chosen $\Delta t$ is either restricted
by accuracy or stability depedning on an error tolerance. The stability bounds
in the experiments are explained by and agree well with the theoretical
results.
