From 83d713f4cad84e22fd9b3b23c83765bcd2913bee Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 12:36:30 -0400
Subject: [PATCH 01/20] Vendor lizard hash to curve from Signal

---
 curve25519-dalek/src/lizard/LICENSE           |  21 ++
 curve25519-dalek/src/lizard/jacobi_quartic.rs |  72 ++++
 .../src/lizard/lizard_constants.rs            |  49 +++
 .../src/lizard/lizard_ristretto.rs            | 332 ++++++++++++++++++
 curve25519-dalek/src/lizard/mod.rs            |  13 +
 curve25519-dalek/src/lizard/u32_constants.rs  |  46 +++
 curve25519-dalek/src/lizard/u64_constants.rs  |  63 ++++
 7 files changed, 596 insertions(+)
 create mode 100644 curve25519-dalek/src/lizard/LICENSE
 create mode 100644 curve25519-dalek/src/lizard/jacobi_quartic.rs
 create mode 100644 curve25519-dalek/src/lizard/lizard_constants.rs
 create mode 100644 curve25519-dalek/src/lizard/lizard_ristretto.rs
 create mode 100644 curve25519-dalek/src/lizard/mod.rs
 create mode 100644 curve25519-dalek/src/lizard/u32_constants.rs
 create mode 100644 curve25519-dalek/src/lizard/u64_constants.rs

diff --git a/curve25519-dalek/src/lizard/LICENSE b/curve25519-dalek/src/lizard/LICENSE
new file mode 100644
index 000000000..2df7c9688
--- /dev/null
+++ b/curve25519-dalek/src/lizard/LICENSE
@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) 2019 Bas Westerbaan
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
diff --git a/curve25519-dalek/src/lizard/jacobi_quartic.rs b/curve25519-dalek/src/lizard/jacobi_quartic.rs
new file mode 100644
index 000000000..af64b4580
--- /dev/null
+++ b/curve25519-dalek/src/lizard/jacobi_quartic.rs
@@ -0,0 +1,72 @@
+//! Helper functions for use with Lizard
+
+#![allow(non_snake_case)]
+
+use subtle::Choice;
+use subtle::ConditionallyNegatable;
+use subtle::ConditionallySelectable;
+use subtle::ConstantTimeEq;
+
+use super::lizard_constants;
+use crate::constants;
+
+use crate::field::FieldElement;
+
+/// Represents a point (s,t) on the the Jacobi quartic associated
+/// to the Edwards curve.
+#[derive(Copy, Clone)]
+#[allow(missing_docs)]
+pub struct JacobiPoint {
+    pub S: FieldElement,
+    pub T: FieldElement,
+}
+
+impl JacobiPoint {
+    /// Elligator2 is defined in two steps: first a field element is converted
+    /// to a point (s,t) on the Jacobi quartic associated to the Edwards curve.
+    /// Then this point is mapped to a point on the Edwards curve.
+    /// This function computes a field element that is mapped to a given (s,t)
+    /// with Elligator2 if it exists.
+    pub(crate) fn elligator_inv(&self) -> (Choice, FieldElement) {
+        let mut out = FieldElement::ZERO;
+
+        // Special case: s = 0.  If s is zero, either t = 1 or t = -1.
+        // If t=1, then sqrt(i*d) is the preimage.  Otherwise it's 0.
+        let s_is_zero = self.S.is_zero();
+        let t_equals_one = self.T.ct_eq(&FieldElement::ONE);
+        out.conditional_assign(&lizard_constants::SQRT_ID, t_equals_one);
+        let mut ret = s_is_zero;
+        let mut done = s_is_zero;
+
+        // a := (t+1) (d+1)/(d-1)
+        let a = &(&self.T + &FieldElement::ONE) * &lizard_constants::DP1_OVER_DM1;
+        let a2 = a.square();
+
+        // y := 1/sqrt(i (s^4 - a^2)).
+        let s2 = self.S.square();
+        let s4 = s2.square();
+        let invSqY = &(&s4 - &a2) * &constants::SQRT_M1;
+
+        // There is no preimage if the square root of i*(s^4-a^2) does not exist.
+        let (sq, y) = invSqY.invsqrt();
+        ret |= sq;
+        done |= !sq;
+
+        // x := (a + sign(s)*s^2) y
+        let mut pms2 = s2;
+        pms2.conditional_negate(self.S.is_negative());
+        let mut x = &(&a + &pms2) * &y;
+        let x_is_negative = x.is_negative();
+        x.conditional_negate(x_is_negative);
+        out.conditional_assign(&x, !done);
+
+        (ret, out)
+    }
+
+    pub(crate) fn dual(&self) -> JacobiPoint {
+        JacobiPoint {
+            S: -(&self.S),
+            T: -(&self.T),
+        }
+    }
+}
diff --git a/curve25519-dalek/src/lizard/lizard_constants.rs b/curve25519-dalek/src/lizard/lizard_constants.rs
new file mode 100644
index 000000000..bed03d633
--- /dev/null
+++ b/curve25519-dalek/src/lizard/lizard_constants.rs
@@ -0,0 +1,49 @@
+//! Constants for use in Lizard
+//!
+//! Could be moved into backend/serial/u??/constants.rs
+
+#[cfg(curve25519_dalek_bits = "64")]
+pub(crate) use super::u64_constants::*;
+
+#[cfg(curve25519_dalek_bits = "32")]
+pub(crate) use super::u32_constants::*;
+
+// ------------------------------------------------------------------------
+// Tests
+// ------------------------------------------------------------------------
+
+#[cfg(test)]
+mod test {
+
+    use super::*;
+    use crate::constants;
+    use crate::field::FieldElement;
+
+    #[test]
+    fn test_lizard_constants() {
+        let (_, sqrt_id) = FieldElement::sqrt_ratio_i(
+            &(&constants::SQRT_M1 * &constants::EDWARDS_D),
+            &FieldElement::ONE,
+        );
+        assert_eq!(sqrt_id, SQRT_ID);
+
+        assert_eq!(
+            &(&constants::EDWARDS_D + &FieldElement::ONE)
+                * &(&constants::EDWARDS_D - &FieldElement::ONE).invert(),
+            DP1_OVER_DM1
+        );
+
+        assert_eq!(
+            MDOUBLE_INVSQRT_A_MINUS_D,
+            -&(&constants::INVSQRT_A_MINUS_D + &constants::INVSQRT_A_MINUS_D)
+        );
+
+        assert_eq!(
+            MIDOUBLE_INVSQRT_A_MINUS_D,
+            &MDOUBLE_INVSQRT_A_MINUS_D * &constants::SQRT_M1
+        );
+
+        let (_, invsqrt_one_plus_d) = (&constants::EDWARDS_D + &FieldElement::ONE).invsqrt();
+        assert_eq!(-&invsqrt_one_plus_d, MINVSQRT_ONE_PLUS_D);
+    }
+}
diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
new file mode 100644
index 000000000..13b2e0f16
--- /dev/null
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -0,0 +1,332 @@
+//! Defines additional methods on RistrettoPoint for Lizard
+
+#![allow(non_snake_case)]
+
+use digest::generic_array::typenum::U32;
+use digest::Digest;
+
+use crate::constants;
+use crate::field::FieldElement;
+
+use subtle::Choice;
+use subtle::ConditionallySelectable;
+use subtle::ConstantTimeEq;
+
+use crate::edwards::EdwardsPoint;
+
+use super::jacobi_quartic::JacobiPoint;
+use super::lizard_constants;
+
+use crate::ristretto::RistrettoPoint;
+#[allow(unused_imports)]
+use core::prelude::*;
+
+impl RistrettoPoint {
+    /// Directly encode 253 bits as a RistrettoPoint, using Elligator
+    pub fn from_uniform_bytes_single_elligator(bytes: &[u8; 32]) -> RistrettoPoint {
+        RistrettoPoint::elligator_ristretto_flavor(&FieldElement::from_bytes(bytes))
+    }
+
+    /// Encode 16 bytes of data to a RistrettoPoint, using the Lizard method
+    pub fn lizard_encode<D: Digest>(data: &[u8; 16]) -> RistrettoPoint
+    where
+        D: Digest<OutputSize = U32>,
+    {
+        let mut fe_bytes: [u8; 32] = Default::default();
+
+        let digest = D::digest(data);
+        fe_bytes[0..32].copy_from_slice(digest.as_slice());
+        fe_bytes[8..24].copy_from_slice(data);
+        fe_bytes[0] &= 254; // make positive since Elligator on r and -r is the same
+        fe_bytes[31] &= 63;
+        let fe = FieldElement::from_bytes(&fe_bytes);
+        RistrettoPoint::elligator_ristretto_flavor(&fe)
+    }
+
+    /// Decode 16 bytes of data from a RistrettoPoint, using the Lizard method
+    pub fn lizard_decode<D: Digest>(&self) -> Option<[u8; 16]>
+    where
+        D: Digest<OutputSize = U32>,
+    {
+        let mut result: [u8; 16] = Default::default();
+        let mut h: [u8; 32] = Default::default();
+        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
+        let mut n_found = 0;
+        for (j, fe_j) in fes.iter().enumerate() {
+            let mut ok = Choice::from((mask >> j) & 1);
+            let buf2 = fe_j.as_bytes(); // array
+            h.copy_from_slice(&D::digest(&buf2[8..24])); // array
+            h[8..24].copy_from_slice(&buf2[8..24]);
+            h[0] &= 254;
+            h[31] &= 63;
+            ok &= h.ct_eq(&buf2);
+            for i in 0..16 {
+                result[i] = u8::conditional_select(&result[i], &buf2[8 + i], ok);
+            }
+            n_found += ok.unwrap_u8();
+        }
+        if n_found == 1 {
+            Some(result)
+        } else {
+            None
+        }
+    }
+
+    /// Directly encode 253 bits as a RistrettoPoint, using Elligator
+    pub fn encode_253_bits(data: &[u8; 32]) -> Option<RistrettoPoint> {
+        if data.len() != 32 {
+            return None;
+        }
+
+        let fe = FieldElement::from_bytes(data);
+        let p = RistrettoPoint::elligator_ristretto_flavor(&fe);
+        Some(p)
+    }
+
+    /// Directly decode a RistrettoPoint as 253 bits, using Elligator
+    pub fn decode_253_bits(&self) -> (u8, [[u8; 32]; 8]) {
+        let mut ret = [[0u8; 32]; 8];
+        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
+
+        for j in 0..8 {
+            ret[j] = fes[j].as_bytes();
+        }
+        (mask, ret)
+    }
+
+    /// Return the coset self + E[4], for debugging.
+    pub fn xcoset4(&self) -> [EdwardsPoint; 4] {
+        [
+            self.0,
+            self.0 + constants::EIGHT_TORSION[2],
+            self.0 + constants::EIGHT_TORSION[4],
+            self.0 + constants::EIGHT_TORSION[6],
+        ]
+    }
+
+    /// Computes the at most 8 positive FieldElements f such that
+    /// self == elligator_ristretto_flavor(f).
+    /// Assumes self is even.
+    ///
+    /// Returns a bitmask of which elements in fes are set.
+    pub fn elligator_ristretto_flavor_inverse(&self) -> (u8, [FieldElement; 8]) {
+        // Elligator2 computes a Point from a FieldElement in two steps: first
+        // it computes a (s,t) on the Jacobi quartic and then computes the
+        // corresponding even point on the Edwards curve.
+        //
+        // We invert in three steps.  Any Ristretto point has four representatives
+        // as even Edwards points.  For each of those even Edwards points,
+        // there are two points on the Jacobi quartic that map to it.
+        // Each of those eight points on the Jacobi quartic might have an
+        // Elligator2 preimage.
+        //
+        // Essentially we first loop over the four representatives of our point,
+        // then for each of them consider both points on the Jacobi quartic and
+        // check whether they have an inverse under Elligator2.  We take the
+        // following shortcut though.
+        //
+        // We can compute two Jacobi quartic points for (x,y) and (-x,-y)
+        // at the same time.  The four Jacobi quartic points are two of
+        // such pairs.
+
+        let mut mask: u8 = 0;
+        let jcs = self.to_jacobi_quartic_ristretto();
+        let mut ret = [FieldElement::ONE; 8];
+
+        for i in 0..4 {
+            let (ok, fe) = jcs[i].elligator_inv();
+            let mut tmp: u8 = 0;
+            ret[2 * i] = fe;
+            tmp.conditional_assign(&1, ok);
+            mask |= tmp << (2 * i);
+
+            let jc = jcs[i].dual();
+            let (ok, fe) = jc.elligator_inv();
+            let mut tmp: u8 = 0;
+            ret[2 * i + 1] = fe;
+            tmp.conditional_assign(&1, ok);
+            mask |= tmp << (2 * i + 1);
+        }
+
+        (mask, ret)
+    }
+
+    /// Find a point on the Jacobi quartic associated to each of the four
+    /// points Ristretto equivalent to p.
+    ///
+    /// There is one exception: for (0,-1) there is no point on the quartic and
+    /// so we repeat one on the quartic equivalent to (0,1).
+    fn to_jacobi_quartic_ristretto(self) -> [JacobiPoint; 4] {
+        let x2 = self.0.X.square(); // X^2
+        let y2 = self.0.Y.square(); // Y^2
+        let y4 = y2.square(); // Y^4
+        let z2 = self.0.Z.square(); // Z^2
+        let z_min_y = &self.0.Z - &self.0.Y; // Z - Y
+        let z_pl_y = &self.0.Z + &self.0.Y; // Z + Y
+        let z2_min_y2 = &z2 - &y2; // Z^2 - Y^2
+
+        // gamma := 1/sqrt( Y^4 X^2 (Z^2 - Y^2) )
+        let (_, gamma) = (&(&y4 * &x2) * &z2_min_y2).invsqrt();
+
+        let den = &gamma * &y2;
+
+        let s_over_x = &den * &z_min_y;
+        let sp_over_xp = &den * &z_pl_y;
+
+        let s0 = &s_over_x * &self.0.X;
+        let s1 = &(-(&sp_over_xp)) * &self.0.X;
+
+        // t_0 := -2/sqrt(-d-1) * Z * sOverX
+        // t_1 := -2/sqrt(-d-1) * Z * spOverXp
+        let tmp = &lizard_constants::MDOUBLE_INVSQRT_A_MINUS_D * &self.0.Z;
+        let mut t0 = &tmp * &s_over_x;
+        let mut t1 = &tmp * &sp_over_xp;
+
+        // den := -1/sqrt(1+d) (Y^2 - Z^2) gamma
+        let den = &(&(-(&z2_min_y2)) * &lizard_constants::MINVSQRT_ONE_PLUS_D) * &gamma;
+
+        // Same as before but with the substitution (X, Y, Z) = (Y, X, i*Z)
+        let iz = &constants::SQRT_M1 * &self.0.Z; // iZ
+        let iz_min_x = &iz - &self.0.X; // iZ - X
+        let iz_pl_x = &iz + &self.0.X; // iZ + X
+
+        let s_over_y = &den * &iz_min_x;
+        let sp_over_yp = &den * &iz_pl_x;
+
+        let mut s2 = &s_over_y * &self.0.Y;
+        let mut s3 = &(-(&sp_over_yp)) * &self.0.Y;
+
+        // t_2 := -2/sqrt(-d-1) * i*Z * sOverY
+        // t_3 := -2/sqrt(-d-1) * i*Z * spOverYp
+        let tmp = &lizard_constants::MDOUBLE_INVSQRT_A_MINUS_D * &iz;
+        let mut t2 = &tmp * &s_over_y;
+        let mut t3 = &tmp * &sp_over_yp;
+
+        // Special case: X=0 or Y=0.  Then return
+        //
+        //  (0,1)   (1,-2i/sqrt(-d-1)   (-1,-2i/sqrt(-d-1))
+        //
+        // Note that if X=0 or Y=0, then s_i = t_i = 0.
+        let x_or_y_is_zero = self.0.X.is_zero() | self.0.Y.is_zero();
+        t0.conditional_assign(&FieldElement::ONE, x_or_y_is_zero);
+        t1.conditional_assign(&FieldElement::ONE, x_or_y_is_zero);
+        t2.conditional_assign(
+            &lizard_constants::MIDOUBLE_INVSQRT_A_MINUS_D,
+            x_or_y_is_zero,
+        );
+        t3.conditional_assign(
+            &lizard_constants::MIDOUBLE_INVSQRT_A_MINUS_D,
+            x_or_y_is_zero,
+        );
+        s2.conditional_assign(&FieldElement::ONE, x_or_y_is_zero);
+        s3.conditional_assign(&(-(&FieldElement::ONE)), x_or_y_is_zero);
+
+        [
+            JacobiPoint { S: s0, T: t0 },
+            JacobiPoint { S: s1, T: t1 },
+            JacobiPoint { S: s2, T: t2 },
+            JacobiPoint { S: s3, T: t3 },
+        ]
+    }
+}
+
+// ------------------------------------------------------------------------
+// Tests
+// ------------------------------------------------------------------------
+
+#[cfg(test)]
+mod test {
+
+    use sha2;
+
+    use self::sha2::Sha256;
+    use super::*;
+    use crate::ristretto::CompressedRistretto;
+    use rand_core::RngCore;
+    #[cfg(feature = "rand")]
+    use rand_os::OsRng;
+
+    fn test_lizard_encode_helper(data: &[u8; 16], result: &[u8; 32]) {
+        let p = RistrettoPoint::lizard_encode::<Sha256>(data);
+        let p_bytes = p.compress().to_bytes();
+        assert!(&p_bytes == result);
+        let p = CompressedRistretto::from_slice(&p_bytes)
+            .unwrap()
+            .decompress()
+            .unwrap();
+        let data_out = p.lizard_decode::<Sha256>().unwrap();
+        assert!(&data_out == data);
+    }
+
+    #[test]
+    fn test_lizard_encode() {
+        test_lizard_encode_helper(
+            &[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+            &[
+                0xf0, 0xb7, 0xe3, 0x44, 0x84, 0xf7, 0x4c, 0xf0, 0xf, 0x15, 0x2, 0x4b, 0x73, 0x85,
+                0x39, 0x73, 0x86, 0x46, 0xbb, 0xbe, 0x1e, 0x9b, 0xc7, 0x50, 0x9a, 0x67, 0x68, 0x15,
+                0x22, 0x7e, 0x77, 0x4f,
+            ],
+        );
+
+        test_lizard_encode_helper(
+            &[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1],
+            &[
+                0xcc, 0x92, 0xe8, 0x1f, 0x58, 0x5a, 0xfc, 0x5c, 0xaa, 0xc8, 0x86, 0x60, 0xd8, 0xd1,
+                0x7e, 0x90, 0x25, 0xa4, 0x44, 0x89, 0xa3, 0x63, 0x4, 0x21, 0x23, 0xf6, 0xaf, 0x7,
+                0x2, 0x15, 0x6e, 0x65,
+            ],
+        );
+
+        test_lizard_encode_helper(
+            &[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15],
+            &[
+                0xc8, 0x30, 0x57, 0x3f, 0x8a, 0x8e, 0x77, 0x78, 0x67, 0x1f, 0x76, 0xcd, 0xc7, 0x96,
+                0xdc, 0xa, 0x23, 0x5c, 0xf1, 0x77, 0xf1, 0x97, 0xd9, 0xfc, 0xba, 0x6, 0xe8, 0x4e,
+                0x96, 0x24, 0x74, 0x44,
+            ],
+        );
+    }
+
+    #[test]
+    fn test_elligator_inv() {
+        let mut rng = rand::thread_rng();
+
+        for i in 0..100 {
+            let mut fe_bytes = [0u8; 32];
+
+            if i == 0 {
+                // Test for first corner-case: fe = 0
+                fe_bytes = [0u8; 32];
+            } else if i == 1 {
+                // Test for second corner-case: fe = +sqrt(i*d)
+                fe_bytes = [
+                    168, 27, 92, 74, 203, 42, 48, 117, 170, 109, 234, 14, 45, 169, 188, 205, 21,
+                    110, 235, 115, 153, 84, 52, 117, 151, 235, 123, 244, 88, 85, 179, 5,
+                ];
+            } else {
+                // For the rest, just generate a random field element to test.
+                rng.fill_bytes(&mut fe_bytes);
+            }
+            fe_bytes[0] &= 254; // positive
+            fe_bytes[31] &= 127; // < 2^255-19
+            let fe = FieldElement::from_bytes(&fe_bytes);
+
+            let pt = RistrettoPoint::elligator_ristretto_flavor(&fe);
+            for pt2 in &pt.xcoset4() {
+                let (mask, fes) = RistrettoPoint(*pt2).elligator_ristretto_flavor_inverse();
+
+                let mut found = false;
+                for (j, fe_j) in fes.iter().enumerate() {
+                    if mask & (1 << j) != 0 {
+                        assert_eq!(RistrettoPoint::elligator_ristretto_flavor(fe_j), pt);
+                        if *fe_j == fe {
+                            found = true;
+                        }
+                    }
+                }
+                assert!(found);
+            }
+        }
+    }
+}
diff --git a/curve25519-dalek/src/lizard/mod.rs b/curve25519-dalek/src/lizard/mod.rs
new file mode 100644
index 000000000..df1ac4506
--- /dev/null
+++ b/curve25519-dalek/src/lizard/mod.rs
@@ -0,0 +1,13 @@
+//! The Lizard method for encoding/decoding 16 bytes into Ristretto points.
+
+#![allow(non_snake_case)]
+
+#[cfg(curve25519_dalek_bits = "32")]
+mod u32_constants;
+
+#[cfg(curve25519_dalek_bits = "64")]
+mod u64_constants;
+
+pub mod jacobi_quartic;
+pub mod lizard_constants;
+pub mod lizard_ristretto;
diff --git a/curve25519-dalek/src/lizard/u32_constants.rs b/curve25519-dalek/src/lizard/u32_constants.rs
new file mode 100644
index 000000000..22720feba
--- /dev/null
+++ b/curve25519-dalek/src/lizard/u32_constants.rs
@@ -0,0 +1,46 @@
+use cfg_if::cfg_if;
+
+cfg_if! {
+    if #[cfg(curve25519_dalek_backend = "fiat")] {
+        pub use crate::backend::serial::fiat_u32::field::FieldElement2625;
+
+        const fn field_element(element: [u32; 10]) -> FieldElement2625 {
+            FieldElement2625(fiat_crypto::curve25519_32::fiat_25519_tight_field_element(element))
+        }
+    } else {
+        pub use crate::backend::serial::u32::field::FieldElement2625;
+
+        const fn field_element(element: [u32; 10]) -> FieldElement2625 {
+            FieldElement2625(element)
+        }
+    }
+}
+
+/// `= sqrt(i*d)`, where `i = +sqrt(-1)` and `d` is the Edwards curve parameter.
+pub const SQRT_ID: FieldElement2625 = field_element([
+    39590824, 701138, 28659366, 23623507, 53932708, 32206357, 36326585, 24309414, 26167230, 1494357,
+]);
+
+/// `= (d+1)/(d-1)`, where `d` is the Edwards curve parameter.
+pub const DP1_OVER_DM1: FieldElement2625 = field_element([
+    58833708, 32184294, 62457071, 26110240, 19032991, 27203620, 7122892, 18068959, 51019405,
+    3776288,
+]);
+
+/// `= -2/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
+pub const MDOUBLE_INVSQRT_A_MINUS_D: FieldElement2625 = field_element([
+    54885894, 25242303, 55597453, 9067496, 51808079, 33312638, 25456129, 14121551, 54921728,
+    3972023,
+]);
+
+/// `= -2i/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters
+/// and `i = +sqrt(-1)`.
+pub const MIDOUBLE_INVSQRT_A_MINUS_D: FieldElement2625 = field_element([
+    58178520, 23970840, 26444491, 29801899, 41064376, 743696, 2900628, 27920316, 41968995, 5270573,
+]);
+
+/// `= -1/sqrt(1+d)`, where `d` is the Edwards curve parameters.
+pub const MINVSQRT_ONE_PLUS_D: FieldElement2625 = field_element([
+    38019585, 4791795, 20332186, 18653482, 46576675, 33182583, 65658549, 2817057, 12569934,
+    30919145,
+]);
diff --git a/curve25519-dalek/src/lizard/u64_constants.rs b/curve25519-dalek/src/lizard/u64_constants.rs
new file mode 100644
index 000000000..e034a3890
--- /dev/null
+++ b/curve25519-dalek/src/lizard/u64_constants.rs
@@ -0,0 +1,63 @@
+use cfg_if::cfg_if;
+
+cfg_if! {
+    if #[cfg(curve25519_dalek_backend = "fiat")] {
+        pub use crate::backend::serial::fiat_u64::field::FieldElement51;
+
+        const fn field_element(element: [u64; 5]) -> FieldElement51 {
+            FieldElement51(fiat_crypto::curve25519_64::fiat_25519_tight_field_element(element))
+        }
+    } else {
+        pub use crate::backend::serial::u64::field::FieldElement51;
+
+        const fn field_element(element: [u64; 5]) -> FieldElement51 {
+            FieldElement51(element)
+        }
+    }
+}
+
+/// `= sqrt(i*d)`, where `i = +sqrt(-1)` and `d` is the Edwards curve parameter.
+pub const SQRT_ID: FieldElement51 = field_element([
+    2298852427963285,
+    3837146560810661,
+    4413131899466403,
+    3883177008057528,
+    2352084440532925,
+]);
+
+/// `= (d+1)/(d-1)`, where `d` is the Edwards curve parameter.
+pub const DP1_OVER_DM1: FieldElement51 = field_element([
+    2159851467815724,
+    1752228607624431,
+    1825604053920671,
+    1212587319275468,
+    253422448836237,
+]);
+
+/// `= -2/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters.
+pub const MDOUBLE_INVSQRT_A_MINUS_D: FieldElement51 = field_element([
+    1693982333959686,
+    608509411481997,
+    2235573344831311,
+    947681270984193,
+    266558006233600,
+]);
+
+/// `= -2i/sqrt(a-d)`, where `a = -1 (mod p)`, `d` are the Edwards curve parameters
+/// and `i = +sqrt(-1)`.
+pub const MIDOUBLE_INVSQRT_A_MINUS_D: FieldElement51 = field_element([
+    1608655899704280,
+    1999971613377227,
+    49908634785720,
+    1873700692181652,
+    353702208628067,
+]);
+
+/// `= -1/sqrt(1+d)`, where `d` is the Edwards curve parameters.
+pub const MINVSQRT_ONE_PLUS_D: FieldElement51 = field_element([
+    321571956990465,
+    1251814006996634,
+    2226845496292387,
+    189049560751797,
+    2074948709371214,
+]);

