diff --git a/scripts/polyrith_sage.py b/scripts/polyrith_sage.py
index 23f5d0d4a09d5..993a41b65f347 100644
--- a/scripts/polyrith_sage.py
+++ b/scripts/polyrith_sage.py
@@ -32,6 +32,26 @@ def create_query(type: str, n_vars: int, eq_list, goal_type):
 '''
     return query
 
+def create_query_radical(type: str, n_vars: int, eq_list, goal_type):
+    """ Create a query to invoke Sage's `MPolynomial_libsingular.lift`. See
+    https://github.com/sagemath/sage/blob/f8df80820dc7321dc9b18c9644c3b8315999670b/src/sage/rings/polynomial/multi_polynomial_libsingular.pyx#L4472-L4518
+    for a description of this method. """
+    var_list = [f"var{i}" for i in range(n_vars)] + ['aux']
+    query = f'''
+import json
+P = {type_str(type)}{var_list}
+[{", ".join(var_list)}] = P.gens()
+p = P({goal_type})
+gens = {eq_list} + [1 - p*aux]
+I = ideal(gens)
+coeffs = P(1).lift(I)
+power = max(cf.degree(aux) for cf in coeffs)
+coeffs = [P(cf.subs(aux = 1/p)*p^power) for cf in coeffs[:int(-1)]]
+js = {{'power': power, 'coeffs': [polynomial_to_string(c) for c in coeffs]}}
+print(json.dumps(js))
+'''
+    return query
+
 class EvaluationError(Exception):
     def __init__(self, ename, evalue, message='Error in Sage communication'):
         self.ename = ename
@@ -68,7 +88,7 @@ def main():
       error_value: Optional[str] }
     ```
     '''
-    command = create_query(sys.argv[2], int(sys.argv[3]), sys.argv[4], sys.argv[5])
+    command = create_query_radical(sys.argv[2], int(sys.argv[3]), sys.argv[4], sys.argv[5])
     final_query = polynomial_formatting_functions + "\n" + command
     if sys.argv[1] == 'tt': # trace dry run enabled
         output = dict(success=True, trace=command)
diff --git a/src/tactic/polyrith.lean b/src/tactic/polyrith.lean
index 17c6ad9872e9c..7f6354e9fa682 100644
--- a/src/tactic/polyrith.lean
+++ b/src/tactic/polyrith.lean
@@ -300,12 +300,18 @@ meta instance : non_null_json_serializable poly :=
     | sum.inr p := pure p
     end}
 
