It is well known that Navier-Stokes equations are not valid for those high-Knudsen and high-Mach flows, in which the local thermodynamically non-equilibrium effects are dominant. To extend the non-equilibrium describing the ability of macroscopic equations, Nonlinear Coupled Constitutive Relation (NCCR) model was developed from Eu's generalized hydrodynamic equa-tions to substitute linear Newton's law of viscosity and Fourier's law of heat conduction in conser-vation laws. In the NCCR model, how to solve the decomposed constitutive equations with reasonable computational cost is a key ingredient of this scheme. In this paper, an analytic method is proposed firstly. Compared to the iterative procedure in the conventional NCCR model, the ana-lytic method not only obtains exact roots of the decomposed constitutive polynomials, but also pre-serves the nonlinear constitutive relations in the original framework of NCCR methods. Numerical tests to assess the efficiency and accuracy of the proposed method are conducted for argon shock structures, Couette flows, two-dimensional hypersonic flows over a cylinder and three-dimensional supersonic flows over a three-dimensional sphere. These superior advantages of the current method are expected to render itself a powerful tool for simulating the hypersonic rarefied flows and microscale flows of high Knudsen number for engineering applications.