From 24625f0f5597a74a69373556404f06c481e08025 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 12:39:39 -0400
Subject: [PATCH 02/20] Fmt

---
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 8 ++------
 1 file changed, 2 insertions(+), 6 deletions(-)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 13b2e0f16..2e5e86d07 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -2,8 +2,8 @@
 
 #![allow(non_snake_case)]
 
-use digest::generic_array::typenum::U32;
 use digest::Digest;
+use digest::generic_array::typenum::U32;
 
 use crate::constants;
 use crate::field::FieldElement;
@@ -65,11 +65,7 @@ impl RistrettoPoint {
             }
             n_found += ok.unwrap_u8();
         }
-        if n_found == 1 {
-            Some(result)
-        } else {
-            None
-        }
+        if n_found == 1 { Some(result) } else { None }
     }
 
     /// Directly encode 253 bits as a RistrettoPoint, using Elligator

From 1caee9bb2e2173624ad59f43be486432b424d6eb Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 12:53:34 -0400
Subject: [PATCH 03/20] Make lizard compile; tests pass

---
 curve25519-dalek/Cargo.toml                     |  1 +
 curve25519-dalek/src/lib.rs                     |  3 +++
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 10 ++++------
 3 files changed, 8 insertions(+), 6 deletions(-)

diff --git a/curve25519-dalek/Cargo.toml b/curve25519-dalek/Cargo.toml
index 1fbd72d00..d78d0ce43 100644
--- a/curve25519-dalek/Cargo.toml
+++ b/curve25519-dalek/Cargo.toml
@@ -82,6 +82,7 @@ legacy_compatibility = []
 group = ["dep:group", "rand_core"]
 group-bits = ["group", "ff/bits"]
 digest = ["dep:digest"]
+lizard = []
 
 [target.'cfg(all(not(curve25519_dalek_backend = "fiat"), not(curve25519_dalek_backend = "serial"), target_arch = "x86_64"))'.dependencies]
 curve25519-dalek-derive = "0.1"
diff --git a/curve25519-dalek/src/lib.rs b/curve25519-dalek/src/lib.rs
index 24e0fa5b8..623918168 100644
--- a/curve25519-dalek/src/lib.rs
+++ b/curve25519-dalek/src/lib.rs
@@ -90,6 +90,9 @@ pub(crate) mod backend;
 // Generic code for window lookups
 pub(crate) mod window;
 
+#[cfg(feature = "lizard")]
+pub mod lizard;
+
 pub use crate::{
     edwards::EdwardsPoint, montgomery::MontgomeryPoint, ristretto::RistrettoPoint, scalar::Scalar,
 };
diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 2e5e86d07..a212a159f 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -3,7 +3,7 @@
 #![allow(non_snake_case)]
 
 use digest::Digest;
-use digest::generic_array::typenum::U32;
+use digest::consts::U32;
 
 use crate::constants;
 use crate::field::FieldElement;
