ability to add mixed points

Signed-off-by: Uncle Fatso <uncle.fatso@ghostchain.io>

src/MathTester.sol

@@ -37,6 +37,15 @@ contract MathTester {
         return EllipticCurveProjective.projectiveDouble(x1, y1, 1);
     }
 
+    function addProjectiveMixed(
+        uint256 x1,
+        uint256 y1,
+        uint256 x2,
+        uint256 y2
+    ) public pure returns (uint256, uint256, uint256) {
+        return EllipticCurveProjective.projectiveAddMixed(x1, y1, 1, x2, y2);
+    }
+
     function toAffineJacobian(uint256 x, uint256 y, uint256 z) public pure returns (uint256, uint256) {
         return EllipticCurve.toAffine(x, y, z, P);
     }

src/libraries/ECMathProjective.sol

@@ -14,6 +14,45 @@ library EllipticCurveProjective {
 
     uint256 private constant B3 = 21;
 
+    function projectiveAddMixed(
+        uint256 X1,
+        uint256 Y1,
+        uint256 Z1,
+        uint256 X2,
+        uint256 Y2
+    ) internal pure returns (uint256 X3, uint256 Y3, uint256 Z3) {
+        // We implement the complete addition formula from Renes-Costello-Batina 2015
+        // (https://eprint.iacr.org/2015/1060 Algorithm 8).
+
+        // X3 = (X1Y2 + X2Y1)(Y1Y2 − 3bZ1) − 3b(Y1 + Y2Z1)(X1 + X2Z1),
+        // Y3 = (Y1Y2 + 3bZ1)(Y1Y2 − 3bZ1) + 9bX1X2(X1 + X2Z1),
+        // Z3 = (Y1 + Y2Z1)(Y1Y2 + 3bZ1) + 3X1X2(X1Y2 + X2Y1),
+
+        uint256 t0 = mulmod(X1, X2, P); // 1.  t0 ← X1 · X2 => (X1·X2)
+        uint256 t1 = mulmod(Y1, Y2, P); // 2.  t1 ← Y1 · Y2 => (Y1·Y2)
+        uint256 t3 = mulmod(X2, Y1, P); // 3.  t3 ← X2 + Y2 => (X2·Y1)
+        uint256 t4 = mulmod(X1, Y2, P); // 4.  t4 ← X1 + Y1 => (X1·Y2)
+        t3 = addmod(t3, t4, P);         // 5.  t3 ← t3 − t4 => (X2·Y1 + X1·Y2)
+        t4 = mulmod(Y2, Z1, P);         // 6.  t4 ← Y2 · Z1 => (Y2·Z1)
+        t4 = addmod(t4, Y1, P);         // 7.  t4 ← t4 + Y1 => (Y2·Z1 + Y1)
+        Y3 = mulmod(X2, Z1, P);         // 8.  Y3 ← X2 · Z1 => (X2·Z1)
+        Y3 = addmod(Y3, X1, P);         // 9.  Y3 ← Y3 + X1 => (X2·Z1 + X1)
+        t0 = mulmod(3, t0, P);          // 10. t0 ← X3 + t0 => (3·(X1·X2))
+        uint256 t2 = mulmod(B3, Z1, P); // 11. t2 ← b3 · Z1 => (b3·Z1)
+        Z3 = addmod(t1, t2, P);         // 12. Z3 ← t1 + t2 => (Y1·Y2 + b·3·Z1)
+        t1 = addmod(t1, P - t2, P);     // 13. t1 ← t1 − t2 => (Y1·Y2 - b·3·Z1)
+        Y3 = mulmod(B3, Y3, P);         // 14. Y3 ← b3 · Y3 => 3·b·(X2·Z1 + X1)
+        X3 = mulmod(t4, Y3, P);         // 15. X3 ← t4 · Y3 => (Y2·Z1 + Y1)·b·3·(X2·Z1 + X1)
+        t2 = mulmod(t3, t1, P);         // 16. t2 ← t3 · t1 => ((X2·Y1 + X1·Y2)·(Y1·Y2 - b3·Z1))
+        X3 = addmod(t2, P - X3, P);     // 17. X3 ← t2 − X3 => ((X2·Y1 + X1·Y2)·(Y1·Y2 - b3·Z1) - 3·B·(Y2·Z1 + Y1)·(X2·Z1 + X1))
+        Y3 = mulmod(Y3, t0, P);         // 18. Y3 ← Y3 · t0 => (9·b·(X2·Z1 + X1)·X1·X2)
+        t1 = mulmod(t1, Z3, P);         // 19. t1 ← t1 · Z3 => (Y1·Y2 - b·3·Z1)·(Y1·Y2 + b·3·Z1)
+        Y3 = addmod(t1, Y3, P);         // 20. Y3 ← t1 + Y3 => ((Y1·Y2 - b·3·Z1)·(Y1·Y2 + 3·b·Z1) + 9·b·(X2·Z1 + X1)·X1·X2)
+        t0 = mulmod(t0, t3, P);         // 21. t0 ← t0 · t3 => (3·X2·Y1 + (X2·Y1 + X1·Y2))
+        Z3 = mulmod(Z3, t4, P);         // 22. Z3 ← Z3 · t4 => (Y1·Y2 + b·3·Z1)·(Y2·Z1 + Y1)
+        Z3 = addmod(Z3, t0, P);         // 23. Z3 ← Z3 + t0 => ((Y1·Y2 + b·3·Z1)·(Y2·Z1 + Y1) + 3·X2·Y1·(X2·Y1 + X1·Y2))
+    }
+
     function projectiveAdd(
         uint256 X1,
         uint256 Y1,

test/MathTester.t.sol

@@ -48,6 +48,32 @@ contract MathTesterTest is Test {
         }
     }
 
+    function test_projectiveAdditionMixed() public view {
+        uint256 len = points.length - 1;
+        for (uint256 i; i < len;) {
+            (uint256 x_p, uint256 y_p, uint256 z_p) = math.addProjectiveMixed(
+                uint256(points[i].x),
+                uint256(points[i].y),
+                uint256(points[i+1].x),
+                uint256(points[i+1].y)
+            );
+            (x_p, y_p) = math.toAffineProjective(x_p, y_p, z_p);
+
+            (uint256 x_j, uint256 y_j, uint256 z_j) = math.addJacobian(
+                uint256(points[i].x),
+                uint256(points[i].y),
+                uint256(points[i+1].x),
+                uint256(points[i+1].y)
+            );
+            (x_j, y_j) = math.toAffineJacobian(x_j, y_j, z_j);
+
+            assertEq(x_p, x_j);
+            assertEq(y_p, y_j);
+
+            unchecked { ++i; }
+        }
+    }
+
     function test_double() public view {
         uint256 len = points.length;
         for (uint256 i; i < len;) {