The material point method (MPM) is a numerical method for simulating ground deformations, that combinesthe benefits of a fixed grid with movingmaterial points, which carry the history dependent information of thecontinuum. During each time step, this information is projected to the grid, on which the equations of motionare solved and a corresponding acceleration field is constructed. Then, the material points undergo motioncorresponding to the acceleration field and the particle properties are updated accordingly. MPM uses a finiteelement method to solve the equations of motion over a triangular grid, but in current versions, the basisfunctions used are only piece-wise linear, which results in low-order spatial convergence for function reconstruction.Furthermore, the piece-wise linear basis functions have discontinuous first derivatives, which causeso called grid-crossing errors when amaterial point crosses between elements, where these discontinuities areoccur. To solve these problems for 2D MPM on a triangular grid, this study investigates the possibility to usehigher-order basis functions with the properties of C1-continuity and non-negativity. B-spline basis functionswould be very suitable for this purpose and earlier studies for 1D have already shown the potential of B-splineMPM. This thesis shows that it is also possible to use B-splines on a 2D triangular grid, which has led to the useof quadratic Powell-Sabin (PS) spline basis functions, which are piece-wise quadratic, C1-continuous, nonnegativebasis functions on triangulations. In this thesis, we have implemented these basis functions inMPMand investigated the performance of PS-spline MPM. Although Constructing and evaluating PS-splines hasproven more cumbersome than the piece-wise linear basis functions classically used in MPM, but PS-splineMPM shows higher-order spatial convergence for function reconstruction and completely erases grid crossingerrors. For these reasons, PS-splineMPM would be very suitable for application in engineering in which largenumbers of basis functions are necessary: due to the higher-order spatial convergence of PS-splineMPM, thenumber of required basis functions could be significantly reduced.Before implementation in engineering application can take place, however, two issues still appear in PSsplineMPM, which should be further researched. The first issue is that PS-spline MPM shows instabilitieswhen elements become nearly empty, as the mass matrix becomes ill-conditioned when this occurs. A way tomitigate this problems is by lumping the mass matrix, however, this leads to the second issue. When the massmatrix is lumped, the accuracy of the solution drops significantly, far more than when classical MPM withpiece-wise linear basis functions is lumped. If these issues are resolved in further research, PS-spline MPMmay prove very suitable for engineering applications.