@@ -54,7 +54,7 @@ impl RistrettoPoint {
         let mut n_found = 0;
         for (j, fe_j) in fes.iter().enumerate() {
             let mut ok = Choice::from((mask >> j) & 1);
-            let buf2 = fe_j.as_bytes(); // array
+            let buf2 = fe_j.to_bytes(); // array
             h.copy_from_slice(&D::digest(&buf2[8..24])); // array
             h[8..24].copy_from_slice(&buf2[8..24]);
             h[0] &= 254;
@@ -85,7 +85,7 @@ impl RistrettoPoint {
         let (mask, fes) = self.elligator_ristretto_flavor_inverse();
 
         for j in 0..8 {
-            ret[j] = fes[j].as_bytes();
+            ret[j] = fes[j].to_bytes();
         }
         (mask, ret)
     }
@@ -239,8 +239,6 @@ mod test {
     use super::*;
     use crate::ristretto::CompressedRistretto;
     use rand_core::RngCore;
-    #[cfg(feature = "rand")]
-    use rand_os::OsRng;
 
     fn test_lizard_encode_helper(data: &[u8; 16], result: &[u8; 32]) {
         let p = RistrettoPoint::lizard_encode::<Sha256>(data);
@@ -286,7 +284,7 @@ mod test {
 
     #[test]
     fn test_elligator_inv() {
-        let mut rng = rand::thread_rng();
+        let mut rng = rand::rng();
 
         for i in 0..100 {
             let mut fe_bytes = [0u8; 32];

From 1b19e868618ad246c102d6a2e30014b3937e536e Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 12:54:01 -0400
Subject: [PATCH 04/20] Move lizard license to acknowledgements

---
 curve25519-dalek/ACKNOWLEDGEMENTS.md | 63 ++++++++++++++++++++++++++++
 curve25519-dalek/LICENSE             | 36 ----------------
 curve25519-dalek/src/lizard/LICENSE  | 21 ----------
 3 files changed, 63 insertions(+), 57 deletions(-)
 create mode 100644 curve25519-dalek/ACKNOWLEDGEMENTS.md
 delete mode 100644 curve25519-dalek/src/lizard/LICENSE

diff --git a/curve25519-dalek/ACKNOWLEDGEMENTS.md b/curve25519-dalek/ACKNOWLEDGEMENTS.md
new file mode 100644
index 000000000..49e4fde96
--- /dev/null
+++ b/curve25519-dalek/ACKNOWLEDGEMENTS.md
@@ -0,0 +1,63 @@
+# Go ed25519
+
+Portions of curve25519-dalek were originally derived from Adam Langley's
+Go ed25519 implementation, found at <https://github.com/agl/ed25519/>,
+under the following licence:
+
+```
+Copyright (c) 2012 The Go Authors. All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are
+met:
+
+   * Redistributions of source code must retain the above copyright
+notice, this list of conditions and the following disclaimer.
+   * Redistributions in binary form must reproduce the above
+copyright notice, this list of conditions and the following disclaimer
+in the documentation and/or other materials provided with the
+distribution.
+   * Neither the name of Google Inc. nor the names of its
+contributors may be used to endorse or promote products derived from
+this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
+OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+```
+
+# Lizard
+
+The `src/lizard` directory was copied from [Signal's curve25519-dalek repo]( https://github.com/signalapp/curve25519-dalek/tree/7c6d34756355a3566a704da84dce7b1c039a6572). Its license is copied below
+
+```
+MIT License
+
+Copyright (c) 2019 Bas Westerbaan
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
+```
diff --git a/curve25519-dalek/LICENSE b/curve25519-dalek/LICENSE
index ff3475745..987ee5311 100644
--- a/curve25519-dalek/LICENSE
+++ b/curve25519-dalek/LICENSE
@@ -26,40 +26,4 @@ TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
-SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 
-
-========================================================================
-
-Portions of curve25519-dalek were originally derived from Adam Langley's
-Go ed25519 implementation, found at <https://github.com/agl/ed25519/>,
-under the following licence:
-
-========================================================================
-
-Copyright (c) 2012 The Go Authors. All rights reserved.
-
-Redistribution and use in source and binary forms, with or without
-modification, are permitted provided that the following conditions are
-met:
-
-   * Redistributions of source code must retain the above copyright
-notice, this list of conditions and the following disclaimer.
-   * Redistributions in binary form must reproduce the above
-copyright notice, this list of conditions and the following disclaimer
-in the documentation and/or other materials provided with the
-distribution.
-   * Neither the name of Google Inc. nor the names of its
-contributors may be used to endorse or promote products derived from
-this software without specific prior written permission.
-
-THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
-IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
-TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
-PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
-OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
-EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
-PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
-PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
-LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
-NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
diff --git a/curve25519-dalek/src/lizard/LICENSE b/curve25519-dalek/src/lizard/LICENSE
deleted file mode 100644
index 2df7c9688..000000000
--- a/curve25519-dalek/src/lizard/LICENSE
+++ /dev/null
@@ -1,21 +0,0 @@
-MIT License
-
-Copyright (c) 2019 Bas Westerbaan
-
-Permission is hereby granted, free of charge, to any person obtaining a copy
-of this software and associated documentation files (the "Software"), to deal
-in the Software without restriction, including without limitation the rights
-to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
-copies of the Software, and to permit persons to whom the Software is
-furnished to do so, subject to the following conditions:
-
-The above copyright notice and this permission notice shall be included in all
-copies or substantial portions of the Software.
-
-THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
-SOFTWARE.

From 2b821c4031c9a3f0c6124cd2be97cc0f06548974 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 17:30:49 -0400
Subject: [PATCH 05/20] Add invalid Lizard encoding test

---
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 15 +++++++++++++++
 1 file changed, 15 insertions(+)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index a212a159f..8f5c98047 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -282,6 +282,21 @@ mod test {
         );
     }
 
+    // Tests that lizard_decode of a random point is None
+    // TODO: what's the false positive rate on this?
+    #[test]
+    fn test_lizard_invalid() {
+        let mut rng = rand::rng();
+        for _ in 0..100 {
+            let pt = RistrettoPoint::random(&mut rng);
+            assert!(
+                pt.lizard_decode::<Sha256>().is_none(),
+                "random point {:02x?} is a valid Lizard encoding",
+                pt.compress().to_bytes()
+            )
+        }
+    }
+
     #[test]
     fn test_elligator_inv() {
         let mut rng = rand::rng();

From c3e13b950b650540adb456e6e1489dd68ce49263 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 17:31:16 -0400
Subject: [PATCH 06/20] Fix nostd build

---
 curve25519-dalek/Cargo.toml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/curve25519-dalek/Cargo.toml b/curve25519-dalek/Cargo.toml
index d78d0ce43..592c7a306 100644
--- a/curve25519-dalek/Cargo.toml
+++ b/curve25519-dalek/Cargo.toml
@@ -82,7 +82,7 @@ legacy_compatibility = []
 group = ["dep:group", "rand_core"]
 group-bits = ["group", "ff/bits"]
 digest = ["dep:digest"]
-lizard = []
+lizard = ["digest"]
 
 [target.'cfg(all(not(curve25519_dalek_backend = "fiat"), not(curve25519_dalek_backend = "serial"), target_arch = "x86_64"))'.dependencies]
 curve25519-dalek-derive = "0.1"

From a71a6f6962361dc7c44d08f2ed1bae03fc27f94c Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Mon, 11 Aug 2025 17:33:04 -0400
Subject: [PATCH 07/20] Remove unnecessary pragmas and imports; improve docs;
 add personal notes

---
 .../src/lizard/lizard_ristretto.rs            | 96 +++++++------------
 curve25519-dalek/src/lizard/mod.rs            |  4 +-
 2 files changed, 34 insertions(+), 66 deletions(-)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 8f5c98047..c18839d31 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -1,8 +1,7 @@
 //! Defines additional methods on RistrettoPoint for Lizard
 
-#![allow(non_snake_case)]
-
 use digest::Digest;
+use digest::HashMarker;
 use digest::consts::U32;
 
 use crate::constants;
@@ -12,25 +11,18 @@ use subtle::Choice;
 use subtle::ConditionallySelectable;
 use subtle::ConstantTimeEq;
 
-use crate::edwards::EdwardsPoint;
-
 use super::jacobi_quartic::JacobiPoint;
 use super::lizard_constants;
 
 use crate::ristretto::RistrettoPoint;
-#[allow(unused_imports)]
-use core::prelude::*;
 
 impl RistrettoPoint {
-    /// Directly encode 253 bits as a RistrettoPoint, using Elligator
-    pub fn from_uniform_bytes_single_elligator(bytes: &[u8; 32]) -> RistrettoPoint {
-        RistrettoPoint::elligator_ristretto_flavor(&FieldElement::from_bytes(bytes))
-    }
-
-    /// Encode 16 bytes of data to a RistrettoPoint, using the Lizard method
+    /// Encode 16 bytes of data to a RistrettoPoint, using the Lizard method. The secure hash
+    /// function `D` is used internally to produce a unique encoding for `data. Use SHA-256 if
+    /// otherwise unsure.
     pub fn lizard_encode<D: Digest>(data: &[u8; 16]) -> RistrettoPoint
     where
-        D: Digest<OutputSize = U32>,
+        D: Digest<OutputSize = U32> + HashMarker,
     {
         let mut fe_bytes: [u8; 32] = Default::default();
 
@@ -43,10 +35,13 @@ impl RistrettoPoint {
         RistrettoPoint::elligator_ristretto_flavor(&fe)
     }
 
-    /// Decode 16 bytes of data from a RistrettoPoint, using the Lizard method
+    /// Decode 16 bytes of data from a RistrettoPoint, using the Lizard method. Returns `None` if
+    /// this point was not generated using Lizard.
+    // TODO: This also fails if more than one inverse exists. How likely is this? Should we test for
+    // this in the encoding step?
     pub fn lizard_decode<D: Digest>(&self) -> Option<[u8; 16]>
     where
-        D: Digest<OutputSize = U32>,
+        D: Digest<OutputSize = U32> + HashMarker,
     {
         let mut result: [u8; 16] = Default::default();
         let mut h: [u8; 32] = Default::default();
@@ -68,44 +63,12 @@ impl RistrettoPoint {
         if n_found == 1 { Some(result) } else { None }
     }
 
-    /// Directly encode 253 bits as a RistrettoPoint, using Elligator
-    pub fn encode_253_bits(data: &[u8; 32]) -> Option<RistrettoPoint> {
-        if data.len() != 32 {
-            return None;
-        }
-
-        let fe = FieldElement::from_bytes(data);
-        let p = RistrettoPoint::elligator_ristretto_flavor(&fe);
-        Some(p)
-    }
-
-    /// Directly decode a RistrettoPoint as 253 bits, using Elligator
-    pub fn decode_253_bits(&self) -> (u8, [[u8; 32]; 8]) {
-        let mut ret = [[0u8; 32]; 8];
-        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
-
-        for j in 0..8 {
-            ret[j] = fes[j].to_bytes();
-        }
-        (mask, ret)
-    }
-
-    /// Return the coset self + E[4], for debugging.
-    pub fn xcoset4(&self) -> [EdwardsPoint; 4] {
-        [
-            self.0,
-            self.0 + constants::EIGHT_TORSION[2],
-            self.0 + constants::EIGHT_TORSION[4],
-            self.0 + constants::EIGHT_TORSION[6],
-        ]
-    }
-
     /// Computes the at most 8 positive FieldElements f such that
-    /// self == elligator_ristretto_flavor(f).
+    /// `self == RistrettoPoint::elligator_ristretto_flavor(f)`.
     /// Assumes self is even.
     ///
     /// Returns a bitmask of which elements in fes are set.
-    pub fn elligator_ristretto_flavor_inverse(&self) -> (u8, [FieldElement; 8]) {
+    fn elligator_ristretto_flavor_inverse(&self) -> (u8, [FieldElement; 8]) {
         // Elligator2 computes a Point from a FieldElement in two steps: first
         // it computes a (s,t) on the Jacobi quartic and then computes the
         // corresponding even point on the Edwards curve.
@@ -226,25 +189,30 @@ impl RistrettoPoint {
     }
 }
 
-// ------------------------------------------------------------------------
-// Tests
-// ------------------------------------------------------------------------
-
 #[cfg(test)]
 mod test {
-
-    use sha2;
-
-    use self::sha2::Sha256;
     use super::*;
-    use crate::ristretto::CompressedRistretto;
+    use crate::{edwards::EdwardsPoint, ristretto::CompressedRistretto};
     use rand_core::RngCore;
+    use sha2::Sha256;
+
+    /// Return the coset self + E\[4\]
+    fn xcoset4(pt: &RistrettoPoint) -> [EdwardsPoint; 4] {
+        [
+            pt.0,
+            pt.0 + constants::EIGHT_TORSION[2],
+            pt.0 + constants::EIGHT_TORSION[4],
+            pt.0 + constants::EIGHT_TORSION[6],
+        ]
+    }
 
-    fn test_lizard_encode_helper(data: &[u8; 16], result: &[u8; 32]) {
+    /// Checks
+    /// `lizard_decode(lizard_encode(data)) == lizard_decode(expected_pt_bytes) == data`
+    fn test_lizard_encode_helper(data: &[u8; 16], expected_pt_bytes: &[u8; 32]) {
         let p = RistrettoPoint::lizard_encode::<Sha256>(data);
-        let p_bytes = p.compress().to_bytes();
-        assert!(&p_bytes == result);
-        let p = CompressedRistretto::from_slice(&p_bytes)
+        let pt_bytes = p.compress().to_bytes();
+        assert!(&pt_bytes == expected_pt_bytes);
+        let p = CompressedRistretto::from_slice(&pt_bytes)
             .unwrap()
             .decompress()
             .unwrap();
@@ -322,8 +290,8 @@ mod test {
             let fe = FieldElement::from_bytes(&fe_bytes);
 
             let pt = RistrettoPoint::elligator_ristretto_flavor(&fe);
-            for pt2 in &pt.xcoset4() {
-                let (mask, fes) = RistrettoPoint(*pt2).elligator_ristretto_flavor_inverse();
+            for pt2 in xcoset4(&pt) {
+                let (mask, fes) = RistrettoPoint(pt2).elligator_ristretto_flavor_inverse();
 
                 let mut found = false;
                 for (j, fe_j) in fes.iter().enumerate() {
diff --git a/curve25519-dalek/src/lizard/mod.rs b/curve25519-dalek/src/lizard/mod.rs
index df1ac4506..45b4aee08 100644
--- a/curve25519-dalek/src/lizard/mod.rs
+++ b/curve25519-dalek/src/lizard/mod.rs
@@ -8,6 +8,6 @@ mod u32_constants;
 #[cfg(curve25519_dalek_bits = "64")]
 mod u64_constants;
 
-pub mod jacobi_quartic;
-pub mod lizard_constants;
+mod jacobi_quartic;
+mod lizard_constants;
 pub mod lizard_ristretto;

From 502d71b15b10e7a8313df0d5686b723ec2dc96da Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Wed, 13 Aug 2025 20:54:40 -0400
Subject: [PATCH 08/20] Make clippy happy; remove old TODOs

---
 curve25519-dalek/build.rs                       | 2 +-
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 8 +++-----
 2 files changed, 4 insertions(+), 6 deletions(-)

diff --git a/curve25519-dalek/build.rs b/curve25519-dalek/build.rs
index 0460b8c13..b8a3fb023 100644
--- a/curve25519-dalek/build.rs
+++ b/curve25519-dalek/build.rs
@@ -17,7 +17,7 @@ impl std::fmt::Display for DalekBits {
             DalekBits::Dalek32 => "32",
             DalekBits::Dalek64 => "64",
         };
-        write!(f, "{}", w_bits)
+        write!(f, "{w_bits}")
     }
 }
 
diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index c18839d31..6f7a99ca2 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -20,7 +20,7 @@ impl RistrettoPoint {
     /// Encode 16 bytes of data to a RistrettoPoint, using the Lizard method. The secure hash
     /// function `D` is used internally to produce a unique encoding for `data. Use SHA-256 if
     /// otherwise unsure.
-    pub fn lizard_encode<D: Digest>(data: &[u8; 16]) -> RistrettoPoint
+    pub fn lizard_encode<D>(data: &[u8; 16]) -> RistrettoPoint
     where
         D: Digest<OutputSize = U32> + HashMarker,
     {
@@ -37,9 +37,7 @@ impl RistrettoPoint {
 
     /// Decode 16 bytes of data from a RistrettoPoint, using the Lizard method. Returns `None` if
     /// this point was not generated using Lizard.
-    // TODO: This also fails if more than one inverse exists. How likely is this? Should we test for
-    // this in the encoding step?
-    pub fn lizard_decode<D: Digest>(&self) -> Option<[u8; 16]>
+    pub fn lizard_decode<D>(&self) -> Option<[u8; 16]>
     where
         D: Digest<OutputSize = U32> + HashMarker,
     {
@@ -60,6 +58,7 @@ impl RistrettoPoint {
             }
             n_found += ok.unwrap_u8();
         }
+        // Per the README, the likelihood that n_found == 2 is something like 2^-122
         if n_found == 1 { Some(result) } else { None }
     }
 
@@ -251,7 +250,6 @@ mod test {
     }
 
     // Tests that lizard_decode of a random point is None
-    // TODO: what's the false positive rate on this?
     #[test]
     fn test_lizard_invalid() {
         let mut rng = rand::rng();

From 4d6800bb41c6da5e26c2c6ece44a7f0c94a23812 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Wed, 13 Aug 2025 21:06:18 -0400
Subject: [PATCH 09/20] Added lizard details to readme and changelog

---
 curve25519-dalek/CHANGELOG.md         | 4 ++++
 curve25519-dalek/README.md            | 1 +
 curve25519-dalek/src/lizard/README.md | 1 +
 3 files changed, 6 insertions(+)
 create mode 100644 curve25519-dalek/src/lizard/README.md

diff --git a/curve25519-dalek/CHANGELOG.md b/curve25519-dalek/CHANGELOG.md
index b126651ad..bc54382ab 100644
--- a/curve25519-dalek/CHANGELOG.md
+++ b/curve25519-dalek/CHANGELOG.md
@@ -5,6 +5,10 @@ major series.
 
 ## 5.x series
 
+## Unreleased
+
+* Add Lizard bytes-to-point injection for Ristretto. Gated under `lizard`.
+
 ## 5.0.0-pre.0
 
 * Update edition to 2024
diff --git a/curve25519-dalek/README.md b/curve25519-dalek/README.md
index 1316dbac8..765cd6da8 100644
--- a/curve25519-dalek/README.md
+++ b/curve25519-dalek/README.md
@@ -57,6 +57,7 @@ curve25519-dalek = ">= 5.0, < 5.2"
 | `legacy_compatibility`|       | Enables `Scalar::from_bits`, which allows the user to build unreduced scalars whose arithmetic is broken. Do not use this unless you know what you're doing. |
 | `group`            |          | Enables external `group` and `ff` crate traits. |
 | `group-bits`       |          | Enables `group` and impls `ff::PrimeFieldBits` for `Scalar`.  |
+| `lizard`           |          | Enables the [Lizard](src/lizard/README.md) bytestring-to-point injection for `RistrettoPoint`. Specifically enables the methods `lizard_encode` and `lizard_decode`. |
 
 To disable the default features when using `curve25519-dalek` as a dependency,
 add `default-features = false` to the dependency in your `Cargo.toml`. To
diff --git a/curve25519-dalek/src/lizard/README.md b/curve25519-dalek/src/lizard/README.md
new file mode 100644
index 000000000..3b94f9157
--- /dev/null
+++ b/curve25519-dalek/src/lizard/README.md
@@ -0,0 +1 @@
+Placeholder

From 158ea3ea0f51f5f60193c4ecd7d755821328504b Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Thu, 14 Aug 2025 11:56:11 -0400
Subject: [PATCH 10/20] Upgrade ristretto.sage to Python3

---
 curve25519-dalek/vendor/ristretto.sage | 46 +++++++++++++-------------
 1 file changed, 23 insertions(+), 23 deletions(-)

diff --git a/curve25519-dalek/vendor/ristretto.sage b/curve25519-dalek/vendor/ristretto.sage
index 085fef6d1..fc9d87906 100644
--- a/curve25519-dalek/vendor/ristretto.sage
+++ b/curve25519-dalek/vendor/ristretto.sage
@@ -6,7 +6,7 @@ class SpecException(Exception): pass
 def lobit(x): return int(x) & 1
 def hibit(x): return lobit(2*x)
 def negative(x): return lobit(x)
-def enc_le(x,n): return bytearray([int(x)>>(8*i) & 0xFF for i in xrange(n)])
+def enc_le(x,n): return bytearray([int(x)>>(8*i) & 0xFF for i in range(n)])
 def dec_le(x): return sum(b<<(8*i) for i,b in enumerate(x))
 def randombytes(n): return bytearray([randint(0,255) for _ in range(n)])
 
@@ -463,8 +463,8 @@ class Decaf_1_1_Point(QuotientEdwardsPoint):
                         if negative(sr) != toggle_r: sr = -sr
                         ret = self.gfToBytes(sr)
                         if self.elligator(ret) != self and self.elligator(ret) != -self:
-                            print "WRONG!",[toggle_rotation,toggle_altx,toggle_s]
-                        if self.elligator(ret) == -self and self != -self: print "Negated!",[toggle_rotation,toggle_altx,toggle_s]
+                            print("WRONG!",[toggle_rotation,toggle_altx,toggle_s])
+                        if self.elligator(ret) == -self and self != -self: print("Negated!",[toggle_rotation,toggle_altx,toggle_s])
                         rets.append(bytes(ret))
         return rets
 
@@ -590,7 +590,7 @@ class Decaf_1_1_Point(QuotientEdwardsPoint):
     def elligatorInverseBruteForce(self):
         """Invert Elligator using SAGE's polynomial solver"""
         a,d = self.a,self.d
-        R.<r0> = self.F[]
+        R = self.F["r0"]
         r = self.qnr * r0^2
         den = (d*r-(d-a))*((d-a)*r-d)
         n1 = (r+1)*(a-2*d)/den
@@ -602,7 +602,7 @@ class Decaf_1_1_Point(QuotientEdwardsPoint):
             y = (1-a*s2) / t
 
             selfT = self
-            for i in xrange(self.cofactor/2):
+            for i in range(self.cofactor/2):
                 xT,yT = selfT
                 polyX = xT^2-x2
                 polyY = yT-y
@@ -721,7 +721,7 @@ class IsoEd25519Point(Decaf_1_1_Point):
 class TestFailedException(Exception): pass
 
 def test(cls,n):
-    print "Testing curve %s" % cls.__name__
+    print("Testing curve %s" % cls.__name__)
     
     specials = [1]
     ii = cls.F(-1)
@@ -744,7 +744,7 @@ def test(cls,n):
     
     P = cls.base()
     Q = cls()
-    for i in xrange(n):
+    for i in range(n):
         #print binascii.hexlify(Q.encode())
         QE = Q.encode()
         QQ = cls.decode(QE)
@@ -766,7 +766,7 @@ def test(cls,n):
                     raise TestFailedException("s -> 1/s should work for cofactor 4")
         
         QT = Q
-        for h in xrange(cls.cofactor):
+        for h in range(cls.cofactor):
             QT = QT.torque()
             if QT.encode() != QE:
                 raise TestFailedException("Can't torque %s,%d" % (str(Q),h+1))
@@ -782,19 +782,19 @@ def test(cls,n):
         Q = Q1
         
 def testElligator(cls,n):
-    print "Testing elligator on %s" % cls.__name__
-    for i in xrange(n):
+    print("Testing elligator on %s" % cls.__name__)
+    for i in range(n):
         r = randombytes(cls.encLen)
         P = cls.elligator(r)
         if hasattr(P,"invertElligator"):
             iv = P.invertElligator()
             modr = bytes(cls.gfToBytes(cls.bytesToGf(r,mustBeProper=False,maskHiBits=True)))
             iv2 = P.torque().invertElligator()
-            if modr not in iv: print "Failed to invert Elligator!"
+            if modr not in iv: print("Failed to invert Elligator!")
             if len(iv) != len(set(iv)):
-                print "Elligator inverses not unique!", len(set(iv)), len(iv)
+                print("Elligator inverses not unique!", len(set(iv)), len(iv))
             if iv != iv2:
-                print "Elligator is untorqueable!"
+                print("Elligator is untorqueable!")
                 #print [binascii.hexlify(j) for j in iv]
                 #print [binascii.hexlify(j) for j in iv2]
                 #break
@@ -802,7 +802,7 @@ def testElligator(cls,n):
             pass # TODO
 
 def gangtest(classes,n):
-    print "Gang test",[cls.__name__ for cls in classes]
+    print("Gang test",[cls.__name__ for cls in classes])
     specials = [1]
     ii = classes[0].F(-1)
     while is_square(ii):
@@ -810,27 +810,27 @@ def gangtest(classes,n):
         ii = sqrt(ii)
     specials.append(ii)
     
-    for i in xrange(n):
+    for i in range(n):
         rets = [bytes((cls.base()*i).encode()) for cls in classes]
         if len(set(rets)) != 1:
-            print "Divergence in encode at %d" % i
+            print("Divergence in encode at %d" % i)
             for c,ret in zip(classes,rets):
-                print c,binascii.hexlify(ret)
-            print
+                print(c,binascii.hexlify(ret))
+            print()
         
         if i < len(specials): r0 = enc_le(specials[i],classes[0].encLen)
         else: r0 = randombytes(classes[0].encLen)
         
         rets = [bytes((cls.elligator(r0)*i).encode()) for cls in classes]
         if len(set(rets)) != 1:
-            print "Divergence in elligator at %d" % i
+            print("Divergence in elligator at %d" % i)
             for c,ret in zip(classes,rets):
-                print c,binascii.hexlify(ret)
-            print
+                print(c,binascii.hexlify(ret))
+            print()
 
 def testDoubleAndEncode(cls,n):
-    print "Testing doubleAndEncode on %s" % cls.__name__
-    for i in xrange(n):
+    print("Testing doubleAndEncode on %s" % cls.__name__)
+    for i in range(n):
         r1 = randombytes(cls.encLen)
         r2 = randombytes(cls.encLen)
         u = cls.elligator(r1) + cls.elligator(r2)

From 18fa45d2fd3cde1ef208523434e519af0de77565 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Thu, 14 Aug 2025 12:06:16 -0400
Subject: [PATCH 11/20] Fix runtime error in ristretto.sage

---
 curve25519-dalek/vendor/ristretto.sage | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/curve25519-dalek/vendor/ristretto.sage b/curve25519-dalek/vendor/ristretto.sage
index fc9d87906..6f06a990a 100644
--- a/curve25519-dalek/vendor/ristretto.sage
+++ b/curve25519-dalek/vendor/ristretto.sage
@@ -22,17 +22,17 @@ def optimized_version_of(spec):
             try: opt_ans = f(self,*args,**kwargs),None
             except Exception as e: opt_ans = None,e
             if spec_ans[1] is None and opt_ans[1] is not None:
-                raise
+                raise opt_ans[1]
                 #raise SpecException("Mismatch in %s: spec returned %s but opt threw %s"
                 #    % (f.__name__,str(spec_ans[0]),str(opt_ans[1])))
             if spec_ans[1] is not None and opt_ans[1] is None:
-                raise
+                raise spec_ans[1]
                 #raise SpecException("Mismatch in %s: spec threw %s but opt returned %s"
                 #    % (f.__name__,str(spec_ans[1]),str(opt_ans[0])))
             if spec_ans[0] != opt_ans[0]:
                 raise SpecException("Mismatch in %s: %s != %s"
                     % (f.__name__,pr(spec_ans[0]),pr(opt_ans[0])))
-            if opt_ans[1] is not None: raise
+            if opt_ans[1] is not None: raise opt_ans[1]
             else: return opt_ans[0]
         wrapper.__name__ = f.__name__
         return wrapper
@@ -133,7 +133,7 @@ class QuotientEdwardsPoint(object):
         s = dec_le(bytes)
         if mustBeProper and s >= cls.F.order():
             raise InvalidEncodingException("%d out of range!" % s)
-        bitlen = int(ceil(log(cls.F.order())/log(2)))
+        bitlen = cls.F.order().bit_length()
         if maskHiBits: s &= 2^bitlen-1
         s = cls.F(s)
         if mustBePositive and negative(s):

From babf4723f618cf2f50495372c58d2463e9f3b782 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Thu, 14 Aug 2025 13:15:30 -0400
Subject: [PATCH 12/20] Add lizard test vector generation to ristretto.sage

---
 curve25519-dalek/vendor/ristretto.sage | 19 +++++++++++++++++++
 1 file changed, 19 insertions(+)

diff --git a/curve25519-dalek/vendor/ristretto.sage b/curve25519-dalek/vendor/ristretto.sage
index 6f06a990a..2148378e4 100644
--- a/curve25519-dalek/vendor/ristretto.sage
+++ b/curve25519-dalek/vendor/ristretto.sage
@@ -1,4 +1,5 @@
 import binascii
+from hashlib import sha256
 class InvalidEncodingException(Exception): pass
 class NotOnCurveException(Exception): pass
 class SpecException(Exception): pass
@@ -836,6 +837,23 @@ def testDoubleAndEncode(cls,n):
         u = cls.elligator(r1) + cls.elligator(r2)
         u.doubleAndEncode()
 
+# Prints test vectors for the Lizard encoding function. Output is the encoding of a compressed
+# Ristretto point
+def testLizard():
+    # 16-byte strings, in hex
+    inputs = ["00000000000000000000000000000000", "01010101010101010101010101010101"]
+
+    for payload in map(binascii.unhexlify, inputs):
+        # Do the lizard encoding of the field element
+        data = bytearray(sha256(payload).digest())
+        data[8:24] = payload
+        data[0] &= 0b11111110
+        data[31] &= 0b00111111
+
+        # Encode to Ristretto (Ed25519Point is actually Ristretto), and print the byte repr
+        pt = Ed25519Point.elligator(data)
+        print("lizard({}) = {}".format(binascii.hexlify(payload), binascii.hexlify(pt.encode())))
+
 testDoubleAndEncode(Ed25519Point,100)
 testDoubleAndEncode(NegEd25519Point,100)
 testDoubleAndEncode(IsoEd25519Point,100)
@@ -853,5 +871,6 @@ testDoubleAndEncode(TwistedEd448GoldilocksPoint,100)
 #testElligator(IsoEd448Point,100)
 #testElligator(Ed448GoldilocksPoint,100)
 #testElligator(TwistedEd448GoldilocksPoint,100)
+#testLizard()
 #gangtest([IsoEd448Point,TwistedEd448GoldilocksPoint,Ed448GoldilocksPoint],100)
 #gangtest([Ed25519Point,IsoEd25519Point],100)

From 9423b5a47c8e497891817b695d453270b5dd6016 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Thu, 14 Aug 2025 15:01:26 -0400
Subject: [PATCH 13/20] Add and clarify lizard test vectors

---
 .../src/lizard/lizard_ristretto.rs            | 59 ++++++++++---------
 curve25519-dalek/vendor/ristretto.sage        |  7 ++-
 2 files changed, 36 insertions(+), 30 deletions(-)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 6f7a99ca2..696303491 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -207,10 +207,13 @@ mod test {
 
     /// Checks
     /// `lizard_decode(lizard_encode(data)) == lizard_decode(expected_pt_bytes) == data`
-    fn test_lizard_encode_helper(data: &[u8; 16], expected_pt_bytes: &[u8; 32]) {
-        let p = RistrettoPoint::lizard_encode::<Sha256>(data);
+    fn test_lizard_encode_helper(data: &[u8], expected_pt_bytes: &[u8]) {
+        assert_eq!(data.len(), 16);
+        assert_eq!(expected_pt_bytes.len(), 32);
+
+        let p = RistrettoPoint::lizard_encode::<Sha256>(data.try_into().unwrap());
         let pt_bytes = p.compress().to_bytes();
-        assert!(&pt_bytes == expected_pt_bytes);
+        assert!(pt_bytes == expected_pt_bytes);
         let p = CompressedRistretto::from_slice(&pt_bytes)
             .unwrap()
             .decompress()
@@ -221,32 +224,30 @@ mod test {
 
     #[test]
     fn test_lizard_encode() {
-        test_lizard_encode_helper(
-            &[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
-            &[
-                0xf0, 0xb7, 0xe3, 0x44, 0x84, 0xf7, 0x4c, 0xf0, 0xf, 0x15, 0x2, 0x4b, 0x73, 0x85,
-                0x39, 0x73, 0x86, 0x46, 0xbb, 0xbe, 0x1e, 0x9b, 0xc7, 0x50, 0x9a, 0x67, 0x68, 0x15,
-                0x22, 0x7e, 0x77, 0x4f,
-            ],
-        );
-
-        test_lizard_encode_helper(
-            &[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1],
-            &[
-                0xcc, 0x92, 0xe8, 0x1f, 0x58, 0x5a, 0xfc, 0x5c, 0xaa, 0xc8, 0x86, 0x60, 0xd8, 0xd1,
-                0x7e, 0x90, 0x25, 0xa4, 0x44, 0x89, 0xa3, 0x63, 0x4, 0x21, 0x23, 0xf6, 0xaf, 0x7,
-                0x2, 0x15, 0x6e, 0x65,
-            ],
-        );
-
-        test_lizard_encode_helper(
-            &[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15],
-            &[
-                0xc8, 0x30, 0x57, 0x3f, 0x8a, 0x8e, 0x77, 0x78, 0x67, 0x1f, 0x76, 0xcd, 0xc7, 0x96,
-                0xdc, 0xa, 0x23, 0x5c, 0xf1, 0x77, 0xf1, 0x97, 0xd9, 0xfc, 0xba, 0x6, 0xe8, 0x4e,
-                0x96, 0x24, 0x74, 0x44,
-            ],
-        );
+        // Test vectors are of the form (x, y) where y is the compressed encoding of the Ristretto
+        // point given by lizard_encode(x).
+        // These values come from the testLizard() function in vendor/ristretto.sage
+        let test_vectors = [
+            (
+                "00000000000000000000000000000000",
+                "f0b7e34484f74cf00f15024b738539738646bbbe1e9bc7509a676815227e774f",
+            ),
+            (
+                "01010101010101010101010101010101",
+                "cc92e81f585afc5caac88660d8d17e9025a44489a363042123f6af0702156e65",
+            ),
+            (
+                "000102030405060708090a0b0c0d0e0f",
+                "c830573f8a8e7778671f76cdc796dc0a235cf177f197d9fcba06e84e96247444",
+            ),
+            (
+                "dddddddddddddddddddddddddddddddd",
+                "ccb60554c081841037f821fa827b6a5bc2531f80e2647f1a858611f4ccfe3056",
+            ),
+        ];
+        for tv in test_vectors {
+            test_lizard_encode_helper(&hex::decode(tv.0).unwrap(), &hex::decode(tv.1).unwrap());
+        }
     }
 
     // Tests that lizard_decode of a random point is None
diff --git a/curve25519-dalek/vendor/ristretto.sage b/curve25519-dalek/vendor/ristretto.sage
index 2148378e4..e518c0ae0 100644
--- a/curve25519-dalek/vendor/ristretto.sage
+++ b/curve25519-dalek/vendor/ristretto.sage
@@ -841,7 +841,12 @@ def testDoubleAndEncode(cls,n):
 # Ristretto point
 def testLizard():
     # 16-byte strings, in hex
-    inputs = ["00000000000000000000000000000000", "01010101010101010101010101010101"]
+    inputs = [
+        "00000000000000000000000000000000",
+        "01010101010101010101010101010101",
+        "000102030405060708090a0b0c0d0e0f",
+        "dddddddddddddddddddddddddddddddd",
+    ]
 
     for payload in map(binascii.unhexlify, inputs):
         # Do the lizard encoding of the field element

From e7bba8e7b9ef4ec3d14d1d19b5f4bdc8b615b2f4 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Thu, 14 Aug 2025 17:15:58 -0400
Subject: [PATCH 14/20] Small cleanup

---
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 696303491..01e5d76a8 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -29,7 +29,7 @@ impl RistrettoPoint {
         let digest = D::digest(data);
         fe_bytes[0..32].copy_from_slice(digest.as_slice());
         fe_bytes[8..24].copy_from_slice(data);
-        fe_bytes[0] &= 254; // make positive since Elligator on r and -r is the same
+        fe_bytes[0] &= 254;
         fe_bytes[31] &= 63;
         let fe = FieldElement::from_bytes(&fe_bytes);
         RistrettoPoint::elligator_ristretto_flavor(&fe)
@@ -47,8 +47,8 @@ impl RistrettoPoint {
         let mut n_found = 0;
         for (j, fe_j) in fes.iter().enumerate() {
             let mut ok = Choice::from((mask >> j) & 1);
-            let buf2 = fe_j.to_bytes(); // array
-            h.copy_from_slice(&D::digest(&buf2[8..24])); // array
+            let buf2 = fe_j.to_bytes();
+            h.copy_from_slice(&D::digest(&buf2[8..24]));
             h[8..24].copy_from_slice(&buf2[8..24]);
             h[0] &= 254;
             h[31] &= 63;

From 7a73437841c381f1fa00d99448de96b469cce986 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Tue, 2 Sep 2025 23:26:32 -0400
Subject: [PATCH 15/20] Rename Lizard e^-1 function to e_inv

---
 curve25519-dalek/src/lizard/jacobi_quartic.rs   | 9 +++++----
 curve25519-dalek/src/lizard/lizard_ristretto.rs | 4 ++--
 2 files changed, 7 insertions(+), 6 deletions(-)

diff --git a/curve25519-dalek/src/lizard/jacobi_quartic.rs b/curve25519-dalek/src/lizard/jacobi_quartic.rs
index af64b4580..34d6189f5 100644
--- a/curve25519-dalek/src/lizard/jacobi_quartic.rs
+++ b/curve25519-dalek/src/lizard/jacobi_quartic.rs
@@ -22,12 +22,12 @@ pub struct JacobiPoint {
 }
 
 impl JacobiPoint {
-    /// Elligator2 is defined in two steps: first a field element is converted
+    /// Elligator2 is defined in two steps: first a function e maps a field element
     /// to a point (s,t) on the Jacobi quartic associated to the Edwards curve.
     /// Then this point is mapped to a point on the Edwards curve.
-    /// This function computes a field element that is mapped to a given (s,t)
-    /// with Elligator2 if it exists.
-    pub(crate) fn elligator_inv(&self) -> (Choice, FieldElement) {
+    /// This function computes a field element that is mapped by e to a given (s,t),
+    /// if it exists.
+    pub(crate) fn e_inv(&self) -> (Choice, FieldElement) {
         let mut out = FieldElement::ZERO;
 
         // Special case: s = 0.  If s is zero, either t = 1 or t = -1.
@@ -56,6 +56,7 @@ impl JacobiPoint {
         let mut pms2 = s2;
         pms2.conditional_negate(self.S.is_negative());
         let mut x = &(&a + &pms2) * &y;
+        // Always pick the positive solution
         let x_is_negative = x.is_negative();
         x.conditional_negate(x_is_negative);
         out.conditional_assign(&x, !done);
diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 01e5d76a8..5c73b19c3 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -92,14 +92,14 @@ impl RistrettoPoint {
         let mut ret = [FieldElement::ONE; 8];
 
         for i in 0..4 {
-            let (ok, fe) = jcs[i].elligator_inv();
+            let (ok, fe) = jcs[i].e_inv();
             let mut tmp: u8 = 0;
             ret[2 * i] = fe;
             tmp.conditional_assign(&1, ok);
             mask |= tmp << (2 * i);
 
             let jc = jcs[i].dual();
-            let (ok, fe) = jc.elligator_inv();
+            let (ok, fe) = jc.e_inv();
             let mut tmp: u8 = 0;
             ret[2 * i + 1] = fe;
             tmp.conditional_assign(&1, ok);

From 2ce67cc17d1b9bebdff05bc9b39da526f60d2a44 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Wed, 3 Sep 2025 00:39:43 -0400
Subject: [PATCH 16/20] Added map_to_curve and its inverse

---
 .../src/lizard/lizard_ristretto.rs            | 76 +++++++++++++++++--
 curve25519-dalek/src/ristretto.rs             |  6 +-
 2 files changed, 74 insertions(+), 8 deletions(-)

diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 5c73b19c3..943e13519 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -29,8 +29,8 @@ impl RistrettoPoint {
         let digest = D::digest(data);
         fe_bytes[0..32].copy_from_slice(digest.as_slice());
         fe_bytes[8..24].copy_from_slice(data);
-        fe_bytes[0] &= 254;
-        fe_bytes[31] &= 63;
+        fe_bytes[0] &= 0b11111110;
+        fe_bytes[31] &= 0b00111111;
         let fe = FieldElement::from_bytes(&fe_bytes);
         RistrettoPoint::elligator_ristretto_flavor(&fe)
     }
@@ -50,8 +50,8 @@ impl RistrettoPoint {
             let buf2 = fe_j.to_bytes();
             h.copy_from_slice(&D::digest(&buf2[8..24]));
             h[8..24].copy_from_slice(&buf2[8..24]);
-            h[0] &= 254;
-            h[31] &= 63;
+            h[0] &= 0b11111110;
+            h[31] &= 0b00111111;
             ok &= h.ct_eq(&buf2);
             for i in 0..16 {
                 result[i] = u8::conditional_select(&result[i], &buf2[8 + i], ok);
@@ -66,7 +66,8 @@ impl RistrettoPoint {
     /// `self == RistrettoPoint::elligator_ristretto_flavor(f)`.
     /// Assumes self is even.
     ///
-    /// Returns a bitmask of which elements in fes are set.
+    /// First return value is a bitmask of which elements are valid inverses (lowest bit corresponds
+    /// to the 0-th element).
     fn elligator_ristretto_flavor_inverse(&self) -> (u8, [FieldElement; 8]) {
         // Elligator2 computes a Point from a FieldElement in two steps: first
         // it computes a (s,t) on the Jacobi quartic and then computes the
@@ -186,6 +187,39 @@ impl RistrettoPoint {
             JacobiPoint { S: s3, T: t3 },
         ]
     }
+
+    /// Interprets the given bytestring as a positive field element and computes the Ristretto
+    /// Elligator map. Note this clears the bottom bit and top two bits of `bytes`.
+    ///
+    /// # Warning
+    ///
+    /// This function does not produce cryptographically random-looking Ristretto points. Use
+    /// [`Self::hash_from_bytes`] for that. DO NOT USE THIS FUNCTION unless you really know what
+    /// you're doing.
+    pub fn map_to_curve(mut bytes: [u8; 32]) -> RistrettoPoint {
+        // We only have a meaningful inverse if we give Elligator a point in its domain, ie a
+        // positive (meaning low bit 0) field element. Mask off the top two bits to ensure it's less
+        // than the modulus, and the bottom bit for evenness.
+        bytes[0] &= 0b11111110;
+        bytes[31] &= 0b00111111;
+
+        let fe = FieldElement::from_bytes(&bytes);
+        RistrettoPoint::elligator_ristretto_flavor(&fe)
+    }
+
+    /// Computes the possible bytestrings that could have produced this point via
+    /// [`Self::map_to_curve`].
+    /// First return value is a bitmask of which elements are valid inverses (lowest bit corresponds
+    /// to the 0-th element).
+    pub fn map_to_curve_inverse(&self) -> (u8, [[u8; 32]; 8]) {
+        let mut ret = [[0u8; 32]; 8];
+        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
+
+        for j in 0..8 {
+            ret[j] = fes[j].to_bytes();
+        }
+        (mask, ret)
+    }
 }
 
 #[cfg(test)]
@@ -305,4 +339,36 @@ mod test {
             }
         }
     }
+
+    // Tests that map_to_curve_inverse ○ map_to_curve is the identity
+    #[test]
+    fn test_map_to_curve_inverse() {
+        let mut rng = rand::rng();
+
+        for _ in 0..100 {
+            let mut input = [0u8; 32];
+            rng.fill_bytes(&mut input);
+
+            // Map to Ristretto and invert it
+            let pt = RistrettoPoint::map_to_curve(input);
+            let (bitmask, inverses) = pt.map_to_curve_inverse();
+
+            // map_to_curve masks the bottom bit and top two bits of `input`
+            let mut expected_inverse = input;
+            expected_inverse[31] &= 0b00111111;
+            expected_inverse[0] &= 0b11111110;
+
+            // Check that one of the valid inverses matches the input
+            let mut found = false;
+            for (i, inv) in inverses.into_iter().enumerate() {
+                let is_valid = bitmask & (1 << i) != 0;
+                if is_valid && inv == expected_inverse {
+                    found = true;
+                }
+            }
+            if !found {
+                panic!("did not find inverse for input {:02x?}", input);
+            }
+        }
+    }
 }
diff --git a/curve25519-dalek/src/ristretto.rs b/curve25519-dalek/src/ristretto.rs
index ca5fa88ec..bc1ccd515 100644
--- a/curve25519-dalek/src/ristretto.rs
+++ b/curve25519-dalek/src/ristretto.rs
@@ -656,9 +656,9 @@ impl RistrettoPoint {
         ]
     }
 
-    /// Computes the Ristretto Elligator map. This is the
-    /// [`MAP`](https://datatracker.ietf.org/doc/html/draft-irtf-cfrg-ristretto255-decaf448-04#section-4.3.4)
-    /// function defined in the Ristretto spec.
+    /// Computes the Ristretto Elligator map for the given field element. This is the second half of
+    /// the [`MAP`](https://www.rfc-editor.org/rfc/rfc9496.html#section-4.3.4-4) function defined in
+    /// the Ristretto spec.
     ///
     /// # Note
     ///

From 1d427571fd001c09a748941a3ea7d6401b372478 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Fri, 5 Sep 2025 16:23:18 -0400
Subject: [PATCH 17/20] Use CtOption instead of bitmasks; other cleanup

---
 curve25519-dalek/src/field.rs                 |   6 +
 curve25519-dalek/src/lizard/jacobi_quartic.rs |   9 +-
 .../src/lizard/lizard_ristretto.rs            | 125 +++++++++---------
 3 files changed, 71 insertions(+), 69 deletions(-)

diff --git a/curve25519-dalek/src/field.rs b/curve25519-dalek/src/field.rs
index db2d9aa4b..3e5631107 100644
--- a/curve25519-dalek/src/field.rs
+++ b/curve25519-dalek/src/field.rs
@@ -98,6 +98,12 @@ impl ConstantTimeEq for FieldElement {
     }
 }
 
+impl Default for FieldElement {
+    fn default() -> Self {
+        FieldElement::ZERO
+    }
+}
+
 impl FieldElement {
     /// Load a `FieldElement` from 64 bytes, by reducing modulo q.
     #[cfg(feature = "digest")]
diff --git a/curve25519-dalek/src/lizard/jacobi_quartic.rs b/curve25519-dalek/src/lizard/jacobi_quartic.rs
index 34d6189f5..897f57074 100644
--- a/curve25519-dalek/src/lizard/jacobi_quartic.rs
+++ b/curve25519-dalek/src/lizard/jacobi_quartic.rs
@@ -6,6 +6,7 @@ use subtle::Choice;
 use subtle::ConditionallyNegatable;
 use subtle::ConditionallySelectable;
 use subtle::ConstantTimeEq;
+use subtle::CtOption;
 
 use super::lizard_constants;
 use crate::constants;
@@ -27,7 +28,7 @@ impl JacobiPoint {
     /// Then this point is mapped to a point on the Edwards curve.
     /// This function computes a field element that is mapped by e to a given (s,t),
     /// if it exists.
-    pub(crate) fn e_inv(&self) -> (Choice, FieldElement) {
+    pub(crate) fn e_inv(&self) -> CtOption<FieldElement> {
         let mut out = FieldElement::ZERO;
 
         // Special case: s = 0.  If s is zero, either t = 1 or t = -1.
@@ -35,7 +36,7 @@ impl JacobiPoint {
         let s_is_zero = self.S.is_zero();
         let t_equals_one = self.T.ct_eq(&FieldElement::ONE);
         out.conditional_assign(&lizard_constants::SQRT_ID, t_equals_one);
-        let mut ret = s_is_zero;
+        let mut is_defined = s_is_zero;
         let mut done = s_is_zero;
 
         // a := (t+1) (d+1)/(d-1)
@@ -49,7 +50,7 @@ impl JacobiPoint {
 
         // There is no preimage if the square root of i*(s^4-a^2) does not exist.
         let (sq, y) = invSqY.invsqrt();
-        ret |= sq;
+        is_defined |= sq;
         done |= !sq;
 
         // x := (a + sign(s)*s^2) y
@@ -61,7 +62,7 @@ impl JacobiPoint {
         x.conditional_negate(x_is_negative);
         out.conditional_assign(&x, !done);
 
-        (ret, out)
+        CtOption::new(out, is_defined)
     }
 
     pub(crate) fn dual(&self) -> JacobiPoint {
diff --git a/curve25519-dalek/src/lizard/lizard_ristretto.rs b/curve25519-dalek/src/lizard/lizard_ristretto.rs
index 943e13519..8456616ce 100644
--- a/curve25519-dalek/src/lizard/lizard_ristretto.rs
+++ b/curve25519-dalek/src/lizard/lizard_ristretto.rs
@@ -1,13 +1,15 @@
 //! Defines additional methods on RistrettoPoint for Lizard
 
-use digest::Digest;
-use digest::HashMarker;
-use digest::consts::U32;
+use digest::{
+    Digest, HashMarker,
+    array::Array,
+    consts::{U8, U32},
+};
+use subtle::CtOption;
 
 use crate::constants;
 use crate::field::FieldElement;
 
-use subtle::Choice;
 use subtle::ConditionallySelectable;
 use subtle::ConstantTimeEq;
 
@@ -42,19 +44,24 @@ impl RistrettoPoint {
         D: Digest<OutputSize = U32> + HashMarker,
     {
         let mut result: [u8; 16] = Default::default();
-        let mut h: [u8; 32] = Default::default();
-        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
+        let fes = self.elligator_ristretto_flavor_inverse();
         let mut n_found = 0;
-        for (j, fe_j) in fes.iter().enumerate() {
-            let mut ok = Choice::from((mask >> j) & 1);
-            let buf2 = fe_j.to_bytes();
-            h.copy_from_slice(&D::digest(&buf2[8..24]));
-            h[8..24].copy_from_slice(&buf2[8..24]);
-            h[0] &= 0b11111110;
-            h[31] &= 0b00111111;
-            ok &= h.ct_eq(&buf2);
+        for fe in fes {
+            let mut ok = fe.is_some();
+            let fe = fe.unwrap_or(FieldElement::ZERO);
+            let bytes = fe.to_bytes();
+
+            let mut expected_bytes: [u8; 32] = Default::default();
+            expected_bytes.copy_from_slice(&D::digest(&bytes[8..24]));
+            expected_bytes[8..24].copy_from_slice(&bytes[8..24]);
+            expected_bytes[0] &= 0b11111110;
+            expected_bytes[31] &= 0b00111111;
+
+            // If we found our inverse (there's at most 1), write the value to our output buffer
+            ok &= expected_bytes.ct_eq(&bytes);
             for i in 0..16 {
-                result[i] = u8::conditional_select(&result[i], &buf2[8 + i], ok);
+                // Copy bytes 8-24, since that's the payload
+                result[i] = u8::conditional_select(&result[i], &bytes[8 + i], ok);
             }
             n_found += ok.unwrap_u8();
         }
@@ -62,13 +69,9 @@ impl RistrettoPoint {
         if n_found == 1 { Some(result) } else { None }
     }
 
-    /// Computes the at most 8 positive FieldElements f such that
-    /// `self == RistrettoPoint::elligator_ristretto_flavor(f)`.
-    /// Assumes self is even.
-    ///
-    /// First return value is a bitmask of which elements are valid inverses (lowest bit corresponds
-    /// to the 0-th element).
-    fn elligator_ristretto_flavor_inverse(&self) -> (u8, [FieldElement; 8]) {
+    /// Computes the at most 8 positive FieldElements f such that `self ==
+    /// RistrettoPoint::elligator_ristretto_flavor(f)`. Assumes self is even.
+    fn elligator_ristretto_flavor_inverse(&self) -> [CtOption<FieldElement>; 8] {
         // Elligator2 computes a Point from a FieldElement in two steps: first
         // it computes a (s,t) on the Jacobi quartic and then computes the
         // corresponding even point on the Edwards curve.
@@ -88,26 +91,13 @@ impl RistrettoPoint {
         // at the same time.  The four Jacobi quartic points are two of
         // such pairs.
 
-        let mut mask: u8 = 0;
         let jcs = self.to_jacobi_quartic_ristretto();
-        let mut ret = [FieldElement::ONE; 8];
-
-        for i in 0..4 {
-            let (ok, fe) = jcs[i].e_inv();
-            let mut tmp: u8 = 0;
-            ret[2 * i] = fe;
-            tmp.conditional_assign(&1, ok);
-            mask |= tmp << (2 * i);
-
-            let jc = jcs[i].dual();
-            let (ok, fe) = jc.e_inv();
-            let mut tmp: u8 = 0;
-            ret[2 * i + 1] = fe;
-            tmp.conditional_assign(&1, ok);
-            mask |= tmp << (2 * i + 1);
-        }
 
-        (mask, ret)
+        // Compute the inverse of the every point and its dual
+        let invs = jcs.iter().flat_map(|jc| [jc.e_inv(), jc.dual().e_inv()]);
+        // This cannot panic because jcs is guaranteed to be size 4, and the above iterator expands
+        // it to size 8
+        Array::<_, U8>::from_iter(invs).0
     }
 
     /// Find a point on the Jacobi quartic associated to each of the four
@@ -209,16 +199,12 @@ impl RistrettoPoint {
 
     /// Computes the possible bytestrings that could have produced this point via
     /// [`Self::map_to_curve`].
-    /// First return value is a bitmask of which elements are valid inverses (lowest bit corresponds
-    /// to the 0-th element).
-    pub fn map_to_curve_inverse(&self) -> (u8, [[u8; 32]; 8]) {
-        let mut ret = [[0u8; 32]; 8];
-        let (mask, fes) = self.elligator_ristretto_flavor_inverse();
-
-        for j in 0..8 {
-            ret[j] = fes[j].to_bytes();
-        }
-        (mask, ret)
+    pub fn map_to_curve_inverse(&self) -> [CtOption<[u8; 32]>; 8] {
+        // Compute the inverses
+        let fes = self.elligator_ristretto_flavor_inverse();
+        // Serialize the field elements
+        let it = fes.map(|fe| fe.map(|f| f.to_bytes()));
+        Array::<_, U8>::from_iter(it).0
     }
 }
 
@@ -257,7 +243,7 @@ mod test {
     }
 
     #[test]
-    fn test_lizard_encode() {
+    fn lizard_encode() {
         // Test vectors are of the form (x, y) where y is the compressed encoding of the Ristretto
         // point given by lizard_encode(x).
         // These values come from the testLizard() function in vendor/ristretto.sage
@@ -286,7 +272,7 @@ mod test {
 
     // Tests that lizard_decode of a random point is None
     #[test]
-    fn test_lizard_invalid() {
+    fn lizard_invalid() {
         let mut rng = rand::rng();
         for _ in 0..100 {
             let pt = RistrettoPoint::random(&mut rng);
@@ -298,8 +284,11 @@ mod test {
         }
     }
 
+    // Test that
+    //   elligator_ristretto_flavor ○ elligator_ristretto_flavor_inverse ○ elligator_ristretto_flavor
+    // is the identity
     #[test]
-    fn test_elligator_inv() {
+    fn elligator_inv() {
         let mut rng = rand::rng();
 
         for i in 0..100 {
@@ -318,19 +307,21 @@ mod test {
                 // For the rest, just generate a random field element to test.
                 rng.fill_bytes(&mut fe_bytes);
             }
-            fe_bytes[0] &= 254; // positive
-            fe_bytes[31] &= 127; // < 2^255-19
+            // Make fe positive (even) and less than the modulus
+            fe_bytes[0] &= 254;
+            fe_bytes[31] &= 127;
             let fe = FieldElement::from_bytes(&fe_bytes);
 
             let pt = RistrettoPoint::elligator_ristretto_flavor(&fe);
             for pt2 in xcoset4(&pt) {
-                let (mask, fes) = RistrettoPoint(pt2).elligator_ristretto_flavor_inverse();
+                let fes = RistrettoPoint(pt2).elligator_ristretto_flavor_inverse();
 
                 let mut found = false;
-                for (j, fe_j) in fes.iter().enumerate() {
-                    if mask & (1 << j) != 0 {
-                        assert_eq!(RistrettoPoint::elligator_ristretto_flavor(fe_j), pt);
-                        if *fe_j == fe {
+                for fe_j in fes {
+                    if fe_j.is_some().into() {
+                        let fe_j = fe_j.unwrap();
+                        assert_eq!(RistrettoPoint::elligator_ristretto_flavor(&fe_j), pt);
+                        if fe_j == fe {
                             found = true;
                         }
                     }
@@ -342,7 +333,7 @@ mod test {
 
     // Tests that map_to_curve_inverse ○ map_to_curve is the identity
     #[test]
-    fn test_map_to_curve_inverse() {
+    fn map_to_curve_inverse() {
         let mut rng = rand::rng();
 
         for _ in 0..100 {
@@ -351,7 +342,7 @@ mod test {
 
             // Map to Ristretto and invert it
             let pt = RistrettoPoint::map_to_curve(input);
-            let (bitmask, inverses) = pt.map_to_curve_inverse();
+            let inverses = pt.map_to_curve_inverse();
 
             // map_to_curve masks the bottom bit and top two bits of `input`
             let mut expected_inverse = input;
@@ -360,9 +351,13 @@ mod test {
 
             // Check that one of the valid inverses matches the input
             let mut found = false;
-            for (i, inv) in inverses.into_iter().enumerate() {
-                let is_valid = bitmask & (1 << i) != 0;
-                if is_valid && inv == expected_inverse {
+            for inv in inverses.into_iter() {
+                if inv.is_some().into() && inv.unwrap() == expected_inverse {
+                    // Per the README, the probability of finding two inverses is ~2^-122
+                    if found == true {
+                        panic!("found two inverses for input {:02x?}", input);
+                    }
+
                     found = true;
                 }
             }

From 9ac71b5d79a4dfc2c689e11877b3064e8895bcf7 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Fri, 5 Sep 2025 22:36:09 -0400
Subject: [PATCH 18/20] Add README to lizard folder

---
 curve25519-dalek/src/lizard/README.md | 516 +++++++++++++++++++++++++-
 1 file changed, 515 insertions(+), 1 deletion(-)

diff --git a/curve25519-dalek/src/lizard/README.md b/curve25519-dalek/src/lizard/README.md
index 3b94f9157..24c1746ca 100644
--- a/curve25519-dalek/src/lizard/README.md
+++ b/curve25519-dalek/src/lizard/README.md
@@ -1 +1,515 @@
-Placeholder
+# Computing Preimages of Ristretto255 Points Under the Elligator2 Map
+
+# Introduction
+
+The [Elligator2][elligator] map provides a way to map field elements to elliptic curve points with a nearly uniform distribution. When adapted for [Ristretto255][ristretto], it maps elements of the field $\mathbb{F}_p$ (where $p = 2^{255} - 19$) to points in the Ristretto255 group. This mapping is not injective; each Ristretto255 point can be the image of up to eight distinct field elements.
+
+Some use cases, however, require an injective map into the Ristretto group with an easy-to-compute inverse. This is handy if, e.g., we want to encrypt data with the [ElGamal](https://en.wikipedia.org/wiki/ElGamal_encryption) encryption scheme, where plaintexts are curve points. This injection is precisely what the ["Lizard"][lizard] encoding provides.
+
+In this document, we describe the Lizard encoding/decoding algorithms. Along the way we will describe an efficient method for computing preimages of the Elligator2 map. These algorithms are implemented in `src/lizard/` in this repo.
+
+# Lizard
+The Lizard encoding is an injective map into the Ristretto255 group $\mathcal{R}$ built on the Elligator2 map, $E_2$.
+
+## Encoding: Make Elligator2 an injection
+
+As mentioned above, $E_2$ is not injective and thus cannot work directly as an encoding function.
+
+Instead we create an injective function from $E_2$ by "tagging" each input $x$ so it is distinguished among the other inverses of $E_2(x)$. This requires restricting the map's domain: rather than mapping $\mathbb{F}_{>0} \to \mathcal{R}$, we define a map $\lbrace 0,1 \rbrace^{128} \to \mathcal{R}$ (i.e., half the input length), as follows:
+1. Given bitstring $b \in \lbrace 0,1 \rbrace^{128}$, compute its hash $h = \mathsf{H}(b)$ where $\mathsf{H}$ is a hash function with a 32-byte digest.
+2. Let `s` be the 32-byte value `h[0:8] || b || h[24:32]`.
+3. Clear the least significant bit: `s[0] &= 254`. This is done so that the field element represented by `s` is "positive" (defined in next section).
+4. Clear the most significant two bits: `s[31] &= 63`. This is done so that the field element represented by `s` does not exceed $p$.
+5. Interpret `s` as a field element (little endian), and output $E_2(s)$
+
+## Decoding: Invert the injection
+
+The tag allows us to decode points $P \in \mathcal{R}$ as follows:
+1. Compute the set of preimages $E_2^{-1}(P)$ using the procedure described in the rest of this document. This gives us a set of up to 8 scalars that are candidates to be the correct inverse.
+2. For each candidate scalar, check whether it adheres to the above form. Namely check that it is positive and has byte representation `h[0:8] || b || h[24:32]`, where $h = \mathsf{H}(b)$, and the appopriate bits are cleared.
+3. If such a candidate is found, this is the decoded value.
+
+The above procedures are exactly what the `lizard_encode` and `lizard_decode` functions do in `src/lizard/lizard_ristretto.rs`.
+
+## Security argument
+
+It is not a-priori clear that the function described is an injection. 
+We repeat [Bas Westerbaan's argument](https://github.com/bwesterb/go-ristretto/blob/7c80df9e61a54045aedb2ccad99d1e636a4b4c90/ristretto.go#L224) for why this method produces an injection with overwhelming probability:
+
+> There are some (and with high probability at most 80) inputs to SetLizard() which cannot be decoded. The chance that you hit such an input is around 1 in 2^122.
+
+> In Lizard there are $256 - 128 - 3 = 125$ check bits to pick out the right preimage among at most eight.  Conservatively assuming there are seven other preimages, the chance that one of them passes the check as well is given by:
+
+> $$
+\begin{align}
+1 - (1 - 2^{-125})^7
+&= 7\cdot 2^{-125} + 21\cdot 2^{-250} - ...
+\\\\&\approx 2^{-125 - 2\log(7)}
+\\\\&= 2^{-122.192\ldots}
+\end{align}
+$$
+
+> Presuming a random hash function, the number of "bad" inputs is binomially distributed with $n=2^{128}$ and $p=2^{-122.192}\ldots$ For such large $n$, the Poisson distribution with $\lambda=np=56$ is a  very good approximation. In fact: the cumulative distribution function (CDF) of the Poission distribution is larger than that of the binomial distribution for $k > \lambda$.[^1]  The value of the former on $k=80$ is larger than 0.999 and so with a probability of 99.9%, there are fewer than 80 bad inputs.
+
+# Inverting $E_2$
+ 
+ Methods for efficiently implementing the Elligator2 map $E_2$ have been carefully documented on the [Ristretto255][ristretto] website, and using this the implementation of Lizard *encoding* is straightforward. Lizard *decoding*, on the other hand requires efficiently inverting $E_2$ and we describe that process here.
+
+## Background and Notation
+
+The computations will involve a family of related curves all defined over the prime field of order $p = 2^{255} - 19$ which we denote by $\mathbb{F}$. We will adopt the definition of [Decaf][decaf] and define an element of $\mathbb{F}$ to be *positive* if the low bit of its least positive representitive is set. We denote the set of positive field elements by $\mathbb{F}_{>0}$.
+
+### Fundamental Curves
+
+**The Twisted Edwards Curve $\mathcal{E}$.**
+Much of our work will take place in the Twisted Edwards Curve $\mathcal{E}_{a,d}$ given by
+
+$$
+ax^2+y^2=1+dx^2y^2.
+$$
+
+where $a = -1$ and $d =-121665/121666$. In what follows we will simply denote this curve by $\mathcal{E}$, but we will refer to $a$ and $d$ explicitly when it helps clarify the relationships between curves and dispel any magic behind the apparently magic numbers.
+
+The addition law on $\mathcal{E}$ is given by $(x_1, y_1) + (x_2,y_2) = (x_3,y_3)$ where:
+
+$$
+x_3 = \frac{x_1y_2 + x_2y_1}{1 + d x_1 y_1 x_2 y_2}, \quad y_3 = \frac{ y_1 y_2 - a x_1 x_2}{1 - d x_1 y_1 x_2 y_2}.
+$$
+
+The 4-torsion subgroup will be of particular interest. It is cyclic and consists of the points:
+
+$$
+\mathcal{E}[4] = \lbrace (0,1), (1,0), (0,-1), (-1, 0) \rbrace.
+$$
+
+The size of $\mathcal{E}$ is $|\mathcal{E}| = 8\ell$ where $\ell = 2^{252}+27742317777372353535851937790883648493$ is prime.
+
+**The Jacobi Quartic $\mathcal{J}$.**
+We also consider the Jacobi quartic $\mathcal{J}_{e,A}$ given by
+
+$$
+t^2=es^4+2As^2+1,
+$$
+
+where $e = (-a)^2 = 1$ and $A = -a - 2ad/(a-d) = 243329$.
+We denote this curve by $\mathcal{J}$.
+
+$\mathcal{J}$ is 2-isogenous to $\mathcal{E}$ through the mapping  $\theta:\mathcal{J} \rightarrow \mathcal{E}$:
+
+$$
+\theta(s,t) = \left(\frac{1}{\sqrt{-d-1}}\cdot\frac{2s}{t}, \frac{1-s^2}{1+s^2} \right).
+$$
+
+The addition law on $\mathcal{J}$ is given by $(s_1, t_1)+(s_2, t_2) = (s_3, t_3)$ where
+
+$$
+s_3 = \frac{s_1 t_2 + s_2 t_1}{1 - e s_{1}^2 s_{2}^2},\quad t_3 = \frac{(t_1 t_2 + 2A s_1 s_2)(1 + e s_1^2s_2^2) + 2e s_1  s_2 (s_1^2 + s_2^2)}{(1 - e s_1^2 s_2^2)^2}
+$$
+
+$\mathcal{J}$ has full 2-torsion, that is, $|\mathcal{J}[2]| = 4$ and given a point $(s,t)\in\mathcal{J}$, its coset for the 2-torsion subgroup of $\mathcal{J}$ is given by
+
+$$
+(s,t) + \mathcal{J}[2] = \lbrace (s,t), (-s,-t), (1/as, -t/as^2), (-1/as, t/as^2) \rbrace
+$$
+
+The size of $\mathcal{J}$ is $|\mathcal{J}| = 8\ell$.
+
+**The Montgomery Curve $\mathcal{M}$.**
+Although we will not use it in this document, it is convenient to know that $\mathcal{E}$ and $\mathcal{J}$ are closely related to Curve25519, the Montgomery curve given by the equation
+
+$$
+  y^{2}=x^{3}+486662x^{2}+x = x^3 + 2\frac{a+d}{a - d}x^2 + x
+$$
+
+In particular $\mathcal{M}$ is birationally equivalent to $\mathcal{E}$ and is 2-isogenous to $\mathcal{J}$.
+The size of $\mathcal{M}$ is $|\mathcal{M}| = 8\ell$ 
+
+### Ristretto
+
+The Ristretto255 group is the group $\mathcal{R} = [2]\mathcal{E}/\mathcal{E}[4]$ of prime order $\ell$. Each element of $\mathcal{R}$ can be represented by an even point $(x,y)\in\mathcal{E}_{a,d}$. Writing the four-element coset out in full:
+
+$$
+(x,y) + \mathcal{E}_{a,d}[4] = \left\lbrace (x,y), (\frac{y}{\sqrt{a}}, -x\sqrt{a}), (-x,-y), (\frac{-y}{\sqrt{a}}, x\sqrt{a}) \right\rbrace.
+$$
+
+With our parameter choice of $a = -1$ this simplifies to
+
+$$
+\lbrace (x,y), (-x,-y), (\mathrm{i} y, \mathrm{i} x), (-\mathrm{i} y, -\mathrm{i} x) \rbrace,
+$$
+
+where $\mathrm{i} = \sqrt{-1}$ is the imaginary unit in $\mathbb{F}$.
+
+## The Components of Elligator2
+
+$E_2$ is built from four functions:
+
+1. A map from positive field elements to a Jacobi quartic: $e: \mathbb{F}_{>0} \rightarrow \mathcal{J}$.
+1. The quotient map: $q_\mathcal{J}: \mathcal{J} \rightarrow \mathcal{J}/\mathcal{J}[2]$
+1. A 2-isogeny from the Jacobi quartic to the even points on the Edwards curve: $\theta: \mathcal{J} \rightarrow [2]\mathcal{E}$. Note that we only know that the image of $\theta$ consists of even points because the torsion group $\mathcal{E}[8]$ is cyclic. See discussion below.
+1. The quotient map into the Ristretto group: $q_\mathcal{E}: [2]\mathcal{E}/\mathcal{E}[2] \rightarrow \mathcal{R}$, where $P \mapsto P + \mathcal{E}[4]/\mathcal{E}[2]$.
+
+Now since $\mathsf{ker}(\theta) \leqslant \mathcal{J}[2]$, [the Decaf paper][decaf] shows that the map
+
+$$
+\hat{\theta}: \frac{\mathcal{J}}{\mathcal{J}[2]} \rightarrow \frac{\theta(\mathcal{J})}{\theta(\mathcal{J}[2])} = \frac{[2]\mathcal{E}}{\mathcal{E}[2]}, \qquad P + \mathcal{J}[2] \mapsto \theta(P) + \mathcal{E}[2]
+$$ 
+
+<div id="eq-reduced-theta" align="right">(1)</div>
+
+
+is an isomorphism.
+
+The full map is $E_2 = q_\mathcal{E} \circ \hat{\theta} \circ q_\mathcal{J} \circ e$. Our goal is to compute the preimage set $E_2^{-1}(P)$ for a given Ristretto point $P \in \mathcal{R}$.
+
+Since $q_\mathcal{E}$ is 2-to-1, $\hat{\theta}$ is 1-to-1, and $q_\mathcal{J}$ is 4-to-1, we expect to find $4 \times 2 = 8$ preimages on the Jacobi quartic $\mathcal{J}$. Each of these 8 points on $\mathcal{J}$ will then correspond to a unique field element via the inverse of $e$. This section explains how to find these 8 field elements.
+
+We now describe how to invert each of these components before presenting an efficient algorithm to compute the full inverse of $E_2$.
+
+### The Quotient Maps and their Inverses
+
+Given a point $(x,y) + \mathcal{E}[2] \in \mathcal{E}/\mathcal{E}[2]$, $q_\mathcal{E}$ maps it to its coset in the quotient group in $\mathcal{E}/\mathcal{E}[4]$:
+
+$$
+q_\mathcal{E}(\lbrace (x,y), (-x,-y)\rbrace) = \lbrace (x,y), (-x,-y), (\mathrm{i} y, \mathrm{i} x), (-\mathrm{i} y, -\mathrm{i} x) \rbrace.
+$$
+
+From which we can see it is a 2-to-1 map, since $\lbrace (x,y), (-x,-y) \rbrace$ and $\lbrace (\mathrm{i} y, \mathrm{i} x), (-\mathrm{i} y, -\mathrm{i} x) \rbrace$ map to the same point.
+
+When its domain is restricted to the even points $[2]\mathcal{E}/\mathcal{E}[2]$ we see that $|[2]\mathcal{E}/\mathcal{E}[2]| = 2\ell$ so its range must have size $\ell$ - the size of $\mathcal{R}$. Thus $q_\mathcal{E}$ maps onto $\mathcal{R}$. $q$ is also simple to invert:
+
+$$
+q_\mathcal{E}^{-1}((x,y) + \mathcal{E}[4]) = \left\lbrace \lbrace(x,y), (-x,-y)\rbrace, \lbrace(\mathrm{i} y, \mathrm{i} x), (-\mathrm{i} y, -\mathrm{i} x)\rbrace \right\rbrace.
+$$
+
+Note the extra braces here. $q_\mathcal{E}^{-1}(P)$ contains two elements of $\mathcal{E}/\mathcal{E}[2]$ and each of these elements is a coset of size two.
+
+This formula is valid for all points in $\mathcal{E}/\mathcal{E}[4]$, but for points that are not in $\mathcal{R} = [2]\mathcal{E}/\mathcal{E}[4]$ their preimages will not be in the range of $\hat{\theta}$, described next, and we will not be able to continue the inversion of $E_2$. This is expected since these are not validly encoded points.
+
+The quotient map $q_\mathcal{J}$ is also simple to describe and invert
+
+$$
+q_\mathcal{J}((s,t)) = (s,t) + \mathcal{J}[2] = \lbrace (s,t), (-s,-t), (-1/s, t/s^2), (1/s, -t/s^2) \rbrace.
+$$
+
+So
+
+$$
+q_\mathcal{J}^{-1}((s,t) + \mathcal{J}[2]) = \lbrace (s,t), (-s,-t), (-1/s, t/s^2), (1/s, -t/s^2) \rbrace.
+$$
+
+
+### The Isogeny $\theta$ and its Inverse
+
+The Elligator2 construction for Ristretto255 uses the Jacobi quartic curve $\mathcal{J}$ defined by the equation $t^2 = s^4 + A s^2 + 1$, which is 2-isogenous to our twisted Edwards curve $\mathcal{E}$. The 2-isogeny $\theta: \mathcal{J} \rightarrow \mathcal{E}$ is given by:
+
+$$
+\theta(s,t) = (x,y) = \left( \frac{2s}{t\sqrt{c}}, \frac{1-s^2}{1+s^2} \right)
+$$
+
+where for a twisted Edwards curve $ax^2+y^2=1+dx^2y^2$, the constant $c = a-d$. For Ristretto255, $a=-1$, so $c = -1-d$. Note that $\theta(s, t) = \theta(-s, -t)$, so this map is inherently 2-to-1 from the points on $\mathcal{J}$ to points on $\mathcal{E}$.
+
+#### The Projected Map, $\hat{\theta}$
+
+By a special case of Tate's Isogeny Theorem,[^2] since the isogenous curves $\mathcal{J}$ and $\mathcal{E}$ are both defined over $\mathbb{F}$, they have the same number of rational points, i.e., they have the same number of points in $\mathbb{F}$. Since $\theta$ is 2-to-1, this means that $\theta(\mathcal{J})$ is a subgroup of $\mathcal{E}$ of index 2. Since $\mathcal{E}$ is cyclic, $[2]\mathcal{E}$ is the only index-2 subgroup and we see that $\theta(\mathcal{J}) = [2]\mathcal{E}$.
+
+We also see that  because 2-torsion points must map to 2-torsion points we have $\theta(\mathcal{J}[2]) \leqslant \mathcal{E}[2]$. Also, since $\mathcal{J}[2]$ has 4 points and $|\mathsf{ker(\theta)}| = 2$, at least one point in $\mathcal{J}[2]$ does not map to zero and $|\theta(\mathcal{J}(2))| \geq 2$. But $|\mathcal{E}[2]| = 2$ so it follows that $\theta(\mathcal{J}[2]) = \mathcal{E}[2]$.
+
+Thus we are justified in talking about the projected map $\hat{\theta}$ defined in [Equation 1](#eq-reduced-theta).
+
+#### Inverting $\theta$ (and $\hat{\theta}$)
+
+Given a point $(x,y) \in \mathcal{E}$, we can invert $\theta$ to find its two preimages on $\mathcal{J}$. From the $y$-coordinate, we have:
+
+$$
+y = \frac{1-s^2}{1+s^2} \implies s^2 = \frac{1-y}{1+y}
+$$
+
+This gives two possible values for $s$, which we denote $s_0$ and $-s_0$. From the $x$-coordinate, we can then find the corresponding $t$ value:
+
+$$
+x = \frac{2s}{t\sqrt{c}} \implies t = \frac{2s}{x\sqrt{c}}
+$$
+
+Therefore, for each Edwards point $(x,y)$:
+
+$$
+\theta^{-1}((x,y)) = \lbrace (s_0,t_0), (-s_0, -t_0) \rbrace
+$$
+
+where $s_0 = \sqrt{\frac{1-y}{1+y}}$ and $t_0 = \frac{2s_0}{x\sqrt{c}}$.
+
+Note that when writing an Edwards point in projective coordinates, $X:Y:Z$ this becomes
+
+$$
+s_0 = \sqrt{\frac{Z-Y}{Z+Y}},\qquad t_0 = \frac{2Z}{X\sqrt{c}}\sqrt{\frac{Z-Y}{Z+Y}}
+$$
+
+ These formulas allow us to write the inverse for $\hat{\theta}$ as well. Note that an equivalence class in $[2]\mathcal{E}/\mathcal{E}[2]$ is given by
+
+ $$
+ [X:Y:Z] + \mathcal{E}[2] = \left\lbrace [X:Y:Z], [-X:-Y:Z]\right\rbrace
+ $$
+
+So
+
+ $$
+ \hat{\theta}^{-1}(\lbrace [X:Y:Z], [-X:-Y:Z] \rbrace) = \left\lbrace (s_0, t_0), (-s_0, -t_0), \left(\frac{1}{s_0}, -\frac{t_0}{s_0^2}\right), \left(-\frac{1}{s_0}, \frac{t_0}{s_0^2}\right) \right\rbrace  = (s_0,t_0) + \mathcal{J}[2]
+ $$
+
+ Note that this inverse will exist if and only if $[X:Y:Z] \in [2]\mathcal{E}$.
+
+### $e$ and its Inverse
+<div id="subsec-e-and-inverse"></div>
+
+
+The Elligator2 map from $\mathbb{F}$ to $\mathcal{J}$ is computed as follows. On input $r_0$, let $r = \mathrm{i}r_0^2$ - either $0$ or a non-residue. Then compute
+
+$$
+\begin{align}
+(s,t) &= \left( + \sqrt{ \frac{ (r+1)(d - 1)(d + 1) }{ (d r +1)(r + d) } } , \frac{ (r-1)(d - 1)^2 }{ (d r + 1)(r + d) } - 1 \right)
+\end{align}
+$$
+
+or
+
+$$
+\begin{align}
+(s,t) &= \left( - \sqrt{ \frac{ r(r+1)(d - 1)(d + 1) }{ (d r + 1)(r + d) } } , \frac{ -r(r-1)(-1 + d)^2 }{ (d r + 1)(r + d) } - 1 \right)
+\end{align}
+$$
+
+depending on which square root exists, preferring the second when $r = 0$ and both are square.
+
+#### The case $s > 0$
+
+We want to invert this and will begin with the first case, where $s > 0$. First note that
+
+$$
+\frac{t+1}{s^2} = \frac{d-1}{d+1}\frac{r-1}{r+1}
+$$
+
+<div id="eq-a-over-s2" align="right">(2)</div>
+
+We introduce the value $\alpha = \frac{d+1}{d-1}(t+1)$, which is denoted `a` in the source code. This gives
+
+$$
+\frac{\alpha}{s^2} = \frac{r-1}{r+1}
+$$
+
+from which we solve for $r$:
+
+$$
+r = \frac{s^2 + \alpha}{s^2 - \alpha}
+$$
+
+Thus our inverse will be
+
+$$
+r_0 = \sqrt{-\mathrm{i}\frac{s^2 + \alpha}{s^2 - \alpha}}
+$$
+
+<div id="eq-r0case1" align="right">(3)</div>
+
+#### The case $s\leq 0$
+
+When considering the second case, $s \leq 0$, note that [Equation 2](#eq-a-over-s2) just negates the right hand side
+
+$$
+\frac{t+1}{s^2} = -\frac{d-1}{d+1}\frac{r-1}{r+1}
+$$
+
+leading to the formula
+
+$$
+r =  \frac{s^2 - \alpha}{s^2 + \alpha}
+$$
+
+Thus our inverse will be
+
+$$
+r_0 = \sqrt{-\mathrm{i}\frac{s^2 - \alpha}{s^2 + \alpha}}
+$$
+
+<div id="eq-r0case2" align="right">(4)</div>
+
+# Algorithm for Computing All 8 Preimages of a Point
+
+We can combine these computations to invert $E_2$.
+
+First, when inverting $q_\mathcal{E}$, observe that if a point $P\in\mathcal{R}$ has a representative with projective coordinates $[X:Y:Z]$ then its inverse under $q_\mathcal{E}$ is given by
+
+$$
+q_\mathcal{E}^{-1}(P) = \lbrace \lbrace [X:Y:Z], [-X:-Y:Z] \rbrace, \lbrace [Y:X:\mathrm{i}Z], [-Y:-X:\mathrm{i}Z] \rbrace \rbrace.
+$$
+
+Above we showed how, given a point $X:Y:Z \in [2]\mathcal{E}/\mathcal{E}[2]$ we can compute its preimage 
+
+$$
+\hat{\theta}^{-1}([X:Y:Z] + \mathcal{E}[2]) = \left\lbrace (s_0, t_0), (-s_0, -t_0), \left(\frac{1}{s_0}, -\frac{t_0}{s_0^2}\right), \left(-\frac{1}{s_0}, \frac{t_0}{s_0^2}\right) \right\rbrace
+$$
+
+We can repeat this process to compute preimage of $[Y:X:\mathrm{i}Z] + \mathcal{E}[2]$, the other point in $q_\mathcal{E}^{-1}(P)$:
+
+$$
+\hat{\theta}^{-1}([Y:X:\mathrm{i}Z] + \mathcal{E}[2]) = \left\lbrace (s_2, t_2), (-s_2, -t_2), \left(\frac{1}{s_2}, -\frac{t_2}{s_2^2}\right), \left(-\frac{1}{s_2}, \frac{t_2}{s_2^2}\right) \right\rbrace
+$$
+
+
+This already yields an algorithm to compute inverses of $E_2$ and complete the Lizard decoding, but we can make some significant optimizations.
+
+## Step 1: Find Four Non-Dual Jacobi Points That Map to a Ristretto Point
+<div id="subsec-step1"></div>
+
+Based on the calculations above, once we compute the the four points $(s_i, t_i)$ as follows
+
+$$
+\begin{align}
+    s_0 = \sqrt{\frac{Z-Y}{Z+Y}},&\quad& t_0 = \frac{2Z}{X\sqrt{c}}\sqrt{\frac{Z-Y}{Z+Y}} \\
+    s_1 = \frac{1}{s_0},&\quad& t_1 = -\frac{t_0}{s_0^2} \\
+    s_2 = \sqrt{\frac{\mathrm{i}Z-X}{\mathrm{i}Z+X}},&\quad& t_0 = \frac{2\mathrm{i}Z}{Y\sqrt{c}}\sqrt{\frac{\mathrm{i}Z-X}{\mathrm{i}Z+X}} \\
+    s_3 = \frac{1}{s_2},&\quad& t_3 = -\frac{t_2}{s_2^2} \\
+\end{align}
+$$
+
+we easily get all eight preimages of the Ristretto point $[X:Y:Z] + \mathcal{E}[4]$ by taking their "dual" points, where the *dual* of $(s,t)$ is $(-s,-t)$: 
+$\lbrace \pm(s_i, t_i) \rbrace$.
+
+This computation is performed by the function `to_jacobi_quartic_ristretto` in `src/lizard/lizard_ristretto.rs`. We describe that computation now.
+
+
+### Sharing Computation
+
+We can compute $s_0$ and  $s_1$ more efficiently together by precomputing $1/\sqrt{Z^2 - Y^2}$, then multiplying by $(Z-Y)$ and $(Z+Y)$ respectively. Similarly we can compute $s_2$ and $s_3$ together by precomputing $1/\sqrt{X^2 + Z^2}$ then multiplying by $(\mathrm{i}Z-X)$ and $(\mathrm{i}Z+X)$ respectively.
+
+It turns out that we can perform both of these precomputations while only computing one square root. In particular, compute the following value:
+
+$$
+\gamma = \frac{1}{\sqrt{Y^4X^2(Z^2-Y^2)}}
+$$
+
+Then we can easily compute $s_0$ and $s_1$ with inexpensive field operations:
+
+$$
+s_0 = \gamma Y^2 (Z-Y) X, \quad s_1 = \gamma Y^2 (Z+Y) X
+$$
+
+Now $\gamma Y^2 = 1/\sqrt{X^2(Z^2 - Y^2)}$ and, when we are computing $s_2, s_3$ the analog of this value is $-i/\sqrt{Y^2(X^2 + Z^2)}$. So we can compute
+
+$$
+\begin{align}
+    \sqrt{\frac{-1}{Y^2(X^2 + Z^2)}} &=& \sqrt{\frac{Z^2 - Y^2}{(X^2 + Z^2)(Y^2 - Z^2)}} \\
+         &=& \sqrt{\frac{Z^2 - Y^2}{X^2Y^2 + dX^2Y^2}} \\
+         &=& \frac{1}{\sqrt{d+1}}\sqrt{\frac{Z^2-Y^2}{X^2Y^2}} \\
+         &=&  \frac{1}{\sqrt{d+1}} \gamma (Z^2 - Y^2)
+\end{align}
+$$
+
+Where the equality in the second line comes from using the curve equation for $\mathcal{E}$: $Y^2Z^2 - X^2Z^2 = Z^4 + dX^2Y^2$. The quantity $1/\sqrt{d+1}$ is a constant that can be precomputed, so $\sqrt{-1/(Y^2(X^2 + Z^2))}$ - and hence $(s_2,t_2)$ and $(s_3,t_3)$ - can be computed from $\gamma$ without computing an additional square root.
+
+### Step by Step Through the Rust Implementation
+
+With that background, we can understand the code in `src/lizard/lizard_ristretto.rs`. The following steps can be found in the `RistrettoPoint::to_jacobi_quartic_ristretto()` function:
+
+1. $\gamma \leftarrow 1/\sqrt{Y^4X^2(Z^2 - Y^2)}$
+1. $D \leftarrow \gamma Y^2$
+1. $\mathsf{s\\_over\\_x} = D(Z-Y)$, $\mathsf{sp\\_over\\_xp} = D(Z+Y)$
+1. $s_0 = \mathsf{s\\_over\\_x}\cdot X$, $s_1 = \mathsf{sp\\_over\\_xp}\cdot X$. Note that
+
+$$
+s_0 = D(Z-Y)X = \frac{X(Z-Y)}{\sqrt{X^2(Z^2 - Y^2)}} = \sqrt{\frac{Z-Y}{Z+Y}}
+$$
+
+as required.
+
+5. $\mathsf{tmp} \leftarrow \frac{-2}{\sqrt{-d-1}}Z$
+5. $t_0 \leftarrow \mathsf{tmp}\cdot \mathsf{s\\_over\\_x}$, $t_1 \leftarrow \mathsf{tmp}\cdot \mathsf{sp\\_over\\_xp}$
+5. $(s_0,t_0)$ and $(-s_0, -t_0)$ are preimages of $X:Y:Z$ while $(s_1,t_1)$ and $(-s_1,-t_1)$ are preimages of $-X:-Y:Z$.
+5. $D \leftarrow (Z^2 - Y^2)\gamma \cdot -1/\sqrt{-d-1}$ (This equals $\frac{1}{\sqrt{-d-1}}\sqrt{-\frac{Z^2 - Y^2}{Y^4X^2}}$)
+5. $\mathsf{s\\_over\\_y} \gets D(\mathrm{i}Z - X)$
+5. $\mathsf{sp\\_over\\_yp} \leftarrow D (\mathrm{i}Z + X)$
+5. $s_2 \leftarrow \mathsf{s\\_over\\_y} \cdot Y$, $s_3 \leftarrow \mathsf{sp\\_over\\_yp} \cdot Y$
+5. $\mathsf{tmp} \leftarrow \frac{-2}{\sqrt{-d-1}}\mathrm{i}Z$
+5. $t_2 \leftarrow \mathsf{tmp}\cdot \mathsf{s\\_over\\_y}$, $t_3 \leftarrow \mathsf{tmp}\cdot \mathsf{sp\\_over\\_yp}$
+5. $(s_2,t_2)$ and $(-s_2, -t_2)$ are preimages of $Y:X:\mathrm{i}Z$ while $(s_3,t_3)$ and $(-s_3,-t_3)$ are preimages of $-Y:-X:\mathrm{i}Z$.
+
+Altogether this gives us eight candidate Jacobi Quartic points:
+
+$$
+(s_0, t_0),(s_1, t_1),(s_2, t_2),(s_3, t_3),
+(-s_0, -t_0),(-s_1, -t_1),(-s_2, -t_2),(-s_3, -t_3)
+$$
+
+Note the `to_jacobi_quartic_ristretto()` function only returns the first 4 points. The remaining ones are derived in `RistrettoPoint::elligator_ristretto_flavor_inverse()` in the same file.
+
+The conditional assignments in `to_jacobi_quartic_ristretto()` are explained below.
+
+### Handling Special Cases
+
+The formulas above assume non-zero denominators. The special cases occur when the input point is a 2-torsion or 4-torsion point.
+
+If the input point is the identity element, its representatives are $(0,1)$, $(0,-1)$, $(\mathrm{i}, 0)$, and $(-\mathrm{i}, 0)$.
+
+1. For $(0,1)$, we have $y=1$, so $s^2 = 0 \implies s=0$. This gives the Jacobi point $(0,1)$.
+1. For $(0,-1)$, we have $y=-1$, so $s^2$ is undefined (division by zero). However, this 1. corresponds to $s^2 \to \infty$. In the code we simply repeat the point $(0,1)$ in this 1. This repetition is harmless since this point has no preimage under the Elligator2 map from $\mathbb{F}$ to $\mathcal{J}$ and cannot correspond to a Lizard encoded bitstring.
+1. For $(\mathrm{i}, 0)$, we have $y=0$, so $s^2=1 \implies s=\pm 1$. The preimages are $(\pm 1, \frac{2\mathrm{i}}{\sqrt{-d-1}})$.
+1. For $(-\mathrm{i}, 0)$, we have $y=0$, so $s^2=1 \implies s=\pm 1$. The preimages are $(\pm 1, \mp\frac{2\mathrm{i}}{\sqrt{-d-1}}))$.
+
+
+
+## Step 3: Compute the 8 Field Element Preimages
+
+The computation of $e^{-1}$ [described above](#subsec-e-and-inverse) required breaking the computation into two cases. To facilitate performing these computations efficiently in constant time, we perform both computations at once. This code can be found in `src/lizard/jacobi_quartic.rs` in `JacbobiPoint::e_inv()`.
+
+We begin by computing the value
+
+$$
+y = \frac{1}{\sqrt{\mathrm{i} (s^4 - \alpha^2)}}
+$$
+
+Note that
+1. $\frac{1}{\mathrm{i} (s^4 - \alpha^2)}$ is a quadratic residue if and only if $(s^4 - \alpha^2)$ is not a quadratic residue.
+1. $s^4 - \alpha^2$ is not a quadratic residue if and only if exactly one of $s^2 - \alpha$ and $s^2 + \alpha$ is a quadratic residue.
+1. exactly one of $s^2 - \alpha$ and $s^2 + \alpha$ is a quadratic residue if and only if both $\frac{s^2 - \alpha}{s^2 + \alpha}$ and $\frac{s^2 + \alpha}{s^2 - \alpha}$ are non-residues.
+1. Both $\frac{s^2 - \alpha}{s^2 + \alpha}$ and $\frac{s^2 + \alpha}{s^2 - \alpha}$ are non-residues if and only if $\sqrt{-\mathrm{i}\frac{s^2 + \alpha}{s^2 - \alpha}}$ and $\sqrt{-\mathrm{i}\frac{s^2 - \alpha}{s^2 + \alpha}}$ exist.
+
+Hence we see that $(s,t)$ has a preimage if and only if $y$, as defined above, exists.
+
+With $y$ computed we can now compute the final preimage. Let $\mathsf{sgn}(s)$ denote the sign function. Then we can consolidate Equations [3](#eq-r0case1) and [4](#eq-r0case2) into
+
+$$
+r_0 =  (\alpha + \mathsf{sgn}(s)s^2) y
+$$
+
+### Relation to the Formulas in Decaf
+
+These computations can also be seen in [Decaf][decaf], which also defines the map $e: \mathbb{F}_{>0} \rightarrow \mathcal{J}$, but we need to clarify an important difference in notation. In that work, they use variables $a$ and $d$ as parameters for a different, but isogenous, twisted Edwards curve. If we let $\hat{a}$ and $\hat{d}$ denote the Edwards curve parameters from Decaf we have the relations
+
+$$
+\hat{a} = -a, \quad \hat{d} = \frac{ad}{a-d}
+$$
+
+Note that in this notation we have
+
+$$
+2\hat{d} - \hat{a} =
+\frac{d-1}{d+1}
+$$
+
+hence the formula given in for the inverse of $e$ can be rewritten
+
+$$
+nr_0^2 = r = \frac{(2\hat{d} - \hat{a})s^2 + c(t+1)}{(2\hat{d} - \hat{a})s^2 - c(t+1)} = \frac{cs^2 + (t+1)\frac{d+1}{d-1}}{cs^2 - (t+1)\frac{d+1}{d-1}}
+$$
+
+where $c = \mathsf{sgn}(s)$. There is another detail we need to take care of. The map in Decaf is defined parameterized by a non-square $n\in\mathbb{F}_{p}$. We take $n = \sqrt{-1}$ as the non-residue parameter.
+
+You can see this choice of $n$, in the implementation of `RistrettoPoint::elligator_ristretto_flavor()`.
+
+[^1]: Anderson, T. W. and Samuels, S. M., Some Inequalities Among Binomial and Poisson Probabilities _Proc. Fifth Berkeley Symp. on Math. Statist. and Prob._, Vol. 1. Univ. of Calif. Press (1967).
+[^2]: Tate, J. Endomorphisms of abelian varieties over finite fields. _Invent Math_ 2, 134–144 (1966). https://doi.org/10.1007/BF01404549 [[pdf](https://pazuki.perso.math.cnrs.fr/index_fichiers/Tate66.pdf)]
+
+[elligator]: https://elligator.org/
+[ristretto]: https://ristretto.group/
+[lizard]: https://arxiv.org/abs/1911.02674
+[decaf]: https://eprint.iacr.org/2015/673

From 2d1128db0916950e2b76b792c21d496f05f34194 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Fri, 5 Sep 2025 22:38:14 -0400
Subject: [PATCH 19/20] Make clippy happy

---
 curve25519-dalek/src/lizard/jacobi_quartic.rs | 6 +-----
 1 file changed, 1 insertion(+), 5 deletions(-)

diff --git a/curve25519-dalek/src/lizard/jacobi_quartic.rs b/curve25519-dalek/src/lizard/jacobi_quartic.rs
index 897f57074..1c6c33dbe 100644
--- a/curve25519-dalek/src/lizard/jacobi_quartic.rs
+++ b/curve25519-dalek/src/lizard/jacobi_quartic.rs
@@ -2,11 +2,7 @@
 
 #![allow(non_snake_case)]
 
-use subtle::Choice;
-use subtle::ConditionallyNegatable;
-use subtle::ConditionallySelectable;
-use subtle::ConstantTimeEq;
-use subtle::CtOption;
+use subtle::{ConditionallyNegatable, ConditionallySelectable, ConstantTimeEq, CtOption};
 
 use super::lizard_constants;
 use crate::constants;

From c98df3c60533dd805982b01fde99eaf5b14157b8 Mon Sep 17 00:00:00 2001
From: Michael Rosenberg <mrosenberg@cloudflare.com>
Date: Fri, 5 Sep 2025 22:41:49 -0400
Subject: [PATCH 20/20] Fix typos

---
 curve25519-dalek/src/lizard/README.md | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/curve25519-dalek/src/lizard/README.md b/curve25519-dalek/src/lizard/README.md
index 24c1746ca..6dbe9dd0b 100644
--- a/curve25519-dalek/src/lizard/README.md
+++ b/curve25519-dalek/src/lizard/README.md
@@ -26,7 +26,7 @@ Instead we create an injective function from $E_2$ by "tagging" each input $x$ s
 
 The tag allows us to decode points $P \in \mathcal{R}$ as follows:
 1. Compute the set of preimages $E_2^{-1}(P)$ using the procedure described in the rest of this document. This gives us a set of up to 8 scalars that are candidates to be the correct inverse.
-2. For each candidate scalar, check whether it adheres to the above form. Namely check that it is positive and has byte representation `h[0:8] || b || h[24:32]`, where $h = \mathsf{H}(b)$, and the appopriate bits are cleared.
+2. For each candidate scalar, check whether it adheres to the above form. Namely check that it is positive and has byte representation `h[0:8] || b || h[24:32]`, where $h = \mathsf{H}(b)$, and the appropriate bits are cleared.
 3. If such a candidate is found, this is the decoded value.
 
 The above procedures are exactly what the `lizard_encode` and `lizard_decode` functions do in `src/lizard/lizard_ristretto.rs`.
@@ -57,7 +57,7 @@ $$
 
 ## Background and Notation
 
-The computations will involve a family of related curves all defined over the prime field of order $p = 2^{255} - 19$ which we denote by $\mathbb{F}$. We will adopt the definition of [Decaf][decaf] and define an element of $\mathbb{F}$ to be *positive* if the low bit of its least positive representitive is set. We denote the set of positive field elements by $\mathbb{F}_{>0}$.
+The computations will involve a family of related curves all defined over the prime field of order $p = 2^{255} - 19$ which we denote by $\mathbb{F}$. We will adopt the definition of [Decaf][decaf] and define an element of $\mathbb{F}$ to be *positive* if the low bit of its least positive representative is set. We denote the set of positive field elements by $\mathbb{F}_{>0}$.
 
 ### Fundamental Curves