Summary: ``The correct treatment of boundary conditions is a key step in the development of the SPH method. The SPH community has to face several challenges in this regard---in particular, a primordial aspect for any boundary formulation is to ensure the consistency of the operators in presence of boundaries and free surfaces. A new implementation is proposed, based on the existing numerical boundary integrals formulation. A new kernel expression is developed to compute the Shepard renormalization factor at the boundary purely as a function of the geometry. In order to evaluate this factor, the resulting expression is split into numerical and analytical parts, which allows accurately computing the Shepard factor. The new expression is satisfactorily tested for different planar geometries, showing that problems featuring free surfaces and boundaries are solved. The methodology is also extended to 3-D geometries without great increase in computational cost.''