In 1999, Christopher gave a necessary and sufficient condition for polynomial Liénard centers, which requires coupled functional equations, where the primitive functions of the damping function and the restoring function are involved, to have polynomial solutions. In order to judge whether the coupled functional equations are solvable, in this paper we give an algorithm to compute a Gröbner basis for irreducible decomposition of algebraic varieties so as to find algebraic relations among coefficients of the damping function and the restoring function. We demonstrate the algorithm for polynomial Liénard systems of degree 5, which are divided into 25 cases. We find all conditions of those coefficients for the polynomial Liénard center in 13 cases and prove that the origin is not a center in the other 12 cases.