We study systems of partial differential-algebraic equations (PDAEs) of first order. Classical solutions of the theory of hyperbolic partial differential equation such as discontinuities (shock and contact discontinuities), rarefactions and diffusive traveling waves are extended for variables living on a surface S, which is defined as solution of a set of algebraic equations. We propose here an alternative formulation to study numerically and theoretically the PDAEs by changing the algebraic conditions into partial differential equations with relaxation source terms (PDREs). The solution of such relaxed systems is proved to tend to the surface S, i.e., to satisfy the algebraic equations for long times. We formulate a unified numerical scheme for systems of PDAEs and PDREs. This scheme is naturally parallelizable and has faster convergence. We do not perform a rigorous analysis about the convergence or accuracy for the method, the evidence of its effectiveness is presented by means of simulations for physical and synthetical problems.