// use Karatsuba's identity here, but I suspect that loses more in
// the shifting, bit-twiddling, and dependency chains that it gains
// in saving a multiplication which otherwise pipelines well.
- // q0 = // (u_0; u_1)
- // q1 = // (v_0; v_1)
+ // q0 = // (u_1; u_0)
+ // q1 = // (v_1; v_0)
vmull.p64 q2, d1, d2 // u_1 v_0
vmull.p64 q3, d0, d3 // u_0 v_1
- vmull.p64 q8, d1, d3 // (x_3; t_1) = u_1 v_1
- vmull.p64 q9, d0, d2 // (t_0; x_0) = u_0 v_0
+ vmull.p64 q8, d1, d3 // (t_1; x_3) = u_1 v_1
+ vmull.p64 q9, d0, d2 // (x_0; t_0) = u_0 v_0
// Arrange the pieces to form a double-precision polynomial.
- veor q2, q2, q3 // (m_1; m_0) = u_0 v_1 + u_1 v_0
+ veor q2, q2, q3 // (m_0; m_1) = u_0 v_1 + u_1 v_0
veor d17, d17, d4 // x_2 = t_1 + m_1
veor d18, d18, d5 // x_1 = t_0 + m_0
- // q8 = // (x_3; x_2)
- // q9 = // (x_1; x_0)
+ // q8 = // (x_2; x_3)
+ // q9 = // (x_0; x_1)
// One-and-a-half problems remain.
//
// This is an inconvenient size. There's nothing for it but to do
// four multiplications, as if for the 128-bit case.
- // q0 = // (u_0 + u_1 t^32; u_2)
- // q1 = // (v_0 + v_1 t^32; v_2)
+ // q0 = // (u_2; u_0 + u_1 t^32)
+ // q1 = // (v_2; v_0 + v_1 t^32)
vmull.p64 q8, d1, d2 // u_2 (v_0 + v_1 t^32) = e_0
vmull.p64 q9, d0, d3 // v_2 (u_0 + u_1 t^32) = e_1
- vmull.p64 q3, d1, d3 // u_2 v_2 t^64 = d = (0; d)
+ vmull.p64 q3, d1, d3 // u_2 v_2 t^64 = d = (d; 0)
vmull.p64 q0, d0, d2 // u_0 v_0 + (u_0 v_1 + u_1 v_0) t^32
// + u_1 v_1 t^64 = f
veor q11, q11, q13 // b = u_1 v_2 + u_2 v_1
// Piece the product together.
- veor d17, d17, d22 // q8 = // (x_5; x_4)
+ veor d17, d17, d22 // q8 = // (x_4; x_5)
veor d18, d18, d23
- veor d19, d19, d24 // q9 = // (x_3; x_2)
- veor d20, d20, d25 // q10 = // (x_1; x_0)
+ veor d19, d19, d24 // q9 = // (x_2; x_3)
+ veor d20, d20, d25 // q10 = // (x_0; x_1)
// Next, the reduction. Our polynomial this time is p(x) = t^192 +
// t^7 + t^2 + t + 1. Yes, the magic numbers are the same as the
// 128-bit case. I don't know why.
// First, shift the high bits down.
- // q8 = // (y_5; y_4)
- // q9 = // (y_3; y_2)
- // q10 = // (y_1; y_0)
- vshl.u64 q11, q8, #63 // (y_5; y_4) b_i for t
+ // q8 = // (y_4; y_5)
+ // q9 = // (y_2; y_3)
+ // q10 = // (y_0; y_1)
+ vshl.u64 q11, q8, #63 // (y_4; y_5) b_i for t
vshl.u64 d28, d18, #63 // y_3 b_i for t
- vshl.u64 q12, q8, #62 // (y_5; y_4) b_i for t^2
+ vshl.u64 q12, q8, #62 // (y_4; y_5) b_i for t^2
vshl.u64 d29, d18, #62 // y_3 b_i for t^2
- vshl.u64 q13, q8, #57 // (y_5; y_4) b_i for t^7
+ vshl.u64 q13, q8, #57 // (y_4; y_5) b_i for t^7
vshl.u64 d30, d18, #57 // y_3 b_i for t^7
veor q11, q11, q12 // mix them all together
veor d28, d28, d29
// And finally shift the low bits up. Also, switch the order of the
// pieces for output.
- // q8 = // (y'_5; y'_4)
- // q9 = // (y'_3; y'_2)
- // q10 = // (y'_1; y'_0)
- vshr.u64 q11, q8, #1 // (y_5; y_4) a_i for t
+ // q8 = // (y'_4; y'_5)
+ // q9 = // (y'_2; y'_3)
+ // q10 = // (y'_0; y'_1)
+ vshr.u64 q11, q8, #1 // (y_4; y_5) a_i for t
vshr.u64 d28, d18, #1 // y'_3 a_i for t
- vshr.u64 q12, q8, #2 // (y_5; y_4) a_i for t^2
+ vshr.u64 q12, q8, #2 // (y_4; y_5) a_i for t^2
vshr.u64 d29, d18, #2 // y'_3 a_i for t^2
- vshr.u64 q13, q8, #7 // (y_5; y_4) a_i for t^7
+ vshr.u64 q13, q8, #7 // (y_4; y_5) a_i for t^7
vshr.u64 d30, d18, #7 // y'_3 a_i for t^7
veor q8, q8, q11
veor d18, d18, d28
// 128-bit multiplications already, and Karatsuba is too annoying
// there, so there'll be 12 multiplications altogether, rather than
// the 16 we'd have if we did this the naïve way.
- // q0 = // u_0 = (u_00; u_01)
- // q1 = // u_1 = (u_10; u_11)
- // q2 = // v_0 = (v_00; v_01)
- // q3 = // v_1 = (v_10; v_11)
+ // q0 = // u_0 = (u_01; u_00)
+ // q1 = // u_1 = (u_11; u_10)
+ // q2 = // v_0 = (v_01; v_00)
+ // q3 = // v_1 = (v_11; v_10)
- veor q8, q0, q1 // u_* = (u_00 + u_10; u_01 + u_11)
- veor q9, q2, q3 // v_* = (v_00 + v_10; v_01 + v_11)
+ veor q8, q0, q1 // u_* = (u_01 + u_11; u_00 + u_10)
+ veor q9, q2, q3 // v_* = (v_01 + v_11; v_00 + v_10)
// Start by building the cross product, q = u_* v_*.
vmull.p64 q14, d16, d19 // u_*0 v_*1
// The polynomial this time is p(t) = t^256 + t^10 + t^5 + t^2 + 1.
// First, shift the high bits down.
- // q8 = // (y_7; y_6)
- // q9 = // (y_5; y_4)
- // q10 = // (y_3; y_2)
- // q11 = // (y_1; y_0)
- vshl.u64 q0, q8, #62 // (y_7; y_6) b_i for t^2
- vshl.u64 q12, q9, #62 // (y_5; y_4) b_i for t^2
- vshl.u64 q1, q8, #59 // (y_7; y_6) b_i for t^5
- vshl.u64 q13, q9, #59 // (y_5; y_4) b_i for t^5
- vshl.u64 q2, q8, #54 // (y_7; y_6) b_i for t^10
- vshl.u64 q14, q9, #54 // (y_5; y_4) b_i for t^10
+ // q8 = // (y_6; y_7)
+ // q9 = // (y_4; y_5)
+ // q10 = // (y_2; y_3)
+ // q11 = // (y_0; y_1)
+ vshl.u64 q0, q8, #62 // (y_6; y_7) b_i for t^2
+ vshl.u64 q12, q9, #62 // (y_4; y_5) b_i for t^2
+ vshl.u64 q1, q8, #59 // (y_6; y_7) b_i for t^5
+ vshl.u64 q13, q9, #59 // (y_4; y_5) b_i for t^5
+ vshl.u64 q2, q8, #54 // (y_6; y_7) b_i for t^10
+ vshl.u64 q14, q9, #54 // (y_4; y_5) b_i for t^10
veor q0, q0, q1 // mix the contributions together
veor q12, q12, q13
veor q0, q0, q2
// And then shift the low bits up. Also, switch the order of the
// pieces for output.
- // q8 = // (y'_7; y'_6)
- // q9 = // (y'_5; y'_4)
- // q10 = // (y'_3; y'_2)
- // q11 = // (y'_1; y'_0)
- vshr.u64 q0, q8, #2 // (y_7; y_6) a_i for t^2
- vshr.u64 q12, q9, #2 // (y_5; y'_4) a_i for t^2
- vshr.u64 q1, q8, #5 // (y_7; y_6) a_i for t^5
- vshr.u64 q13, q9, #5 // (y_5; y_4) a_i for t^5
- vshr.u64 q2, q8, #10 // (y_7; y_6) a_i for t^10
- vshr.u64 q14, q9, #10 // (y_5; y_4) a_i for t^10
+ // q8 = // (y'_6; y'_7)
+ // q9 = // (y'_4; y'_5)
+ // q10 = // (y'_2; y'_3)
+ // q11 = // (y'_0; y'_1)
+ vshr.u64 q0, q8, #2 // (y_6; y_7) a_i for t^2
+ vshr.u64 q12, q9, #2 // (y'_4; y_5) a_i for t^2
+ vshr.u64 q1, q8, #5 // (y_6; y_7) a_i for t^5
+ vshr.u64 q13, q9, #5 // (y_4; y_5) a_i for t^5
+ vshr.u64 q2, q8, #10 // (y_6; y_7) a_i for t^10
+ vshr.u64 q14, q9, #10 // (y_4; y_5) a_i for t^10
veor q8, q8, q0 // mix the contributions together
veor q1, q1, q2
~mdw / catacomb / blobdiff