diff --git a/src/fp.rs b/src/fp.rs
index 687b82ead..127dd62b2 100644
--- a/src/fp.rs
+++ b/src/fp.rs
@@ -256,6 +256,24 @@ impl Fp {
         CtOption::new(sqrt, sqrt.square().ct_eq(self))
     }
 
+    #[inline]
+    /// Computes the multiplicative inverse of this field
+    /// element, returning None in the case that this element
+    /// is zero.
+    pub fn invert(&self) -> CtOption<Self> {
+        // Exponentiate by p - 2
+        let t = self.pow_vartime(&[
+            0xb9feffffffffaaa9,
+            0x1eabfffeb153ffff,
+            0x6730d2a0f6b0f624,
+            0x64774b84f38512bf,
+            0x4b1ba7b6434bacd7,
+            0x1a0111ea397fe69a,
+        ]);
+
+        CtOption::new(t, !self.is_zero())
+    }
+
     #[inline]
     const fn subtract_p(&self) -> Fp {
         let (r0, borrow) = sbb(self.0[0], MODULUS[0], 0);
@@ -748,3 +766,26 @@ fn test_sqrt() {
         ])
     );
 }
+
+#[test]
+fn test_inversion() {
+    let a = Fp([
+        0x43b43a5078ac2076,
+        0x1ce0763046f8962b,
+        0x724a5276486d735c,
+        0x6f05c2a6282d48fd,
+        0x2095bd5bb4ca9331,
+        0x3b35b3894b0f7da,
+    ]);
+    let b = Fp([
+        0x69ecd7040952148f,
+        0x985ccc2022190f55,
+        0xe19bba36a9ad2f41,
+        0x19bb16c95219dbd8,
+        0x14dcacfdfb478693,
+        0x115ff58afff9a8e1,
+    ]);
+
+    assert_eq!(a.invert().unwrap(), b);
+    assert!(Fp::zero().invert().is_none().unwrap_u8() == 1);
+}
diff --git a/src/fp2.rs b/src/fp2.rs
index b088a2681..8725b6525 100644
--- a/src/fp2.rs
+++ b/src/fp2.rs
@@ -242,6 +242,30 @@ impl Fp2 {
         })
     }
 
+    /// Computes the multiplicative inverse of this field
+    /// element, returning None in the case that this element
+    /// is zero.
+    pub fn invert(&self) -> CtOption<Self> {
+        // We wish to find the multiplicative inverse of a nonzero
+        // element a + bu in Fp2. We leverage an identity
+        //
+        // (a + bu)(a - bu) = a^2 + b^2
+        //
+        // which holds because u^2 = -1. This can be rewritten as
+        //
+        // (a + bu)(a - bu)/(a^2 + b^2) = 1
+        //
+        // because a^2 + b^2 = 0 has no nonzero solutions for (a, b).
+        // This gives that (a - bu)/(a^2 + b^2) is the inverse
+        // of (a + bu). Importantly, this can be computing using
+        // only a single inversion in Fp.
+
+        (self.c0.square() + self.c1.square()).invert().map(|t| Fp2 {
+            c0: self.c0 * t,
+            c1: self.c1 * -t,
+        })
+    }
+
     /// Although this is labeled "vartime", it is only
     /// variable time with respect to the exponent. It
     /// is also not exposed in the public API.
@@ -671,3 +695,48 @@ fn test_sqrt() {
         .is_none()
     ));
 }
+
+#[test]
+fn test_inversion() {
+    let a = Fp2 {
+        c0: Fp::from_raw_unchecked([
+            0x1128ecad67549455,
+            0x9e7a1cff3a4ea1a8,
+            0xeb208d51e08bcf27,
+            0xe98ad40811f5fc2b,
+            0x736c3a59232d511d,
+            0x10acd42d29cfcbb6,
+        ]),
+        c1: Fp::from_raw_unchecked([
+            0xd328e37cc2f58d41,
+            0x948df0858a605869,
+            0x6032f9d56f93a573,
+            0x2be483ef3fffdc87,
+            0x30ef61f88f483c2a,
+            0x1333f55a35725be0,
+        ]),
+    };
+
+    let b = Fp2 {
+        c0: Fp::from_raw_unchecked([
+            0x581a1333d4f48a6,
+            0x58242f6ef0748500,
+            0x292c955349e6da5,
+            0xba37721ddd95fcd0,
+            0x70d167903aa5dfc5,
+            0x11895e118b58a9d5,
+        ]),
+        c1: Fp::from_raw_unchecked([
+            0xeda09d2d7a85d17,
+            0x8808e137a7d1a2cf,
+            0x43ae2625c1ff21db,
+            0xf85ac9fdf7a74c64,
+            0x8fccdda5b8da9738,
+            0x8e84f0cb32cd17d,
+        ]),
+    };
+
+    assert_eq!(a.invert().unwrap(), b);
+
+    assert!(Fp2::zero().invert().is_none().unwrap_u8() == 1);
+}