Abstract: The precise prediction of the gradient of flow variables at the center of an element is essential to the computation of inviscid and viscous flux for a 2nd order accurate flow solver based on unstructured grids. The overall order of accuracy of a flow solver is largely determined by the order of accuracy of gradient computations. The orders of accuracy for several gradient computation methods, including cell-based Green-Gauss method, node-based Green-Gauss method and least-squares method, are investigated theoretically and verified through numerical experiments. During numerical tests, the local coordinate is scaled down in proportion to guarantee that the grid topology is strictly unaltered during grid refinement. The same conclusions are drawn by theoretical derivation and numerical tests:for cell-based Green-Gauss method or node-based Green-Gauss method cooperated with inverse distance weightings, 1st order of accuracy is only achieved on certain grids, but more generally the accuracy is less than 1st order. The 1st order of accuracy for gradient computation is guaranteed independent of the grid topology if the gradient is calculated with least-squares method or node-based Green-Gauss method coupled with interpolation weight which is 2nd order accurate in space. The order of accuracy is further proved by numerical test of subsonic inviscid flow around NACA 0012 airfoil.