The spline element method with constraints is a discretization method where the unknowns are expanded as polynomials on each element and Lagrange multipliers are used to enforce the interelement conditions, the boundary conditions and the constraints in numerical solution of partial differential equations. Spaces of piecewise polynomials with global smoothness conditions are known as multivariate splines and have been extensively studied using the Bernstein-Bezier representation of polynomials. It is used here to write the constraints mentioned above as linear equations. In this paper, we illustrate the robustness of this approach on two singular perturbation problems, a fourth order problem and a Stokes-Darcy flow. It is shown that the method converges uniformly in the perturbation parameter.