A shape optimization algorithm is presented that estimates the ice thickness distribution within a three-dimensional, shallow glacier, given a transient surface geometry and a mass-balance distribution. The approach is based on the minimization of the surface topography misfit in the shallow ice approximation by means of a primal-dual procedure. The method's essential novelty is that it uses surface topography and mass-balance data only within the context of a time-dependent problem with evolving surface topography. Moreover, the algorithm is capable of computing some of the model parameters concurrently with the ice thickness distribution. The method is validated on synthetic and real-world data, where the choice of its Tikhonov regularization parameter by means of an L-curve criterion is discussed. (C) 2014 Elsevier Ltd. All rights reserved.