This paper deals with the development of accurate surrogate models for first-principles models constructed for the dissolution of phosphate ore in a phosphoric acid solution. The surrogate models are based on sparse multivariate polynomial interpolation. The main goal is to reduce the computational time of the first-principles model while preserving its properties, namely the outputs’ monotonicity and positivity. The temporal profiles of the concentrations of the different components involved in the liquid phase and in the ore particles, the particle size, and the thickness of the liquid film surrounding the particles are approximated by the surrogate models. The inputs of the models are the particle size distribution, the initial acid concentration, and the hydrodynamic conditions. A design of experiments method is used to produce the required sampling points for the surrogate models and several simulations are conducted in the MATLAB environment. A final comparison of the proposed surrogate models’ performance demonstrates high efficiency, both in terms of accuracy as well as computation time, and highlights the strength of the introduced surrogate modeling methodology.