The accurate and efficient simulation of granular flows, made of complex particle distributions often immersed in a fluid, remains a challenge in the mathematical and computational realm. Historically, different numerical techniques have been developed to deal with this problem either in nonlinear solid mechanics or fluid mechanics. The purpose of this work is to explore a possible path to conciliate these two fields. For this aim, we have developed a displacement-based formulation for free surface viscolastic fluids in a nonlinear solid mechanics framework considering both a weakly-compressible and an incompressible regime. However, the standard Galerkin displacement-based formulation fails in the fully incompressible case, due to the fact that the formulation is unable to evaluate correctly the strain field. For this reason, also we have explored a mixed formulation for simulating incompressible viscoplastic fluid flows. This formulation does not satisfy the well-known inf-sup condition and therefore a stabilization technique is required to obtain robust results. The Material Point Method is identified as a suitable method to bridge the gap between accuracy and correct tackling of large nonlinear deformation phenomena and for solving free surface problems with viscoplastic fluid flows. Therefore, this work presents two main contributions: on the one hand the development of an irreducible and a mixed formulation for free surface viscoplastic fluids in a nonlinear solid mechanics framework and a weakly and fully incompressible regimes. On the other hand, some stabilization techniques based on the Variational MultiScale (VMS) method are employed for solving the mixed formulation. The proposed formulations, with displacement and pressure as primary variables, are tested through several benchmarks.