Algebraic substitution for polynomial simplification

Question

Consider a polynomial $p(x,y)$ and we want to simplify $p(x,y_0)$ where $y_0$ is a root of some other polynomial $q(y)$.

In Maple I would use something like:

simplify(algsubs(q(y0)=0),p(x,y_0))

But using the following syntax in Mathematica doesn't work:

pxy /. {q(y0)=0}

Here is my precise example:

pxy0=-(x^14*(1 - 2*y0)^6) + 32*(1 - 3*y0)^6*(-1 + y0)^8*
  (1 + 2*y0 - 4*y0^2 + 2*y0^3) + 4*x^12*(1 - 2*y0)^4*
  (4 - 33*y0 + 92*y0^2 - 99*y0^3 + 36*y0^4) + 
 4*x^8*(1 - 6*y0 + 11*y0^2 - 6*y0^3)^2*(37 - 391*y0 + 1699*y0^2 - 3910*y0^3 + 
   4964*y0^4 - 3238*y0^5 + 840*y0^6) - 32*x^2*(-1 + y0)^6*(-1 + 3*y0)^5*
  (5 - 31*y0 + 73*y0^2 - 65*y0^3 - 7*y0^4 + 44*y0^5 - 26*y0^6 + 6*y0^7) + 
 8*x^6*(-1 + y0)^4*(-1 + 3*y0)^3*(-29 + 306*y0 - 1409*y0^2 + 3828*y0^3 - 
   6720*y0^4 + 7432*y0^5 - 4572*y0^6 + 1152*y0^7) + 
 8*x^4*(1 - 4*y0 + 3*y0^2)^4*(31 - 276*y0 + 1020*y0^2 - 2172*y0^3 + 
   3226*y0^4 - 3696*y0^5 + 3032*y0^6 - 1452*y0^7 + 288*y0^8) + 
 2*x^10*(1 - 2*y0)^2*(37 - 586*y0 + 3933*y0^2 - 14594*y0^3 + 32722*y0^4 - 
   45380*y0^5 + 38016*y0^6 - 17604*y0^7 + 3456*y0^8)

But the following expression still contains power of $y_0$ more than 3:

pxy0 /. {y0^3->2 y0^2-y0-1/2} //Expand // Simplify

---

EDIT for the first comment:

This command does not make the job:

(243 - 13*x^2 - 564*y0 - 1588*y0^2 + 7550*y0^3 - 9284*y0^4 + 2040*y0^5 + 5232*y0^6 - 5408*y0^7 + 2240*y0^8 - 384*y0^9)
Simplify[%, {y0^3==2*y0^2-y0-x^2 (1-2 y0)/2-1/2}]

Answer 1

PolynomialReduce is useful for this purpose.

Here is the polynomial and the one we use to reduce it.

pxy0 = -(x^14*(1 - 2*y0)^6) + 
   32*(1 - 3*y0)^6*(-1 + y0)^8*(1 + 2*y0 - 4*y0^2 + 2*y0^3) + 
   4*x^12*(1 - 2*y0)^4*(4 - 33*y0 + 92*y0^2 - 99*y0^3 + 36*y0^4) + 
   4*x^8*(1 - 6*y0 + 11*y0^2 - 6*y0^3)^2*(37 - 391*y0 + 1699*y0^2 - 
      3910*y0^3 + 4964*y0^4 - 3238*y0^5 + 840*y0^6) - 
   32*x^2*(-1 + y0)^6*(-1 + 3*y0)^5*(5 - 31*y0 + 73*y0^2 - 65*y0^3 - 
      7*y0^4 + 44*y0^5 - 26*y0^6 + 6*y0^7) + 
   8*x^6*(-1 + y0)^4*(-1 + 3*y0)^3*(-29 + 306*y0 - 1409*y0^2 + 
      3828*y0^3 - 6720*y0^4 + 7432*y0^5 - 4572*y0^6 + 1152*y0^7) + 
   8*x^4*(1 - 4*y0 + 3*y0^2)^4*(31 - 276*y0 + 1020*y0^2 - 2172*y0^3 + 
      3226*y0^4 - 3696*y0^5 + 3032*y0^6 - 1452*y0^7 + 288*y0^8) + 
   2*x^10*(1 - 2*y0)^2*(37 - 586*y0 + 3933*y0^2 - 14594*y0^3 + 
      32722*y0^4 - 45380*y0^5 + 38016*y0^6 - 17604*y0^7 + 3456*y0^8);
qy = y0^3 - (2*y0^2 - y0 - x^2 (1 - 2 y0)/2 - 1/2);

Now do the reduction.

PolynomialReduce[pxy0, qy /. x -> 0, {y0, x}][[2]]

(* Out[26]= 58338 x^2 + 54108 x^4 + 44509 x^6 + 22460 x^8 + 3092 x^10 + 
 454 x^12 - 17 x^14 - 70728 x^2 y0 - 94698 x^4 y0 - 55724 x^6 y0 - 
 17790 x^8 y0 - 10640 x^10 y0 - 480 x^12 y0 + 4 x^14 y0 + 
 28672 x^2 y0^2 + 41002 x^4 y0^2 + 33606 x^6 y0^2 - 510 x^8 y0^2 + 
 5842 x^10 y0^2 + 504 x^12 y0^2 + 20 x^14 y0^2 *)