In computational models of geomorphologic processes Stokes equations have been employed. Multigrid method is added to our solver so that we may construct a solution for these problems. When trying to solve the ellipse differential equation with ill-conditioned matrix, the interpolation method is frequently used to start reducing the iterative and incremental stages due to the point sets inside the matrix merging this same momentum as well as mass formulae and the firmly customisable viscosity due to rheology. Utilizing advantages of the Matlab and the capabilities of the available Graphic Processing Units (GPUs), we accelerate the original Matlab routines using the Compute Unified Device Model (CUDA).