The evolution equations for the particle distribution functions are written in a divergence form applicable in three dimensions. From this set, it is shown that the continuity equation and the equations of motion are satisfied in Cartesian, cylindrical and spherical coordinates for all fluids when additional source terms are added to the equations of evolution in the latter two coordinate systems. If the body forces are present, a new set of source functions is required in each coordinate system and these are described as well. Next, the energy equation is derived by using a separate set of particle distribution functions. Modifications of the relevant equations to be applicable to incompressible fluids is described. The incorporation of boundary conditions and the description of the numerical scheme for the simulation of the flows employing the new approach is given. Validation results obtained through the modelling of a mixed convection flow of a Bingham fluid in a lid-driven square cavity, and the steady flow of a Bingham fluid in a pipe of square cross-section are presented. Next, using the cylindrical coordinate version of the evolution equations, numerical modelling of the steady flow of a Bingham fluid and the Herschel-Bulkley fluid in a pipe of circular cross-section have been performed and compared with the simulation results using the augmented Lagrangian method as well as the analytical solutions for the velocity field and the flow rate. Finally, some comments on the theoretical differences between the present approach and the existing formulations regarding Lattice Boltzmann Equations are offered.