math/mpx-mul4-*.S: Fix up some of the commentary.

  * Fix bogus formatting.

  * Fill in the `...' in the AMD64 version.

  * Explain the common notation and register allocation conventions.

diff --git a/math/mpx-mul4-amd64-sse2.S b/math/mpx-mul4-amd64-sse2.S
index bd8ff2f..0f2dbd3 100644 (file)
--- a/math/mpx-mul4-amd64-sse2.S
+++ b/math/mpx-mul4-amd64-sse2.S
@@ -41,11 +41,11 @@
 /// construct more general variable-length multipliers.
 ///
 /// The basic trick is the same throughout.  In an operand-scanning
 /// construct more general variable-length multipliers.
 ///
 /// The basic trick is the same throughout.  In an operand-scanning

-/// multiplication, the inner multiplication loop multiplies a
-/// multiple-precision operand by a single precision factor, and adds the
-/// result, appropriately shifted, to the result.  A `finely integrated
-/// operand scanning' implementation of Montgomery multiplication also adds
-/// the product of a single-precision `Montgomery factor' and the modulus,
+/// multiplication, the inner multiplication loop multiplies a multiple-
+/// precision operand by a single precision factor, and adds the result,
+/// appropriately shifted, to the result.  A `finely integrated operand
+/// scanning' implementation of Montgomery multiplication also adds the
+/// product of a single-precision `Montgomery factor' and the modulus,

 /// calculated in the same pass.  The more common `coarsely integrated
 /// operand scanning' alternates main multiplication and Montgomery passes,
 /// which requires additional carry propagation.
 /// calculated in the same pass.  The more common `coarsely integrated
 /// operand scanning' alternates main multiplication and Montgomery passes,
 /// which requires additional carry propagation.


@@ -71,22 +71,65 @@
 /// many products together before we must deal with carrying; it also allows
 /// for some calculations to be performed on the above expanded form.
 ///
 /// many products together before we must deal with carrying; it also allows
 /// for some calculations to be performed on the above expanded form.
 ///

-/// ...
+/// We maintain four `carry' registers XMM12--XMM15 accumulating intermediate
+/// results.  The registers' precise roles rotate during the computation; we
+/// name them `c0', `c1', `c2', and `c3'.  Each carry register holds two
+/// 64-bit halves: the register c0, for example, holds c'_0 (low half) and
+/// c''_0 (high half), and represents the value c_0 = c'_0 + c''_0 b; the
+/// carry registers collectively represent the value c_0 + c_1 B + c_2 B^2 +
+/// c_3 B^3.  The `pmuluqdq' instruction acting on a scalar operand
+/// (broadcast across all lanes of its vector) and an operand in the expanded
+/// form above produces a result which can be added directly to the
+/// appropriate carry register.  Following a pass of four multiplications, we
+/// perform some limited carry propagation: let t = c''_0 mod B, and let d =
+/// c'_0 + t b; then we output z = d mod B, add (floor(d/B), floor(c''_0/B))
+/// to c1, and cycle the carry registers around, so that c1 becomes c0, and
+/// the old (implicitly) zeroed c0 becomes c3.

 ///
 ///

