pub fn double(&self) -> Self {
        if self.is_neutral_element() {
            return self.clone();
        }
        let [x1, y1, z1] = self.coordinates();
        //http://www.hyperelliptic.org/EFD/g1p/data/shortw/jacobian-0/doubling/dbl-2009-l

        if E::a() == FieldElement::zero() {
            let a = x1.square(); // A = x1^2
            let b = y1.square(); // B = y1^2
            let c = b.square(); // C = B^2
            let x1_plus_b = x1 + &b; // (X1 + B)
            let x1_plus_b_square = x1_plus_b.square(); // (X1 + B)^2
            let d = (&x1_plus_b_square - &a - &c).double(); // D = 2 * ((X1 + B)^2 - A - C)
            let e = &a.double() + &a; // E = 3 * A
            let f = e.square(); // F = E^2
            let x3 = &f - &d.double(); // X3 = F - 2 * D
            let y3 = &e * (&d - &x3) - &c.double().double().double(); // Y3 = E * (D - X3) - 8 * C
            let z3 = (y1 * z1).double(); // Z3 = 2 * Y1 * Z1

            debug_assert_eq!(
                E::defining_equation_jacobian(&x3, &y3, &z3),
                FieldElement::<E::BaseField>::zero()
            );
            // SAFETY: The values `x_3, y_3, z_3` are computed correctly to be on the curve.
            // The assertion above verifies that the resulting point is valid.
            Self::new([x3, y3, z3]).unwrap()
        } else {
            // http://www.hyperelliptic.org/EFD/g1p/data/shortw/jacobian-0/doubling/dbl-2009-alnr
            // http://www.hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-0.html#doubling-dbl-2009-l
            let xx = x1.square(); // XX = X1^2
            let yy = y1.square(); // YY = Y1^2
            let yyyy = yy.square(); // YYYY = YY^2
            let zz = z1.square(); // ZZ = Z1^2
            let s = ((x1 + &yy).square() - &xx - &yyyy).double(); // S = 2 * ((X1 + YY)^2 - XX - YYYY)
            let m = &xx.double() + &xx + &E::a() * &zz.square(); // M = 3 * XX + a * ZZ^2
            let x3 = m.square() - &s.double(); // X3 = M^2 - 2 * S
            let y3 = m * (&s - &x3) - &yyyy.double().double().double(); // Y3 = M * (S - X3) - 8 * YYYY
            let z3 = (y1 + z1).square() - &yy - &zz; // Z3 = (Y1 + Z1)^2 - YY - ZZ

            debug_assert_eq!(
                E::defining_equation_jacobian(&x3, &y3, &z3),
                FieldElement::<E::BaseField>::zero()
            );
            // SAFETY: The values `x_3, y_3, z_3` are computed correctly to be on the curve.
            // The assertion above verifies that the resulting point is valid.
            Self::new([x3, y3, z3]).unwrap()
        }
    }
Applies the group operation between a point and itself
crates/math/src/elliptic_curve/short_weierstrass/point.rs