This paper focuses on an accurate evaluation of short fibers suspensions models coupling the flow kinematics with the fiber orientation evolution. In coupled models the flow kinematics is usually solved using the finite element method, where the fiber orientation is introduced in the constitutive euqation throught its value in some points (nodes or integration points). In this paper we will compare in a simple steady shear flow, the exact solutions of the extrastresses associated with the fibers' presence with the numerical simulations obtained using both the method of characteristics and the discontinuous Galerkin's method to solve the equation governing the generalized gradient evolution, in order to avoid the introduction of any closure relation. The error introduced if a quadratic closure relation is considered in the constitutive equation will be also quantified.