The melted ice layer between epochs 2009 and 2016 of a part of the Vernagtferner glacier located in the Austrian Alps is modelled dynamically. Three deterministic and two stochastic approaches are applied to compute the induced gravity signal differences. For the deterministic approach, the applied algorithms include the analytical method of right rectangular prism, the line integral analytical solution of general polyhedron and the numerical solution of fully normalized spherical harmonics series expansion. For the stochastic approach, the glacier is approximated as an uncertain general polyhedron of variable shape and as a summation of smaller variable polyhedral masses, where the melted body is expressed in the form of a covariance matrix. The implemented stochastic algorithm, based on the variance propagation law, uses the partial derivatives of the corresponding gravitational functionals with respect to the polyhedral vertex coordinates to derive gravity signal variations implied by the shape changes as described in the covariance matrix. The numerical tests were performed on a set of 20 points which are in fact observation sites of a local gravity survey network. The relative differences of the estimated variations considering the glacier as one polyhedron from the line integral analytical solution are below 7% for gravitational potential, up to 50% for its first and up to 86% for its second order derivatives. These differences between stochastic approach using individual polyhedral elements and analytical solution of a general polyhedron range up to 0.25%, 40% and 60% respectively. Prismatic representation derives relative differences with respect to the analytical general polyhedron up to 53%, 23% and 50% for gravitational potential, its first and second order derivatives, while the differences of harmonic series from the analytical polyhedral solution are 0.35%, 18% and 82% respectively.