diff --git a/src/backend/vector/ifma/edwards.rs b/src/backend/vector/ifma/edwards.rs
index 5339f05..8e65f2f 100644
--- a/src/backend/vector/ifma/edwards.rs
+++ b/src/backend/vector/ifma/edwards.rs
@@ -51,7 +51,50 @@ impl From<ExtendedPoint> for CachedPoint {
         CachedPoint(F51x4Reduced::from(x))
     }
 }
-    
+
+impl ExtendedPoint {
+    pub fn double(&self) -> ExtendedPoint {
+        // Set tmp0 = (X1 Y1 X1 Y1)
+        let mut tmp0 = self.0.shuffle(Shuffle::ABAB);
+
+        // Set tmp1 = (Y1 X1 Y1 X1)
+        let mut tmp1 = tmp0.shuffle(Shuffle::BADC);
+
+        // Set tmp0 = (X1 Y1 Z1 X1+Y1)
+        tmp0 = self.0.blend(&(tmp0 + tmp1), Lanes::D);
+
+        tmp1 = F51x4Reduced::from(tmp0).square();
+        // Now tmp1 = (S1 S2 S3 S4)
+
+        // We want to compute
+        //
+        //    + | S1 | S1 | S1 | S1 |
+        //    + | S2 |    |    | S2 |
+        //    + |    |    | S3 |    |
+        //    + |    |    | S3 |    |
+        //    + |    |16p |16p |16p |
+        //    - |    | S2 | S2 |    |
+        //    - |    |    |    | S4 |
+        //    =======================
+        //        S5   S6   S8   S9
+
+        let zero = F51x4Unreduced::zero();
+
+        let S1_S1_S1_S1 = tmp1.shuffle(Shuffle::AAAA);
+        let S2_S2_S2_S2 = tmp1.shuffle(Shuffle::BBBB);
+
+        let S2_S2_S2_S4 = S2_S2_S2_S2.blend(&tmp1, Lanes::D).negate_lazy();
+
+        tmp0 = S1_S1_S1_S1 + zero.blend(&(tmp1 + tmp1), Lanes::C);
+        tmp0 = tmp0 + zero.blend(&S2_S2_S2_S2, Lanes::AD);
+        tmp0 = tmp0 + zero.blend(&S2_S2_S2_S4, Lanes::BCD);
+
+        let tmp2 = F51x4Reduced::from(tmp0);
+
+        ExtendedPoint(&tmp2.shuffle(Shuffle::DBBD) * &tmp2.shuffle(Shuffle::CACA))
+    }
+}
+
 impl<'a, 'b> Add<&'b CachedPoint> for &'a ExtendedPoint {
     type Output = ExtendedPoint;
 
@@ -134,4 +177,35 @@ mod test {
         let Q = &constants::ED25519_BASEPOINT_TABLE * &Scalar::from(8475983829u64);
         addition_test_helper(P, Q);
     }
+
+    fn doubling_test_helper(P: edwards::EdwardsPoint) {
+        //let R1: edwards::EdwardsPoint = serial_double(P.into()).into();
+        let R2: edwards::EdwardsPoint = ExtendedPoint::from(P).double().into();
+        println!("Testing point doubling:");
+        println!("P = {:?}", P);
+        //println!("(serial) R1 = {:?}", R1);
+        println!("(vector) R2 = {:?}", R2);
+        println!("P + P = {:?}", &P + &P);
+        //assert_eq!(R1.compress(), (&P + &P).compress());
+        assert_eq!(R2.compress(), (&P + &P).compress());
+        println!("OK!\n");
+    }
+
+    #[test]
+    fn vector_doubling_vs_serial_doubling_vs_edwards_extendedpoint() {
+        use constants;
+        use scalar::Scalar;
+
+        println!("Testing [2]id");
+        let P = edwards::EdwardsPoint::identity();
+        doubling_test_helper(P);
+
+        println!("Testing [2]B");
+        let P = constants::ED25519_BASEPOINT_POINT;
+        doubling_test_helper(P);
+
+        println!("Testing [2]([k]B)");
+        let P = &constants::ED25519_BASEPOINT_TABLE * &Scalar::from(8475983829u64);
+        doubling_test_helper(P);
+    }
 }
diff --git a/src/backend/vector/ifma/field.rs b/src/backend/vector/ifma/field.rs
index d344514..40220c5 100644
--- a/src/backend/vector/ifma/field.rs
+++ b/src/backend/vector/ifma/field.rs
@@ -28,20 +28,17 @@ pub struct F51x4Unreduced(pub(crate) [u64x4; 5]);
 #[derive(Copy, Clone, Debug)]
 pub struct F51x4Reduced(pub(crate) [u64x4; 5]);
 
-#[derive(Copy, Clone)]
-pub enum Lanes {
-    D,
-    AB,
-    AC,
-}
-
 #[derive(Copy, Clone)]
 pub enum Shuffle {
     AAAA,
+    BBBB,
     BADC,
     ADDA,
     CBCB,
     ABDC,
+    ABAB,
+    DBBD,
+    CACA,
 }
 
 #[inline(always)]
@@ -51,14 +48,28 @@ fn shuffle_lanes(x: u64x4, control: Shuffle) -> u64x4 {
 
         match control {
             Shuffle::AAAA => perm(x.into_bits(), 0b00_00_00_00).into_bits(),
+            Shuffle::BBBB => perm(x.into_bits(), 0b01_01_01_01).into_bits(),
             Shuffle::BADC => perm(x.into_bits(), 0b10_11_00_01).into_bits(),
             Shuffle::ADDA => perm(x.into_bits(), 0b00_11_11_00).into_bits(),
             Shuffle::CBCB => perm(x.into_bits(), 0b01_10_01_10).into_bits(),
             Shuffle::ABDC => perm(x.into_bits(), 0b10_11_01_00).into_bits(),
+            Shuffle::ABAB => perm(x.into_bits(), 0b01_00_01_00).into_bits(),
+            Shuffle::DBBD => perm(x.into_bits(), 0b11_01_01_11).into_bits(),
+            Shuffle::CACA => perm(x.into_bits(), 0b00_10_00_10).into_bits(),
         }
     }
 }
 
+#[derive(Copy, Clone)]
+pub enum Lanes {
+    D,
+    C,
+    AB,
+    AC,
+    AD,
+    BCD,
+}
+
 #[inline]
 fn blend_lanes(x: u64x4, y: u64x4, control: Lanes) -> u64x4 {
     unsafe {
@@ -66,13 +77,20 @@ fn blend_lanes(x: u64x4, y: u64x4, control: Lanes) -> u64x4 {
 
         match control {
             Lanes::D => blend(x.into_bits(), y.into_bits(), 0b11_00_00_00).into_bits(),
+            Lanes::C => blend(x.into_bits(), y.into_bits(), 0b00_11_00_00).into_bits(),
             Lanes::AB => blend(x.into_bits(), y.into_bits(), 0b00_00_11_11).into_bits(),
             Lanes::AC => blend(x.into_bits(), y.into_bits(), 0b00_11_00_11).into_bits(),
+            Lanes::AD => blend(x.into_bits(), y.into_bits(), 0b11_00_00_11).into_bits(),
+            Lanes::BCD => blend(x.into_bits(), y.into_bits(), 0b11_11_11_00).into_bits(),
         }
     }
 }
 
 impl F51x4Unreduced {
+    pub fn zero() -> F51x4Unreduced {
+        F51x4Unreduced([u64x4::splat(0); 5])
+    }
+
     pub fn new(
         x0: &FieldElement51,
         x1: &FieldElement51,