+/-- A schema for the data reported by the Sage calculation -/
+@[derive [non_null_json_serializable, inhabited]]
+structure sage_json_coeffs_and_power :=
+(coeffs : list poly)
+(power  : ℕ)
+
 /-- A schema for success messages from the python script -/
 @[derive [non_null_json_serializable, inhabited]]
 structure sage_json_success :=
 (success : {b : bool // b = tt})
 (trace : option string := none)
-(data : option (list poly) := none)
+(data : option sage_json_coeffs_and_power := none)
 
 /-- A schema for failure messages from the python script -/
 @[derive [non_null_json_serializable, inhabited]]
@@ -316,7 +322,7 @@ structure sage_json_failure :=
 
 /-- Parse the json output from `scripts/polyrith.py` into either an error message, a list of `poly`
 objects, or `none` if only trace output was requested. -/
-meta def convert_sage_output (j : json) : tactic (option (list poly)) :=
+meta def convert_sage_output (j : json) : tactic (option (ℕ × list poly)) :=
 do
   r : sage_json_success ⊕ sage_json_failure ← decorate_ex "internal json error: "
     -- try the error format first, so that if both fail we get the message from the success parser
@@ -326,7 +332,7 @@ do
       fail!"polyrith failed to retrieve a solution from Sage! {f.error_name}: {f.error_value}"
   | sum.inl s := do
       s.trace.mmap trace,
-      pure s.data
+      pure $ s.data.map $ λ sd, (sd.power, sd.coeffs)
   end
 
 /-!
@@ -484,18 +490,19 @@ a call to `linear_combination`.
 -/
 meta def process_output (eq_names : list expr) (m : list expr) (R : expr) (sage_out : json) :
   tactic format := focus1 $ do
-  some coeffs_as_poly ← convert_sage_output sage_out | fail!"internal error: No output available",
+  some (power, coeffs_as_poly) ← convert_sage_output sage_out | fail!"internal error: No output available",
   coeffs_as_pexpr ← coeffs_as_poly.mmap (poly.to_pexpr m),
   let eq_names_pexpr := eq_names.map to_pexpr,
   coeffs_as_expr ← coeffs_as_pexpr.mmap $ λ e, to_expr ``(%%e : %%R),
-  linear_combo.linear_combination eq_names_pexpr coeffs_as_pexpr,
+  linear_combo.linear_combination eq_names_pexpr coeffs_as_pexpr (some power),
   let components := (eq_names.zip coeffs_as_expr).filter
     $ λ pr, bnot $ pr.2.is_app_of `has_zero.zero,
   expr_string ← components_to_lc_format components,
-  let lc_fmt : format := "linear_combination " ++ format.nest 2 (format.group expr_string),
+  let lc_exp : format := if power = 1 then "" else format!" with exponent {power}",
+  let lc_fmt : format := "linear_combination " ++ format.nest 2 (format.group expr_string ++ lc_exp),
   done <|>
     fail!"polyrith found the following certificate, but it failed to close the goal:\n{lc_fmt}",
-  return $ "linear_combination " ++ format.nest 2 (format.group expr_string)
+  return lc_fmt
 
 /-- Tactic for the special case when no hypotheses are available. -/
 meta def no_hypotheses_case : tactic (option format) :=
diff --git a/test/polyrith.lean b/test/polyrith.lean
index c3b71514e1e42..2ba08073c5130 100644
--- a/test/polyrith.lean
+++ b/test/polyrith.lean
@@ -98,10 +98,13 @@ end tactic
 ## SageCell communcation tests
 -/
 
+-- set_option trace.polyrith true
+
 example (x y : ℚ) (h1 : x*y + 2*x = 1) (h2 : x = y) :
   x*y = -2*y + 1 :=
 begin
-  test_sage_output "{\"data\":[\"(poly.const 1/1)\",\"(poly.const -2/1)\"],\"success\":true}",
+  polyrith,
+  -- test_sage_output "{\"data\":[\"(poly.const 1/1)\",\"(poly.const -2/1)\"],\"success\":true}",
   linear_combination h1 - 2 * h2
 end
 
@@ -436,10 +439,39 @@ by test_polyrith
   "linear_combination h1 - h2"
 
 
--- We comment the following tests so that we don't overwhelm the SageCell API.
-
+/- ### Cases with powers -/
 
+example (x y z : ℚ) (h : x = y) (h2 : x * y = 0) : x + y*z = 0 :=
+by test_polyrith
+  "{\"data\":{\"coeffs\":[\"(poly.add (poly.neg (poly.mul (poly.var 1) (poly.pow (poly.var 2) 2))) (poly.var 0))\",\"(poly.add (poly.add (poly.pow (poly.var 2) 2) (poly.mul (poly.const 2/1) (poly.var 2))) (poly.const 1/1))\"],\"power\":2},\"success\":true}"
+  ["ff",
+  "rat",
+  "3",
+  "[(var0 - var1), ((var0 * var1) - 0)]",
+  "((var0 + (var1 * var2)) - 0)"]
+  "linear_combination (-(y * z ^ 2) + x) * h + (z ^ 2 + 2 * z + 1) * h2 with exponent 2"
+
+.
+example (K : Type) (h_this : 3 - 1 = 2)
+  [field K]
+  [char_zero K]
+  {x y z : K}
+  (h₂ : y ^ 3 + x * (3 * z ^ 2) = 0)
+  (h₁ : x ^ 3 + z * (3 * y ^ 2) = 0)
+  (h₀ : y * (3 * x ^ 2) + z ^ 3 = 0)
+  (h : x ^ 3 * y + y ^ 3 * z + z ^ 3 * x = 0) :
+  x = 0 :=
+by test_polyrith
+  "{\"data\":{\"coeffs\":[\"(poly.mul (poly.mul (poly.const 2/7) (poly.var 1)) (poly.pow (poly.var 2) 2))\",\"(poly.sub (poly.pow (poly.var 0) 3) (poly.mul (poly.mul (poly.const 1/7) (poly.pow (poly.var 1) 2)) (poly.var 2)))\",\"(poly.neg (poly.mul (poly.mul (poly.var 0) (poly.var 1)) (poly.var 2)))\",\"(poly.mul (poly.mul (poly.const 1/7) (poly.var 1)) (poly.var 2))\"],\"power\":6},\"success\":true}"
+  ["ff",
+  "K",
+  "3",
+  "[(((var1 ^ 3) + (var0 * (3 * (var2 ^ 2)))) - 0), (((var0 ^ 3) + (var2 * (3 * (var1 ^ 2)))) - 0), (((var1 * (3 * (var0 ^ 2))) + (var2 ^ 3)) - 0), (((((var0 ^ 3) * var1) + ((var1 ^ 3) * var2)) + ((var2 ^ 3) * var0)) - 0)]",
+  "(var0 - 0)"]
+  "linear_combination 2 / 7 * y * z ^ 2 * h₂ + (x ^ 3 - 1 / 7 * y ^ 2 * z) * h₁ - x * y * z * h₀ + 1 / 7 * y *
+    z * h with exponent 6"
 
+-- We comment the following tests so that we don't overwhelm the SageCell API.
 
 
 /-