diff --git a/src/constants.rs b/src/constants.rs
index 2f57e0e..a43677d 100644
--- a/src/constants.rs
+++ b/src/constants.rs
@@ -1634,7 +1634,8 @@ mod test {
     fn test_sqrt_constants_sign() {
         let one       = FieldElement([ 1,0,0,0,0,0,0,0,0,0]);
         let minus_one = FieldElement([-1,0,0,0,0,0,0,0,0,0]);
-        let invsqrt_m1 = minus_one.invsqrt().unwrap();
+        let (was_nonzero_square, invsqrt_m1) = minus_one.invsqrt();
+        assert_eq!(was_nonzero_square, 1u8);
         let sign_test_sqrt  = &invsqrt_m1 * &constants::SQRT_M1;
         let sign_test_msqrt = &invsqrt_m1 * &constants::MSQRT_M1;
         // XXX it seems we have flipped the sign relative to
diff --git a/src/decaf.rs b/src/decaf.rs
index 715190e..1b777c9 100644
--- a/src/decaf.rs
+++ b/src/decaf.rs
@@ -71,10 +71,12 @@ impl CompressedDecaf {
         let Z = &FieldElement::one() - &ss; // Z = 1+as^2
         let u = &(&Z * &Z) - &(&constants::d4 * &ss); // u = Z^2 - 4ds^2
         let uss = &u * &ss;
-        let mut v = match uss.invsqrt() {
-            Some(v) => v,
-            None => return None,
-        };
+
+        let (uss_is_nonzero_square, mut v) = uss.invsqrt();
+        if (uss_is_nonzero_square | uss.is_zero()) == 0u8 {
+            return None; // us^2 is nonzero nonsquare
+        }
+
         // Now v = 1/sqrt(us^2) if us^2 is a nonzero square, 0 if us^2 is zero.
         let uv = &v * &u;
         if uv.is_negative_decaf() == 1u8 {
@@ -158,9 +160,9 @@ impl DecafPoint {
         let Z_plus_Y  = &self.0.Z + &Y;
         let Z_minus_Y = &self.0.Z - &Y;
         let t = &constants::a_minus_d * &(&Z_plus_Y * &Z_minus_Y);
+        let (t_is_nonzero_square, mut r) = t.invsqrt();
         // t should always be square (why?)
-        // XXX is it safe to use option types here?
-        let mut r = t.invsqrt().unwrap();
+        debug_assert_eq!( t_is_nonzero_square | t.is_zero(), 1u8 );
 
         // Step 2: Compute u = (a-d)r
         let u = &constants::a_minus_d * &r;
diff --git a/src/field.rs b/src/field.rs
index e0a1203..99bc104 100644
--- a/src/field.rs
+++ b/src/field.rs
@@ -215,6 +215,11 @@ impl FieldElement {
         FieldElement([ 1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ])
     }
 
+    /// Construct -1.
+    pub fn minus_one() -> FieldElement {
+        FieldElement([-1, 0, 0, 0, 0, 0, 0, 0, 0, 0 ])
+    }
+
     fn combine_coeffs(input: &[i64;10]) -> FieldElement { //FeCombine
         let mut c = [0i64;10];
         let mut h = input.clone();
@@ -817,18 +822,20 @@ impl FieldElement {
         t21
     }
 
-    /// Try to compute 1/sqrt(self).
+    /// For `self` a nonzero square, compute 1/sqrt(self) in
+    /// constant time.
+    ///
+    /// It would be much better to use an `Option` type here, but
+    /// doing so forces the caller to branch, which we don't want to
+    /// do.  This seems like the least bad solution.
     ///
     /// # Return
     ///
-    /// * If `self` is zero, returns zero.
-    /// * If `self` is square, returns 1/sqrt(self).
-    /// * If `self` is nonsquare, returns `None`.
-    pub fn invsqrt(&self) -> Option<FieldElement> {
-        // We are to compute v as:
-        //     / 1/sqrt(self)  if self is square, nonzero;
-        // v = |            0  if self is zero;
-        //     \     [reject]  if self is nonsquare.
+    /// - `(1u8, 1/sqrt(self))` if `self` is a nonzero square;
+    /// - `(0u8, zero)`         if `self` is zero;
+    /// - `(0u8, garbage)`      if `self` is nonsquare.
+    ///
+    pub fn invsqrt(&self) -> (u8, FieldElement) {
         //
         // Using the same trick as in ed25519 decoding, we merge the
         // inversion, the square root, and the square test as follows.
@@ -841,22 +848,23 @@ impl FieldElement {
         //    1/β = α^(p-1 - (p+3)/8) = α^((7p-11)/8)
         //                            = α^3 * (α^7)^((p-5)/8).
         //
-        // If α is square, then (1/β)^2 = ±(1/α), so that (1/β)^2 α = ±1.
+        // If α is nonzero square, then (1/β)^2   = ±(1/α),
+        //                      so that (1/β)^2 α = ±1.
         let a3 = &self.square() * self;        // α^3
         let a7 =   &a3.square() * self;        // α^7
         let mut v = &a3 * &a7.pow_p58();       // α^(p-1-(p+3)/8)
-        let check =  self * &v.square();       // ±1 if α is square
+        let check =  self * &v.square();       // ±1 if α is nz square
 
-        if v.is_zero() == 1u8 {
-            return Some(v);  // α was zero all along
-        } else if check == FieldElement::one() {
-            return Some(v);  // computed the correct sqrt
-        } else if check == -&FieldElement::one() {
-                             // wrong sign, multiply by sqrt(-1)
-            return Some(&v * &constants::SQRT_M1);
-        } else {
-            return None;     // input was nonsquare
-        }
+        let correct_sign_sqrt = check.ct_eq(&FieldElement::one());
+        let flipped_sign_sqrt = check.ct_eq(&FieldElement::minus_one());
+
+        // If check = -1, we're off by a factor of sqrt(-1).
+        let v_prime = &constants::SQRT_M1 * &v;
+        v.conditional_assign(&v_prime, flipped_sign_sqrt);
+
+        let was_nonzero_square = correct_sign_sqrt | flipped_sign_sqrt;
+
+        (was_nonzero_square, v)
     }
 
     /// chi calculates `self^((p-1)/2)`.