In this paper we discuss several complete flux schemes for advection-diffusion-reaction problems. We consider both scalar equations as well as systems of equations. For the flux approximations in the latter case, we take into account the coupling between the constituent equations. We study conservation laws with discontinuous diffusion matrix/coefficient and show that the (matrix) harmonic average should be employed in the expressions for the numerical fluxes. The vectorial harmonic complete flux schemes are validated for a test problem