snarkvm_curves/templates/short_weierstrass_jacobian/
projective.rs

1// Copyright 2024-2025 Aleo Network Foundation
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, serialize::*};
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.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    #[must_use]
295    fn double(&self) -> Self {
296        let mut tmp = *self;
297        tmp.double_in_place();
298        tmp
299    }
300
301    #[inline]
302    fn double_in_place(&mut self) {
303        if self.is_zero() {
304            return;
305        }
306
307        if P::WEIERSTRASS_A.is_zero() {
308            // A = X1^2
309            let mut a = self.x.square();
310
311            // B = Y1^2
312            let b = self.y.square();
313
314            // C = B^2
315            let mut c = b.square();
316
317            // D = 2*((X1+B)2-A-C)
318            let d = ((self.x + b).square() - a - c).double();
319
320            // E = 3*A
321            let old_a = a;
322            a.double_in_place();
323            let e = old_a + a;
324
325            // F = E^2
326            let f = e.square();
327
328            // Z3 = 2*Y1*Z1
329            self.z *= &self.y;
330            self.z.double_in_place();
331
332            // X3 = F-2*D
333            self.x = f - d.double();
334
335            // Y3 = E*(D-X3)-8*C
336            c.double_in_place();
337            c.double_in_place();
338            c.double_in_place();
339            self.y = (d - self.x) * e - c;
340        } else {
341            // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-dbl-2009-l
342            // XX = X1^2
343            let xx = self.x.square();
344
345            // YY = Y1^2
346            let yy = self.y.square();
347
348            // YYYY = YY^2
349            let mut yyyy = yy.square();
350
351            // ZZ = Z1^2
352            let zz = self.z.square();
353
354            // S = 2*((X1+YY)^2-XX-YYYY)
355            let s = ((self.x + yy).square() - xx - yyyy).double();
356
357            // M = 3*XX+a*ZZ^2
358            let m = xx.double() + xx + P::mul_by_a(&zz.square());
359
360            // T = M^2-2*S
361            let t = m.square() - s.double();
362
363            // X3 = T
364            self.x = t;
365            // Y3 = M*(S-T)-8*YYYY
366            let old_y = self.y;
367            yyyy.double_in_place();
368            yyyy.double_in_place();
369            yyyy.double_in_place();
370            self.y = m * (s - t) - yyyy;
371            // Z3 = (Y1+Z1)^2-YY-ZZ
372            self.z = (old_y + self.z).square() - yy - zz;
373        }
374    }
375
376    #[inline]
377    fn to_affine(&self) -> Affine<P> {
378        (*self).into()
379    }
380}
381
382impl<P: Parameters> Neg for Projective<P> {
383    type Output = Self;
384
385    #[inline]
386    fn neg(self) -> Self {
387        if !self.is_zero() { Self::new(self.x, -self.y, self.z) } else { self }
388    }
389}
390
391impl_add_sub_from_field_ref!(Projective, Parameters);
392
393impl<'a, P: Parameters> Add<&'a Self> for Projective<P> {
394    type Output = Self;
395
396    #[inline]
397    fn add(self, other: &'a Self) -> Self {
398        let mut copy = self;
399        copy += other;
400        copy
401    }
402}
403
404impl<'a, P: Parameters> AddAssign<&'a Self> for Projective<P> {
405    #[allow(clippy::many_single_char_names)]
406    #[allow(clippy::suspicious_op_assign_impl)]
407    fn add_assign(&mut self, other: &'a Self) {
408        if self.is_zero() {
409            *self = *other;
410            return;
411        }
412
413        if other.is_zero() {
414            return;
415        }
416
417        // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#addition-add-2007-bl
418        // Works for all curves.
419
420        // Z1Z1 = Z1^2
421        let z1z1 = self.z.square();
422
423        // Z2Z2 = Z2^2
424        let z2z2 = other.z.square();
425
426        // U1 = X1*Z2Z2
427        let u1 = self.x * z2z2;
428
429        // U2 = X2*Z1Z1
430        let u2 = other.x * z1z1;
431
432        // S1 = Y1*Z2*Z2Z2
433        let s1 = self.y * other.z * z2z2;
434
435        // S2 = Y2*Z1*Z1Z1
436        let s2 = other.y * self.z * z1z1;
437
438        if u1 == u2 && s1 == s2 {
439            // The two points are equal, so we double.
440            self.double_in_place();
441        } else {
442            // If we're adding -a and a together, self.z becomes zero as H becomes zero.
443
444            // H = U2-U1
445            let h = u2 - u1;
446
447            // I = (2*H)^2
448            let i = (h.double()).square();
449
450            // J = H*I
451            let j = h * i;
452
453            // r = 2*(S2-S1)
454            let r = (s2 - s1).double();
455
456            // V = U1*I
457            let v = u1 * i;
458
459            // X3 = r^2 - J - 2*V
460            self.x = r.square() - j - (v.double());
461
462            // Y3 = r*(V - X3) - 2*S1*J
463            self.y = P::BaseField::sum_of_products([r, -s1.double()].iter(), [(v - self.x), j].iter());
464
465            // Z3 = ((Z1+Z2)^2 - Z1Z1 - Z2Z2)*H
466            self.z = ((self.z + other.z).square() - z1z1 - z2z2) * h;
467        }
468    }
469}
470
471impl<'a, P: Parameters> Sub<&'a Self> for Projective<P> {
472    type Output = Self;
473
474    #[inline]
475    fn sub(self, other: &'a Self) -> Self {
476        let mut copy = self;
477        copy -= other;
478        copy
479    }
480}
481
482impl<'a, P: Parameters> SubAssign<&'a Self> for Projective<P> {
483    fn sub_assign(&mut self, other: &'a Self) {
484        *self += &(-(*other));
485    }
486}
487
488impl<P: Parameters> Mul<P::ScalarField> for Projective<P> {
489    type Output = Self;
490
491    /// Performs scalar multiplication of this element.
492    #[allow(clippy::suspicious_arithmetic_impl)]
493    #[inline]
494    fn mul(self, other: P::ScalarField) -> Self {
495        P::mul_projective(self, other)
496    }
497}
498
499impl<P: Parameters> MulAssign<P::ScalarField> for Projective<P> {
500    /// Performs scalar multiplication of this element.
501    fn mul_assign(&mut self, other: P::ScalarField) {
502        *self = *self * other
503    }
504}
505
506/// The affine point X, Y is represented in the Jacobian coordinates with Z = 1.
507impl<P: Parameters> From<Affine<P>> for Projective<P> {
508    #[inline]
509    fn from(p: Affine<P>) -> Projective<P> {
510        if p.is_zero() { Self::zero() } else { Self::new(p.x, p.y, P::BaseField::one()) }
511    }
512}