ring r = 0,(x,y),dp;
  poly f = 9x16-18x13y2-9x12y3+9x10y4-18x11y2+36x8y4
         +18x7y5-18x5y6+9x6y4-18x3y6-9x2y7+9y8;
  // = 9 * (x5-1y2)^2 * (x6-2x3y2-1x2y3+y4)
  factorize(f);
→ [1]:
→    _[1]=9
→    _[2]=x6-2x3y2-x2y3+y4
→    _[3]=-x5+y2
→ [2]:
→    1,1,2
  // returns factors and multiplicities,
  // first factor is a constant.
  poly g = (y4+x8)*(x2+y2);
  factorize(g);
→ [1]:
→    _[1]=1
→    _[2]=x8+y4
→    _[3]=x2+y2
→ [2]:
→    1,1,1
  // The same in characteristic 2:
  ring s =2,(x,y),dp;
  poly g = (y4+x8)*(x2+y2);
  factorize(g);
→ [1]:
→    _[1]=1
→    _[2]=x+y
→    _[3]=x2+y
→ [2]:
→    1,2,4