-/// We maintain four `carry' registers accumulating intermediate results.
-/// The registers' precise roles rotate during the computation; we name them
-/// `c0', `c1', `c2', and `c3'.  Each carry register holds two 64-bit halves:
-/// the register c0, for example, holds c'_0 (low half) and c''_0 (high
-/// half), and represents the value c_0 = c'_0 + c''_0 b; the carry registers
-/// collectively represent the value c_0 + c_1 B + c_2 B^2 + c_3 B^3.  The
-/// `pmuluqdq' instruction acting on a scalar operand (broadcast across all
-/// lanes of its vector) and an operand in the expanded form above produces a
-/// result which can be added directly to the appropriate carry register.
-/// Following a pass of four multiplications, we perform some limited carry
-/// propagation: let t = c''_0 mod B, and let d = c'_0 + t b; then we output
-/// z = d mod B, add (floor(d/B), floor(c''_0/B)) to c1, and cycle the carry
-/// registers around, so that c1 becomes c0, and the old c0 is (implicitly)
-/// zeroed becomes c3.
+/// On 64-bit AMD64, we have a reasonable number of registers: the expanded
+/// operands are kept in registers.  The packed operands are read from memory
+/// into working registers XMM4 and XMM5; XMM0--XMM3 are used for the actual
+/// multiplications; and XMM6 and XMM7 are used for combining the results.
+/// The following conventional argument names and locations are used
+/// throughout.
+///
+/// Arg        Format      Location    Notes
+///
+/// U  packed      [RAX]
+/// X  packed      [RBX]       In Montgomery multiplication, X = N
+/// V  expanded    XMM8/XMM9
+/// Y  expanded    XMM10/XMM11 In Montgomery multiplication, Y = (A + U V) M
+/// M  expanded    (see below) Montgomery factor, M = -N^{-1} (mod B^4)
+/// N                          Modulus, for Montgomery multiplication
+/// A  packed      [RDI]       Destination/accumulator
+/// C  carry       XMM12--XMM15
+///
+/// The calculation is some variant of
+///
+///    A' + C' B^4 <- U V + X Y + A + C
+///
+/// The low-level functions fit into a fairly traditional (finely-integrated)
+/// operand scanning loop over operand pairs (U, X) (indexed by j) and (V, Y)
+/// (indexed by i).
+///
+/// The variants are as follows.
+///
+/// Function   Variant                 Use             i       j
+///
+/// mmul4      A = C = 0, Y = M        Montgomery      0       0
+/// dmul4      A = 0                   Montgomery      0       +
+/// mmla4      C = 0, Y = M            Montgomery      +       0
+/// dmla4      exactly as shown        Montgomery      +       +
+/// mont4      U = C = 0, V = M        Montgomery      any     0
+///
+/// mul4zc     U = V = A = C = 0       Plain           0       0
+/// mul4       U = V = A = 0           Plain           0       +
+/// mla4zc     U = V = C = 0           Plain           +       0
+/// mla4       U = V = 0               Plain           +       +
+///
+/// The `mmul4' and `mmla4' functions are also responsible for calculating
+/// the Montgomery reduction factor Y = (A + U V) M used by the rest of the
+/// inner loop.

 
 ///--------------------------------------------------------------------------
 /// Macro definitions.
 
 ///--------------------------------------------------------------------------
 /// Macro definitions.


@@ -532,8 +575,8 @@ ENDFUNC
 INTFUNC(mmul4)
        // On entry, RDI points to the destination buffer; RAX and RBX point
        // to the packed operands U and N; and XMM8/XMM9 and XMM10/XMM11 hold
 INTFUNC(mmul4)
        // On entry, RDI points to the destination buffer; RAX and RBX point
        // to the packed operands U and N; and XMM8/XMM9 and XMM10/XMM11 hold

-       // the expanded operands V and M.  The stack pointer must be 8 modulo 16
-       // (as usual for AMD64 ABIs).
+       // the expanded operands V and M.  The stack pointer must be 8 modulo
+       // 16 (as usual for AMD64 ABIs).

        //
        // On exit, we store Y = U V M mod B in XMM10/XMM11, and write the
        // low 128 bits of the sum U V + N Y to [RDI], leaving the remaining
        //
        // On exit, we store Y = U V M mod B in XMM10/XMM11, and write the
        // low 128 bits of the sum U V + N Y to [RDI], leaving the remaining


@@ -750,7 +793,7 @@ ENDFUNC
 
 FUNC(mpx_umul4_amd64_sse2)
        // void mpx_umul4_amd64_sse2(mpw *dv, const mpw *av, const mpw *avl,
 
 FUNC(mpx_umul4_amd64_sse2)
        // void mpx_umul4_amd64_sse2(mpw *dv, const mpw *av, const mpw *avl,

-       //                         const mpw *bv, const mpw *bvl);
+       //                           const mpw *bv, const mpw *bvl);

 
        // Establish the arguments and do initial setup.
        //
 
        // Establish the arguments and do initial setup.
        //

diff --git a/math/mpx-mul4-x86-sse2.S b/math/mpx-mul4-x86-sse2.S
index 2f7b5ec..b1072ff 100644 (file)
--- a/math/mpx-mul4-x86-sse2.S
+++ b/math/mpx-mul4-x86-sse2.S
@@ -40,11 +40,11 @@
 /// construct more general variable-length multipliers.
 ///
 /// The basic trick is the same throughout.  In an operand-scanning
 /// construct more general variable-length multipliers.
 ///
 /// The basic trick is the same throughout.  In an operand-scanning

-/// multiplication, the inner multiplication loop multiplies a
-/// multiple-precision operand by a single precision factor, and adds the
-/// result, appropriately shifted, to the result.  A `finely integrated
-/// operand scanning' implementation of Montgomery multiplication also adds
-/// the product of a single-precision `Montgomery factor' and the modulus,
+/// multiplication, the inner multiplication loop multiplies a multiple-
+/// precision operand by a single precision factor, and adds the result,
+/// appropriately shifted, to the result.  A `finely integrated operand
+/// scanning' implementation of Montgomery multiplication also adds the
+/// product of a single-precision `Montgomery factor' and the modulus,

 /// calculated in the same pass.  The more common `coarsely integrated
 /// operand scanning' alternates main multiplication and Montgomery passes,
 /// which requires additional carry propagation.
 /// calculated in the same pass.  The more common `coarsely integrated
 /// operand scanning' alternates main multiplication and Montgomery passes,
 /// which requires additional carry propagation.


