This paper presents a novel approach to model pressure and flow regulating devices in the context of the Newton-Raphson loop method for water distribution network simulation. The proposed approach uses a symmetric matrix for the underlying linear systems, which enables simpler implementation and faster solution, while producing iterations very close to the global gradient algorithm of Epanet. The structure of the matrix is kept unchanged regardless of the operational status of the valves. The paper presents results that validate its formulation, accuracy and speed in various case studies.