snarkvm_curves/templates/short_weierstrass_jacobian/
projective.rs

1// Copyright (c) 2019-2025 Provable Inc.
2// This file is part of the snarkVM library.
3
4// Licensed under the Apache License, Version 2.0 (the "License");
5// you may not use this file except in compliance with the License.
6// You may obtain a copy of the License at:
7
8// http://www.apache.org/licenses/LICENSE-2.0
9
10// Unless required by applicable law or agreed to in writing, software
11// distributed under the License is distributed on an "AS IS" BASIS,
12// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13// See the License for the specific language governing permissions and
14// limitations under the License.
15
16use crate::{
17    templates::short_weierstrass_jacobian::Affine,
18    traits::{AffineCurve, ProjectiveCurve, ShortWeierstrassParameters as Parameters},
19};
20use snarkvm_fields::{Field, One, Zero, impl_add_sub_from_field_ref};
21use snarkvm_utilities::{FromBytes, ToBytes, cfg_iter_mut, rand::Uniform};
22
23use core::{
24    fmt::{Display, Formatter, Result as FmtResult},
25    hash::{Hash, Hasher},
26    ops::{Add, AddAssign, Mul, MulAssign, Neg, Sub, SubAssign},
27};
28use rand::{
29    Rng,
30    distributions::{Distribution, Standard},
31};
32#[cfg(not(feature = "serial"))]
33use rayon::prelude::*;
34use std::io::{Read, Result as IoResult, Write};
35
36#[derive(Copy, Clone, Debug)]
37pub struct Projective<P: Parameters> {
38    pub x: P::BaseField,
39    pub y: P::BaseField,
40    pub z: P::BaseField,
41}
42
43impl<P: Parameters> Projective<P> {
44    pub const fn new(x: P::BaseField, y: P::BaseField, z: P::BaseField) -> Self {
45        Self { x, y, z }
46    }
47}
48
49impl<P: Parameters> Zero for Projective<P> {
50    // The point at infinity is always represented by Z = 0.
51    #[inline]
52    fn zero() -> Self {
53        Self::new(P::BaseField::zero(), P::BaseField::one(), P::BaseField::zero())
54    }
55
56    // The point at infinity is always represented by Z = 0.
57    #[inline]
58    fn is_zero(&self) -> bool {
59        self.z.is_zero()
60    }
61}
62
63impl<P: Parameters> Default for Projective<P> {
64    #[inline]
65    fn default() -> Self {
66        Self::zero()
67    }
68}
69
70impl<P: Parameters> Display for Projective<P> {
71    fn fmt(&self, f: &mut Formatter<'_>) -> FmtResult {
72        write!(f, "{}", self.to_affine())
73    }
74}
75
76impl<P: Parameters> Hash for Projective<P> {
77    fn hash<H: Hasher>(&self, state: &mut H) {
78        self.to_affine().hash(state);
79    }
80}
81
82impl<P: Parameters> Eq for Projective<P> {}
83
84impl<P: Parameters> PartialEq for Projective<P> {
85    fn eq(&self, other: &Self) -> bool {
86        if self.is_zero() {
87            return other.is_zero();
88        }
89
90        if other.is_zero() {
91            return false;
92        }
93
94        // The points (X, Y, Z) and (X', Y', Z')
95        // are equal when (X * Z^2) = (X' * Z'^2)
96        // and (Y * Z^3) = (Y' * Z'^3).
97        let z1 = self.z.square();
98        let z2 = other.z.square();
99
100        !(self.x * z2 != other.x * z1 || self.y * (z2 * other.z) != other.y * (z1 * self.z))
101    }
102}
103
104impl<P: Parameters> PartialEq<Affine<P>> for Projective<P> {
105    fn eq(&self, other: &Affine<P>) -> bool {
106        if self.is_zero() {
107            return other.is_zero();
108        }
109
110        if other.is_zero() {
111            return false;
112        }
113
114        // The points (X, Y, Z) and (X', Y', Z')
115        // are equal when (X * Z^2) = (X' * Z'^2)
116        // and (Y * Z^3) = (Y' * Z'^3).
117        let z1 = self.z.square();
118        (self.x == other.x * z1) & (self.y == other.y * z1 * self.z)
119    }
120}
121
122impl<P: Parameters> Distribution<Projective<P>> for Standard {
123    #[inline]
124    fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> Projective<P> {
125        loop {
126            let x = P::BaseField::rand(rng);
127            let greatest = rng.r#gen();
128
129            if let Some(p) = Affine::from_x_coordinate(x, greatest) {
130                return p.mul_by_cofactor_to_projective();
131            }
132        }
133    }
134}
135
136impl<P: Parameters> ToBytes for Projective<P> {
137    #[inline]
138    fn write_le<W: Write>(&self, mut writer: W) -> IoResult<()> {
139        self.x.write_le(&mut writer)?;
140        self.y.write_le(&mut writer)?;
141        self.z.write_le(writer)
142    }
143}
144
145impl<P: Parameters> FromBytes for Projective<P> {
146    #[inline]
147    fn read_le<R: Read>(mut reader: R) -> IoResult<Self> {
148        let x = P::BaseField::read_le(&mut reader)?;
149        let y = P::BaseField::read_le(&mut reader)?;
150        let z = P::BaseField::read_le(reader)?;
151        Ok(Self::new(x, y, z))
152    }
153}
154
155impl<P: Parameters> ProjectiveCurve for Projective<P> {
156    type Affine = Affine<P>;
157    type BaseField = P::BaseField;
158    type ScalarField = P::ScalarField;
159
160    #[inline]
161    fn prime_subgroup_generator() -> Self {
162        Affine::prime_subgroup_generator().into()
163    }
164
165    #[inline]
166    fn is_normalized(&self) -> bool {
167        self.is_zero() || self.z.is_one()
168    }
169
170    /// TODO (howardwu): This method can likely be sped up.
171    #[inline]
172    fn batch_normalization(v: &mut [Self]) {
173        // Montgomery's Trick and Fast Implementation of Masked AES
174        // Genelle, Prouff and Quisquater
175        // Section 3.2
176
177        // First pass: compute [a, ab, abc, ...]
178        let mut prod = Vec::with_capacity(v.len());
179        let mut tmp = P::BaseField::one();
180        for g in v
181            .iter_mut()
182            // Ignore normalized elements
183            .filter(|g| !g.is_normalized())
184        {
185            tmp.mul_assign(&g.z);
186            prod.push(tmp);
187        }
188
189        // Invert `tmp`.
190        tmp = tmp.inverse().unwrap(); // Guaranteed to be nonzero.
191
192        // Second pass: iterate backwards to compute inverses
193        for (g, s) in v
194            .iter_mut()
195            // Backwards
196            .rev()
197            // Ignore normalized elements
198            .filter(|g| !g.is_normalized())
199            // Backwards, skip last element, fill in one for last term.
200            .zip(
201                prod.into_iter()
202                    .rev()
203                    .skip(1)
204                    .chain(Some(P::BaseField::one())),
205            )
206        {
207            // tmp := tmp * g.z; g.z := tmp * s = 1/z
208            let newtmp = tmp * g.z;
209            g.z = tmp * s;
210            tmp = newtmp;
211        }
212        cfg_iter_mut!(v).filter(|g| !g.is_normalized()).for_each(|g| {
213            // Perform affine transformations
214            let z2 = g.z.square(); // 1/z
215            g.x *= &z2; // x/z^2
216            g.y *= &(z2 * g.z); // y/z^3
217            g.z = P::BaseField::one(); // z = 1
218        });
219    }
220
221    #[allow(clippy::many_single_char_names)]
222    fn add_assign_mixed(&mut self, other: &Self::Affine) {
223        if other.is_zero() {
224            return;
225        }
226
227        if self.is_zero() {
228            self.x = other.x;
229            self.y = other.y;
230            self.z = P::BaseField::one();
231            return;
232        }
233
234        // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-madd-2007-bl
235        // Works for all curves.
236
237        // Z1Z1 = Z1^2
238        let z1z1 = self.z.square();
239
240        // U2 = X2*Z1Z1
241        let u2 = other.x * z1z1;
242
243        // S2 = Y2*Z1*Z1Z1
244        let s2 = (other.y * self.z) * z1z1;
245
246        if self.x == u2 && self.y == s2 {
247            // The two points are equal, so we double.
248            self.double_in_place();
249        } else {
250            // If we're adding -a and a together, self.z becomes zero as H becomes zero.
251
252            // H = U2-X1
253            let mut h = u2;
254            h -= &self.x;
255
256            // HH = H^2
257            let hh = h.square();
258
259            // I = 4*HH
260            let mut i = hh;
261            i.double_in_place();
262            i.double_in_place();
263
264            // J = H*I
265            let mut j = h;
266            j *= &i;
267
268            // r = 2*(S2-Y1)
269            let mut r = s2;
270            r -= &self.y;
271            r.double_in_place();
272
273            // V = X1*I
274            let mut v = self.x;
275            v *= &i;
276
277            // X3 = r^2 - J - 2*V
278            self.x = r.square();
279            self.x -= &j;
280            self.x -= &v.double();
281
282            // Y3 = r*(V-X3)-2*Y1*J
283            self.y = P::BaseField::sum_of_products([r, -self.y.double()].iter(), [(v - self.x), j].iter());
284
285            // Z3 = (Z1+H)^2-Z1Z1-HH
286            self.z += &h;
287            self.z.square_in_place();
288            self.z -= &z1z1;
289            self.z -= &hh;
290        }
291    }
292
293    #[inline]
294    fn double(&self) -> Self {
295        let mut tmp = *self;
296        tmp.double_in_place();
297        tmp
298    }
299
300    #[inline]
301    fn double_in_place(&mut self) {
302        if self.is_zero() {
303            return;
304        }
305
306        if P::WEIERSTRASS_A.is_zero() {
307            // A = X1^2
308            let mut a = self.x.square();
309
310            // B = Y1^2
311            let b = self.y.square();
312
313            // C = B^2
314            let mut c = b.square();
315
316            // D = 2*((X1+B)2-A-C)
317            let d = ((self.x + b).square() - a - c).double();
318
319            // E = 3*A
320            let old_a = a;
321            a.double_in_place();
322            let e = old_a + a;
323
324            // F = E^2
325            let f = e.square();
326
327            // Z3 = 2*Y1*Z1
328            self.z *= &self.y;
329            self.z.double_in_place();
330
331            // X3 = F-2*D
332            self.x = f - d.double();
333
334            // Y3 = E*(D-X3)-8*C
335            c.double_in_place();
336            c.double_in_place();
337            c.double_in_place();
338            self.y = (d - self.x) * e - c;
339        } else {
340            // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-dbl-2009-l
341            // XX = X1^2
342            let xx = self.x.square();
343
344            // YY = Y1^2
345            let yy = self.y.square();
346
347            // YYYY = YY^2
348            let mut yyyy = yy.square();
349
350            // ZZ = Z1^2
351            let zz = self.z.square();
352
353            // S = 2*((X1+YY)^2-XX-YYYY)
354            let s = ((self.x + yy).square() - xx - yyyy).double();
355
356            // M = 3*XX+a*ZZ^2
357            let m = xx.double() + xx + P::mul_by_a(&zz.square());
358
359            // T = M^2-2*S
360            let t = m.square() - s.double();
361
362            // X3 = T
363            self.x = t;
364            // Y3 = M*(S-T)-8*YYYY
365            let old_y = self.y;
366            yyyy.double_in_place();
367            yyyy.double_in_place();
368            yyyy.double_in_place();
369            self.y = m * (s - t) - yyyy;
370            // Z3 = (Y1+Z1)^2-YY-ZZ
371            self.z = (old_y + self.z).square() - yy - zz;
372        }
373    }
374
375    #[inline]
376    fn to_affine(&self) -> Affine<P> {
377        (*self).into()
378    }
379}
380
381impl<P: Parameters> Neg for Projective<P> {
382    type Output = Self;
383
384    #[inline]
385    fn neg(self) -> Self {
386        if !self.is_zero() { Self::new(self.x, -self.y, self.z) } else { self }
387    }
388}
389
390impl_add_sub_from_field_ref!(Projective, Parameters);
391
392impl<'a, P: Parameters> Add<&'a Self> for Projective<P> {
393    type Output = Self;
394
395    #[inline]
396    fn add(self, other: &'a Self) -> Self {
397        let mut copy = self;
398        copy += other;
399        copy
400    }
401}
402
403impl<'a, P: Parameters> AddAssign<&'a Self> for Projective<P> {
404    #[allow(clippy::many_single_char_names)]
405    #[allow(clippy::suspicious_op_assign_impl)]
406    fn add_assign(&mut self, other: &'a Self) {
407        if self.is_zero() {
408            *self = *other;
409            return;
410        }
411
412        if other.is_zero() {
413            return;
414        }
415
416        // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-add-2007-bl
417        // Works for all curves.
418
419        // Z1Z1 = Z1^2
420        let z1z1 = self.z.square();
421
422        // Z2Z2 = Z2^2
423        let z2z2 = other.z.square();
424
425        // U1 = X1*Z2Z2
426        let u1 = self.x * z2z2;
427
428        // U2 = X2*Z1Z1
429        let u2 = other.x * z1z1;
430
431        // S1 = Y1*Z2*Z2Z2
432        let s1 = self.y * other.z * z2z2;
433
434        // S2 = Y2*Z1*Z1Z1
435        let s2 = other.y * self.z * z1z1;
436
437        if u1 == u2 && s1 == s2 {
438            // The two points are equal, so we double.
439            self.double_in_place();
440        } else {
441            // If we're adding -a and a together, self.z becomes zero as H becomes zero.
442
443            // H = U2-U1
444            let h = u2 - u1;
445
446            // I = (2*H)^2
447            let i = (h.double()).square();
448
449            // J = H*I
450            let j = h * i;
451
452            // r = 2*(S2-S1)
453            let r = (s2 - s1).double();
454
455            // V = U1*I
456            let v = u1 * i;
457
458            // X3 = r^2 - J - 2*V
459            self.x = r.square() - j - (v.double());
460
461            // Y3 = r*(V - X3) - 2*S1*J
462            self.y = P::BaseField::sum_of_products([r, -s1.double()].iter(), [(v - self.x), j].iter());
463
464            // Z3 = ((Z1+Z2)^2 - Z1Z1 - Z2Z2)*H
465            self.z = ((self.z + other.z).square() - z1z1 - z2z2) * h;
466        }
467    }
468}
469
470impl<'a, P: Parameters> Sub<&'a Self> for Projective<P> {
471    type Output = Self;
472
473    #[inline]
474    fn sub(self, other: &'a Self) -> Self {
475        let mut copy = self;
476        copy -= other;
477        copy
478    }
479}
480
481impl<'a, P: Parameters> SubAssign<&'a Self> for Projective<P> {
482    fn sub_assign(&mut self, other: &'a Self) {
483        *self += &(-(*other));
484    }
485}
486
487impl<P: Parameters> Mul<P::ScalarField> for Projective<P> {
488    type Output = Self;
489
490    /// Performs scalar multiplication of this element.
491    #[allow(clippy::suspicious_arithmetic_impl)]
492    #[inline]
493    fn mul(self, other: P::ScalarField) -> Self {
494        P::mul_projective(self, other)
495    }
496}
497
498impl<P: Parameters> MulAssign<P::ScalarField> for Projective<P> {
499    /// Performs scalar multiplication of this element.
500    fn mul_assign(&mut self, other: P::ScalarField) {
501        *self = *self * other
502    }
503}
504
505/// The affine point X, Y is represented in the Jacobian coordinates with Z = 1.
506impl<P: Parameters> From<Affine<P>> for Projective<P> {
507    #[inline]
508    fn from(p: Affine<P>) -> Projective<P> {
509        if p.is_zero() { Self::zero() } else { Self::new(p.x, p.y, P::BaseField::one()) }
510    }
511}