@@ -70,23 +70,64 @@
 /// many products together before we must deal with carrying; it also allows
 /// for some calculations to be performed on the above expanded form.
 ///
 /// many products together before we must deal with carrying; it also allows
 /// for some calculations to be performed on the above expanded form.
 ///

+/// We maintain four `carry' registers XMM4--XMM7 accumulating intermediate
+/// results.  The registers' precise roles rotate during the computation; we
+/// name them `c0', `c1', `c2', and `c3'.  Each carry register holds two
+/// 64-bit halves: the register c0, for example, holds c'_0 (low half) and
+/// c''_0 (high half), and represents the value c_0 = c'_0 + c''_0 b; the
+/// carry registers collectively represent the value c_0 + c_1 B + c_2 B^2 +
+/// c_3 B^3.  The `pmuluqd' instruction acting on a scalar operand (broadcast
+/// across all lanes of its vector) and an operand in the expanded form above
+/// produces a result which can be added directly to the appropriate carry
+/// register.  Following a pass of four multiplications, we perform some
+/// limited carry propagation: let t = c''_0 mod B, and let d = c'_0 + t b;
+/// then we output z = d mod B, add (floor(d/B), floor(c''_0/B)) to c1, and
+/// cycle the carry registers around, so that c1 becomes c0, and the old
+/// (implicitly) zeroed c0 becomes c3.
+///

 /// On 32-bit x86, we are register starved: the expanded operands are kept in
 /// On 32-bit x86, we are register starved: the expanded operands are kept in

-/// memory, typically in warm L1 cache.
+/// memory, typically in warm L1 cache.  The packed operands are read from
+/// memory into working registers XMM0--XMM3 and processed immediately.
+/// The following conventional argument names and locations are used
+/// throughout.
+///
+/// Arg        Format      Location    Notes
+///
+/// U  packed      [EAX]
+/// X  packed      [EBX]       In Montgomery multiplication, X = N
+/// V  expanded    [ECX]
+/// Y  expanded    [EDX]       In Montgomery multiplication, Y = (A + U V) M
+/// M  expanded    [ESI]       -N^{-1} (mod B^4)
+/// N                          Modulus, for Montgomery multiplication
+/// A  packed      [EDI]       Destination/accumulator
+/// C  carry       XMM4--XMM7
+///
+/// The calculation is some variant of
+///
+///    A' + C' B^4 <- U V + X Y + A + C
+///
+/// The low-level functions fit into a fairly traditional (finely-integrated)
+/// operand scanning loop over operand pairs (U, X) (indexed by j) and (V, Y)
+/// (indexed by i).
+///
+/// The variants are as follows.
+///
+/// Function   Variant                 Use             i       j
+///
+/// mmul4      A = C = 0               Montgomery      0       0
+/// dmul4      A = 0                   Montgomery      0       +
+/// mmla4      C = 0                   Montgomery      +       0
+/// dmla4      exactly as shown        Montgomery      +       +
+/// mont4      U = V = C = 0           Montgomery      any     0
+///
+/// mul4zc     U = V = A = C = 0       Plain           0       0
+/// mul4       U = V = A = 0           Plain           0       +
+/// mla4zc     U = V = C = 0           Plain           +       0
+/// mla4       U = V = 0               Plain           +       +

 ///
 ///

-/// We maintain four `carry' registers accumulating intermediate results.
-/// The registers' precise roles rotate during the computation; we name them
-/// `c0', `c1', `c2', and `c3'.  Each carry register holds two 64-bit halves:
-/// the register c0, for example, holds c'_0 (low half) and c''_0 (high
-/// half), and represents the value c_0 = c'_0 + c''_0 b; the carry registers
-/// collectively represent the value c_0 + c_1 B + c_2 B^2 + c_3 B^3.  The
-/// `pmuluqd' instruction acting on a scalar operand (broadcast across all
-/// lanes of its vector) and an operand in the expanded form above produces a
-/// result which can be added directly to the appropriate carry register.
-/// Following a pass of four multiplications, we perform some limited carry
-/// propagation: let t = c''_0 mod B, and let d = c'_0 + t b; then we output
-/// z = d mod B, add (floor(d/B), floor(c''_0/B)) to c1, and cycle the carry
-/// registers around, so that c1 becomes c0, and the old c0 is (implicitly)
-/// zeroed becomes c3.
+/// The `mmul4' and `mmla4' functions are also responsible for calculating
+/// the Montgomery reduction factor Y = (A + U V) M used by the rest of the
+/// inner loop.

 
 ///--------------------------------------------------------------------------
 /// Macro definitions.
 
 ///--------------------------------------------------------------------------
 /// Macro definitions.