diff --git a/src/bandersnatch/points/extended.zig b/src/bandersnatch/points/extended.zig
index f7a9538..281b3ba 100644
--- a/src/bandersnatch/points/extended.zig
+++ b/src/bandersnatch/points/extended.zig
@@ -70,34 +70,15 @@ pub fn equal(p: ExtendedPoint, q: ExtendedPoint) bool {
 }
 
 pub fn add(p: ExtendedPoint, q: ExtendedPoint) ExtendedPoint {
-    // See "Twisted Edwards Curves Revisited" (https: // eprint.iacr.org/2008/522.pdf)
-    // by Huseyin Hisil, Kenneth Koon-Ho Wong, Gary Carter, and Ed Dawson
-    // 3.1 Unified Addition in E^e
-
-    const x1 = p.x;
-    const y1 = p.y;
-    const t1 = p.t;
-    const z1 = p.z;
-
-    const x2 = q.x;
-    const y2 = q.y;
-    const t2 = q.t;
-    const z2 = q.z;
-
-    const a = Fp.mul(x1, x2);
-
-    const b = Fp.mul(y1, y2);
-
-    const c = Fp.mul(Bandersnatch.D, Fp.mul(t1, t2));
-
-    const d = Fp.mul(z1, z2);
-
-    const h = Fp.sub(b, Fp.mul(a, Bandersnatch.A));
-
-    const e = Fp.sub(Fp.sub(Fp.mul(Fp.add(x1, y1), Fp.add(x2, y2)), b), a);
-
+    // https://hyperelliptic.org/EFD/g1p/auto-twisted-extended.html#addition-add-2008-hwcd
+    const a = Fp.mul(p.x, q.x);
+    const b = Fp.mul(p.y, q.y);
+    const c = Fp.mul(Bandersnatch.D, Fp.mul(p.t, q.t));
+    const d = Fp.mul(p.z, q.z);
+    const e = Fp.sub(Fp.sub(Fp.mul(Fp.add(p.x, p.y), Fp.add(q.x, q.y)), a), b);
     const f = Fp.sub(d, c);
     const g = Fp.add(d, c);
+    const h = Fp.sub(b, a.neg().mulBy5());
 
     return ExtendedPoint{
         .x = Fp.mul(e, f),
@@ -107,6 +88,10 @@ pub fn add(p: ExtendedPoint, q: ExtendedPoint) ExtendedPoint {
     };
 }
 
+inline fn mulByA(x: Fp) Fp {
+    x.neg().mulBy5();
+}
+
 pub fn sub(p: ExtendedPoint, q: ExtendedPoint) ExtendedPoint {
     const neg_q = q.neg();
     return add(p, neg_q);
diff --git a/src/bench.zig b/src/bench.zig
index d1c52e9..edf40b6 100644
--- a/src/bench.zig
+++ b/src/bench.zig
@@ -93,7 +93,7 @@ fn benchIPAs() !void {
     const PrecomputedWeights = @import("polynomial/precomputed_weights.zig").PrecomputedWeights(crs.DomainSize, crs.Domain);
 
     std.debug.print("Setting up IPA benchmark...\n", .{});
-    const N = 500;
+    const N = 100;
 
     var weights = PrecomputedWeights.init();
     const xcrs = crs.CRS.init();
diff --git a/src/fields/fields.zig b/src/fields/fields.zig
index 6b31934..f68d481 100644
--- a/src/fields/fields.zig
+++ b/src/fields/fields.zig
@@ -99,12 +99,20 @@ fn Field(comptime F: type, comptime mod: u256) type {
             return Self{ .fe = ret };
         }
 
-        pub inline fn mul(self: Self, other: Self) Self {
+        pub fn mul(self: Self, other: Self) Self {
             var ret: F.MontgomeryDomainFieldElement = undefined;
             F.mul(&ret, self.fe, other.fe);
             return Self{ .fe = ret };
         }
 
+        pub fn mulBy5(self: Self) Self {
+            var ret: F.MontgomeryDomainFieldElement = undefined;
+            F.add(&ret, self.fe, self.fe);
+            F.add(&ret, ret, ret);
+            F.add(&ret, ret, self.fe);
+            return Self{ .fe = ret };
+        }
+
         pub fn neg(self: Self) Self {
             var ret: F.MontgomeryDomainFieldElement = undefined;
             F.sub(&ret, baseZero.fe, self.fe);
@@ -123,7 +131,7 @@ fn Field(comptime F: type, comptime mod: u256) type {
             return self.mul(self);
         }
 
-        pub inline fn pow2(self: Self, comptime exponent: u8) Self {
+        pub fn pow2(self: Self, comptime exponent: u8) Self {
             var ret = self;
             inline for (exponent) |_| {
                 ret = ret.mul(ret);
@@ -184,7 +192,7 @@ fn Field(comptime F: type, comptime mod: u256) type {
             return std.mem.eql(u64, &self.fe, &other.fe);
         }
 
-        pub inline fn toInteger(self: Self) u256 {
+        pub fn toInteger(self: Self) u256 {
             var non_mont: F.NonMontgomeryDomainFieldElement = undefined;
             F.fromMontgomery(&non_mont, self.fe);