/* This Source Code Form is subject to the terms of the Mozilla Public

 * License, v. 2.0. If a copy of the MPL was not distributed with this

 * file, You can obtain one at http://mozilla.org/MPL/2.0/. */

/*

 *

 *  This PA-RISC 2.0 function computes the product of two unsigned integers,

 *  and adds the result to a previously computed integer.  The multiplicand

 *  is a 512-bit (64-byte, eight doubleword) unsigned integer, stored in

 *  memory in little-double-wordian order.  The multiplier is an unsigned

 *  64-bit integer.  The previously computed integer to which the product is

 *  added is located in the result ("res") area, and is assumed to be a

 *  576-bit (72-byte, nine doubleword) unsigned integer, stored in memory

 *  in little-double-wordian order.  This value normally will be the result

 *  of a previously computed nine doubleword result.  It is not necessary

 *  to pad the multiplicand with an additional 64-bit zero doubleword.

 *

 *  Multiplicand, multiplier, and addend ideally should be aligned at

 *  16-byte boundaries for best performance.  The code will function

 *  correctly for alignment at eight-byte boundaries which are not 16-byte

 *  boundaries, but the execution may be slightly slower due to even/odd

 *  bank conflicts on PA-RISC 8000 processors.

 *

 *  This function is designed to accept the same calling sequence as Bill

 *  Ackerman's "maxpy_little" function.  The carry from the ninth doubleword

 *  of the result is written to the tenth word of the result, as is done by

 *  Bill Ackerman's function.  The final carry also is returned as an

 *  integer, which may be ignored.  The function prototype may be either

 *  of the following:

 *

 *      void multacc512( int l, chunk* m, const chunk* a, chunk* res );

 *          or

 *      int multacc512( int l, chunk* m, const chunk* a, chunk* res );

 *

 *  where:  "l" originally denoted vector lengths.  This parameter is

 *      ignored.  This function always assumes a multiplicand length of

 *      512 bits (eight doublewords), and addend and result lengths of

 *      576 bits (nine doublewords).

 *

 *      "m" is a pointer to the doubleword multiplier, ideally aligned

 *      on a 16-byte boundary.

 *

 *      "a" is a pointer to the eight-doubleword multiplicand, stored

 *      in little-double-wordian order, and ideally aligned on a 16-byte

 *      boundary.

 *

 *      "res" is a pointer to the nine doubleword addend, and to the

 *      nine-doubleword product computed by this function.  The result

 *      also is stored in little-double-wordian order, and ideally is

 *      aligned on a 16-byte boundary. It is expected that the alignment

 *      of the "res" area may alternate between even/odd doubleword

 *      boundaries for successive calls for 512-bit x 512-bit

 *      multiplications.

 *

 *  The code for this function has been scheduled to use the parallelism

 *  of the PA-RISC 8000 series microprocessors as well as the author was

 *  able.  Comments and/or suggestions for improvement are welcomed.

 *

 *  The code is "64-bit safe".  This means it may be called in either

 *  the 32ILP context or the 64LP context.  All 64-bits of registers are

 *  saved and restored.

 *

 *  This code is self-contained.  It requires no other header files in order

 *  to compile and to be linkable on a PA-RISC 2.0 machine.  Symbolic

 *  definitions for registers and stack offsets are included within this

 *  one source file.

 *

 *  This is a leaf routine.  As such, minimal use is made of the stack area.

 *  Of the 192 bytes allocated, 64 bytes are used for saving/restoring eight

 *  general registers, and 128 bytes are used to move intermediate products

 *  from the floating-point registers to the general registers.  Stack

 *  protocols assure proper alignment of these areas.

 *

 */

/*  ====================================================================*/

/*      symbolic definitions for PA-RISC registers      */

/*      in the MIPS style, avoids lots of case shifts       */

/*      assigments (except t4) preserve register number parity  */

/*  ====================================================================*/

#define zero    %r0         /* permanent zero */

#define t5      %r1         /* temp register, altered by addil */

#define rp      %r2         /* return pointer */

#define s1      %r3         /* callee saves register*/

#define s0      %r4         /* callee saves register*/

#define s3      %r5         /* callee saves register*/

#define s2      %r6         /* callee saves register*/

#define s5      %r7         /* callee saves register*/

#define s4      %r8         /* callee saves register*/

#define s7      %r9         /* callee saves register*/

#define s6      %r10        /* callee saves register*/

#define t1      %r19        /* caller saves register*/

#define t0      %r20        /* caller saves register*/

#define t3      %r21        /* caller saves register*/

#define t2      %r22        /* caller saves register*/

#define a3      %r23        /* fourth argument register, high word */

#define a2      %r24        /* third argument register, low word*/

#define a1      %r25        /* second argument register, high word*/

#define a0      %r26        /* first argument register, low word*/

#define v0      %r28        /* high order return value*/

#define v1      %r29        /* low order return value*/

#define sp      %r30        /* stack pointer*/

#define t4      %r31        /* temporary register   */

#define fa0     %fr4        /* first argument register*/

#define fa1     %fr5        /* second argument register*/

#define fa2     %fr6        /* third argument register*/

#define fa3     %fr7        /* fourth argument register*/

#define fa0r    %fr4R       /* first argument register*/

#define fa1r    %fr5R       /* second argument register*/

#define fa2r    %fr6R       /* third argument register*/

#define fa3r    %fr7R       /* fourth argument register*/

#define ft0     %fr8        /* caller saves register*/

#define ft1     %fr9        /* caller saves register*/

#define ft2     %fr10       /* caller saves register*/

#define ft3     %fr11       /* caller saves register*/

#define ft0r    %fr8R       /* caller saves register*/

#define ft1r    %fr9R       /* caller saves register*/

#define ft2r    %fr10R      /* caller saves register*/

#define ft3r    %fr11R      /* caller saves register*/

#define ft4     %fr22       /* caller saves register*/

#define ft5     %fr23       /* caller saves register*/

#define ft6     %fr24       /* caller saves register*/

#define ft7     %fr25       /* caller saves register*/

#define ft8     %fr26       /* caller saves register*/

#define ft9     %fr27       /* caller saves register*/

#define ft10    %fr28       /* caller saves register*/

#define ft11    %fr29       /* caller saves register*/

#define ft12    %fr30       /* caller saves register*/

#define ft13    %fr31       /* caller saves register*/

#define ft4r    %fr22R      /* caller saves register*/

#define ft5r    %fr23R      /* caller saves register*/

#define ft6r    %fr24R      /* caller saves register*/

#define ft7r    %fr25R      /* caller saves register*/

#define ft8r    %fr26R      /* caller saves register*/

#define ft9r    %fr27R      /* caller saves register*/

#define ft10r   %fr28R      /* caller saves register*/

#define ft11r   %fr29R      /* caller saves register*/

#define ft12r   %fr30R      /* caller saves register*/

#define ft13r   %fr31R      /* caller saves register*/

/*  ================================================================== */

/*      functional definitions for PA-RISC registers           */

/*  ================================================================== */

/*              general registers           */

#define T1      a0          /* temp, (length parameter ignored)             */

#define pM      a1          /* -> 64-bit multiplier                         */

#define T2      a1          /* temp, (after fetching multiplier)            */

#define pA      a2          /* -> multiplicand vector (8 64-bit words)      */

#define T3      a2          /* temp, (after fetching multiplicand)          */

#define pR      a3          /* -> addend vector (8 64-bit doublewords,

                                  result vector (9 64-bit words)            */

#define S0      s0          /* callee saves summand registers               */

#define S1      s1

#define S2      s2

#define S3      s3

#define S4      s4

#define S5      s5

#define S6      s6

#define S7      s7

#define S8      v0          /* caller saves summand registers               */

#define S9      v1

#define S10     t0

#define S11     t1

#define S12     t2

#define S13     t3

#define S14     t4

#define S15     t5

/*              floating-point registers                                    */

#define M       fa0         /* multiplier double word                       */

#define MR      fa0r        /* low order half of multiplier double word     */

#define ML      fa0         /* high order half of multiplier double word    */

#define A0      fa2         /* multiplicand double word 0                   */

#define A0R     fa2r        /* low order half of multiplicand double word   */

#define A0L     fa2         /* high order half of multiplicand double word  */

#define A1      fa3         /* multiplicand double word 1                   */

#define A1R     fa3r        /* low order half of multiplicand double word   */

#define A1L     fa3         /* high order half of multiplicand double word  */

#define A2      ft0         /* multiplicand double word 2                   */

#define A2R     ft0r        /* low order half of multiplicand double word   */

#define A2L     ft0         /* high order half of multiplicand double word  */

#define A3      ft1         /* multiplicand double word 3                   */

#define A3R     ft1r        /* low order half of multiplicand double word   */

#define A3L     ft1         /* high order half of multiplicand double word  */

#define A4      ft2         /* multiplicand double word 4                   */

#define A4R     ft2r        /* low order half of multiplicand double word   */

#define A4L     ft2         /* high order half of multiplicand double word  */

#define A5      ft3         /* multiplicand double word 5                   */

#define A5R     ft3r        /* low order half of multiplicand double word   */

#define A5L     ft3         /* high order half of multiplicand double word  */

#define A6      ft4         /* multiplicand double word 6                   */

#define A6R     ft4r        /* low order half of multiplicand double word   */

#define A6L     ft4         /* high order half of multiplicand double word  */

#define A7      ft5         /* multiplicand double word 7                   */

#define A7R     ft5r        /* low order half of multiplicand double word   */

#define A7L     ft5         /* high order half of multiplicand double word  */

#define P0      ft6         /* product word 0                               */

#define P1      ft7         /* product word 0                               */

#define P2      ft8         /* product word 0                               */

#define P3      ft9         /* product word 0                               */

#define P4      ft10        /* product word 0                               */

#define P5      ft11        /* product word 0                               */

#define P6      ft12        /* product word 0                               */

#define P7      ft13        /* product word 0                               */

/*  ======================================================================  */

/*      symbolic definitions for HP-UX stack offsets                        */

/*      symbolic definitions for memory NOPs                                */

/*  ======================================================================  */

#define ST_SZ       192         /* stack area total size                    */

#define SV0         -192(sp)    /* general register save area               */

#define SV1         -184(sp)

#define SV2         -176(sp)

#define SV3         -168(sp)

#define SV4         -160(sp)

#define SV5         -152(sp)

#define SV6         -144(sp)

#define SV7         -136(sp)

#define XF0         -128(sp)    /* data transfer area                       */

#define XF1         -120(sp)    /* for floating-pt to integer regs          */

#define XF2         -112(sp)

#define XF3         -104(sp)

#define XF4         -96(sp)

#define XF5         -88(sp)

#define XF6         -80(sp)

#define XF7         -72(sp)

#define XF8         -64(sp)

#define XF9         -56(sp)

#define XF10        -48(sp)

#define XF11        -40(sp)

#define XF12        -32(sp)

#define XF13        -24(sp)

#define XF14        -16(sp)

#define XF15        -8(sp)

#define mnop    proberi (sp),3,zero     /* memory NOP                       */

/*  ======================================================================  */

/*      assembler formalities                                               */

/*  ======================================================================  */

#ifdef __LP64__

                .level  2.0W

#else

                .level  2.0

#endif

                .space    $TEXT$

                .subspa   $CODE$

                .align    16

/*  ======================================================================  */

/*      here to compute 64-bit x 512-bit product + 512-bit addend           */

/*  ======================================================================  */

multacc512

        .PROC

        .CALLINFO

        .ENTRY

    fldd    0(pM),M                 ; multiplier double word

    ldo     ST_SZ(sp),sp            ; push stack

    fldd    0(pA),A0                ; multiplicand double word 0

    std     S1,SV1                  ; save s1

    fldd    16(pA),A2               ; multiplicand double word 2

    std     S3,SV3                  ; save s3

    fldd    32(pA),A4               ; multiplicand double word 4

    std     S5,SV5                  ; save s5

    fldd    48(pA),A6               ; multiplicand double word 6

    std     S7,SV7                  ; save s7

    std     S0,SV0                  ; save s0

    fldd    8(pA),A1                ; multiplicand double word 1

    xmpyu   MR,A0L,P0               ; A0 cross 32-bit word products

    xmpyu   ML,A0R,P2

    std     S2,SV2                  ; save s2

    fldd    24(pA),A3               ; multiplicand double word 3

    xmpyu   MR,A2L,P4               ; A2 cross 32-bit word products

    xmpyu   ML,A2R,P6

    std     S4,SV4                  ; save s4

    fldd    40(pA),A5               ; multiplicand double word 5

    std     S6,SV6                  ; save s6

    fldd    56(pA),A7               ; multiplicand double word 7

    fstd    P0,XF0                  ; MR * A0L

    xmpyu   MR,A0R,P0               ; A0 right 32-bit word product

    xmpyu   MR,A1L,P1               ; A1 cross 32-bit word product

    fstd    P2,XF2                  ; ML * A0R

    xmpyu   ML,A0L,P2               ; A0 left 32-bit word product

    xmpyu   ML,A1R,P3               ; A1 cross 32-bit word product

    fstd    P4,XF4                  ; MR * A2L

    xmpyu   MR,A2R,P4               ; A2 right 32-bit word product

    xmpyu   MR,A3L,P5               ; A3 cross 32-bit word product

    fstd    P6,XF6                  ; ML * A2R

    xmpyu   ML,A2L,P6               ; A2 parallel 32-bit word product

    xmpyu   ML,A3R,P7               ; A3 cross 32-bit word product

    ldd     XF0,S0                  ; MR * A0L

    fstd    P1,XF1                  ; MR * A1L

    ldd     XF2,S2                  ; ML * A0R

    fstd    P3,XF3                  ; ML * A1R

    ldd     XF4,S4                  ; MR * A2L

    fstd    P5,XF5                  ; MR * A3L

    xmpyu   MR,A1R,P1               ; A1 parallel 32-bit word products

    xmpyu   ML,A1L,P3

    ldd     XF6,S6                  ; ML * A2R

    fstd    P7,XF7                  ; ML * A3R

    xmpyu   MR,A3R,P5               ; A3 parallel 32-bit word products

    xmpyu   ML,A3L,P7

    fstd    P0,XF0                  ; MR * A0R

    ldd     XF1,S1                  ; MR * A1L

    nop

    add     S0,S2,T1                ; A0 cross product sum

    fstd    P2,XF2                  ; ML * A0L

    ldd     XF3,S3                  ; ML * A1R

    add,dc  zero,zero,S0            ; A0 cross product sum carry

    depd,z  T1,31,32,S2             ; A0 cross product sum << 32

    fstd    P4,XF4                  ; MR * A2R

    ldd     XF5,S5                  ; MR * A3L

    shrpd   S0,T1,32,S0             ; A0 carry | cross product sum >> 32

    add     S4,S6,T3                ; A2 cross product sum

    fstd    P6,XF6                  ; ML * A2L

    ldd     XF7,S7                  ; ML * A3R

    add,dc  zero,zero,S4            ; A2 cross product sum carry

    depd,z  T3,31,32,S6             ; A2 cross product sum << 32

    ldd     XF0,S8                  ; MR * A0R

    fstd    P1,XF1                  ; MR * A1R

    xmpyu   MR,A4L,P0               ; A4 cross 32-bit word product

    xmpyu   MR,A5L,P1               ; A5 cross 32-bit word product

    ldd     XF2,S10                 ; ML * A0L

    fstd    P3,XF3                  ; ML * A1L

    xmpyu   ML,A4R,P2               ; A4 cross 32-bit word product

    xmpyu   ML,A5R,P3               ; A5 cross 32-bit word product

    ldd     XF4,S12                 ; MR * A2R

    fstd    P5,XF5                  ; MR * A3L

    xmpyu   MR,A6L,P4               ; A6 cross 32-bit word product

    xmpyu   MR,A7L,P5               ; A7 cross 32-bit word product

    ldd     XF6,S14                 ; ML * A2L

    fstd    P7,XF7                  ; ML * A3L

    xmpyu   ML,A6R,P6               ; A6 cross 32-bit word product

    xmpyu   ML,A7R,P7               ; A7 cross 32-bit word product

    fstd    P0,XF0                  ; MR * A4L

    ldd     XF1,S9                  ; MR * A1R

    shrpd   S4,T3,32,S4             ; A2 carry | cross product sum >> 32

    add     S1,S3,T1                ; A1 cross product sum

    fstd    P2,XF2                  ; ML * A4R

    ldd     XF3,S11                 ; ML * A1L

    add,dc  zero,zero,S1            ; A1 cross product sum carry

    depd,z  T1,31,32,S3             ; A1 cross product sum << 32

    fstd    P4,XF4                  ; MR * A6L

    ldd     XF5,S13                 ; MR * A3R

    shrpd   S1,T1,32,S1             ; A1 carry | cross product sum >> 32

    add     S5,S7,T3                ; A3 cross product sum

    fstd    P6,XF6                  ; ML * A6R

    ldd     XF7,S15                 ; ML * A3L

    add,dc  zero,zero,S5            ; A3 cross product sum carry

    depd,z  T3,31,32,S7             ; A3 cross product sum << 32

    shrpd   S5,T3,32,S5             ; A3 carry | cross product sum >> 32

    add     S2,S8,S8                ; M * A0 right doubleword, P0 doubleword

    add,dc  S0,S10,S10              ; M * A0 left doubleword

    add     S3,S9,S9                ; M * A1 right doubleword

    add,dc  S1,S11,S11              ; M * A1 left doubleword

    add     S6,S12,S12              ; M * A2 right doubleword

    ldd     24(pR),S3               ; Addend word 3

    fstd    P1,XF1                  ; MR * A5L

    add,dc  S4,S14,S14              ; M * A2 left doubleword

    xmpyu   MR,A5R,P1               ; A5 right 32-bit word product

    ldd     8(pR),S1                ; Addend word 1

    fstd    P3,XF3                  ; ML * A5R

    add     S7,S13,S13              ; M * A3 right doubleword

    xmpyu   ML,A5L,P3               ; A5 left 32-bit word product

    ldd     0(pR),S7                ; Addend word 0

    fstd    P5,XF5                  ; MR * A7L

    add,dc  S5,S15,S15              ; M * A3 left doubleword

    xmpyu   MR,A7R,P5               ; A7 right 32-bit word product

    ldd     16(pR),S5               ; Addend word 2

    fstd    P7,XF7                  ; ML * A7R

    add     S10,S9,S9               ; P1 doubleword

    xmpyu   ML,A7L,P7               ; A7 left 32-bit word products

    ldd     XF0,S0                  ; MR * A4L

    fstd    P1,XF9                  ; MR * A5R

    add,dc  S11,S12,S12             ; P2 doubleword

    xmpyu   MR,A4R,P0               ; A4 right 32-bit word product

    ldd     XF2,S2                  ; ML * A4R

    fstd    P3,XF11                 ; ML * A5L

    add,dc  S14,S13,S13             ; P3 doubleword

    xmpyu   ML,A4L,P2               ; A4 left 32-bit word product

    ldd     XF6,S6                  ; ML * A6R

    fstd    P5,XF13                 ; MR * A7R

    add,dc  zero,S15,T2             ; P4 partial doubleword

    xmpyu   MR,A6R,P4               ; A6 right 32-bit word product

    ldd     XF4,S4                  ; MR * A6L

    fstd    P7,XF15                 ; ML * A7L

    add     S7,S8,S8                ; R0 + P0, new R0 doubleword

    xmpyu   ML,A6L,P6               ; A6 left 32-bit word product

    fstd    P0,XF0                  ; MR * A4R

    ldd     XF7,S7                  ; ML * A7R

    add,dc  S1,S9,S9                ; c + R1 + P1, new R1 doubleword

    fstd    P2,XF2                  ; ML * A4L

    ldd     XF1,S1                  ; MR * A5L

    add,dc  S5,S12,S12              ; c + R2 + P2, new R2 doubleword

    fstd    P4,XF4                  ; MR * A6R

    ldd     XF5,S5                  ; MR * A7L

    add,dc  S3,S13,S13              ; c + R3 + P3, new R3 doubleword

    fstd    P6,XF6                  ; ML * A6L

    ldd     XF3,S3                  ; ML * A5R

    add,dc  zero,T2,T2              ; c + partial P4

    add     S0,S2,T1                ; A4 cross product sum

    std     S8,0(pR)                ; save R0

    add,dc  zero,zero,S0            ; A4 cross product sum carry

    depd,z  T1,31,32,S2             ; A4 cross product sum << 32

    std     S9,8(pR)                ; save R1

    shrpd   S0,T1,32,S0             ; A4 carry | cross product sum >> 32

    add     S4,S6,T3                ; A6 cross product sum

    std     S12,16(pR)              ; save R2

    add,dc  zero,zero,S4            ; A6 cross product sum carry

    depd,z  T3,31,32,S6             ; A6 cross product sum << 32

    std     S13,24(pR)              ; save R3

    shrpd   S4,T3,32,S4             ; A6 carry | cross product sum >> 32

    add     S1,S3,T1                ; A5 cross product sum

    ldd     XF0,S8                  ; MR * A4R

    add,dc  zero,zero,S1            ; A5 cross product sum carry

    depd,z  T1,31,32,S3             ; A5 cross product sum << 32

    ldd     XF2,S10                 ; ML * A4L

    ldd     XF9,S9                  ; MR * A5R

    shrpd   S1,T1,32,S1             ; A5 carry | cross product sum >> 32

    add     S5,S7,T3                ; A7 cross product sum

    ldd     XF4,S12                 ; MR * A6R

    ldd     XF11,S11                ; ML * A5L

    add,dc  zero,zero,S5            ; A7 cross product sum carry

    depd,z  T3,31,32,S7             ; A7 cross product sum << 32

    ldd     XF6,S14                 ; ML * A6L

    ldd     XF13,S13                ; MR * A7R

    shrpd   S5,T3,32,S5             ; A7 carry | cross product sum >> 32

    add     S2,S8,S8                ; M * A4 right doubleword

    ldd     XF15,S15                ; ML * A7L

    add,dc  S0,S10,S10              ; M * A4 left doubleword

    add     S3,S9,S9                ; M * A5 right doubleword

    add,dc  S1,S11,S11              ; M * A5 left doubleword

    add     S6,S12,S12              ; M * A6 right doubleword

    ldd     32(pR),S0               ; Addend word 4

    ldd     40(pR),S1               ; Addend word 5

    add,dc  S4,S14,S14              ; M * A6 left doubleword

    add     S7,S13,S13              ; M * A7 right doubleword

    ldd     48(pR),S2               ; Addend word 6

    ldd     56(pR),S3               ; Addend word 7

    add,dc  S5,S15,S15              ; M * A7 left doubleword

    add     S8,T2,S8                ; P4 doubleword

    ldd     64(pR),S4               ; Addend word 8

    ldd     SV5,s5                  ; restore s5

    add,dc  S10,S9,S9               ; P5 doubleword

    add,dc  S11,S12,S12             ; P6 doubleword

    ldd     SV6,s6                  ; restore s6

    ldd     SV7,s7                  ; restore s7

    add,dc  S14,S13,S13             ; P7 doubleword

    add,dc  zero,S15,S15            ; P8 doubleword

    add     S0,S8,S8                ; new R4 doubleword

    ldd     SV0,s0                  ; restore s0

    std     S8,32(pR)               ; save R4

    add,dc  S1,S9,S9                ; new R5 doubleword

    ldd     SV1,s1                  ; restore s1

    std     S9,40(pR)               ; save R5

    add,dc  S2,S12,S12              ; new R6 doubleword

    ldd     SV2,s2                  ; restore s2

    std     S12,48(pR)              ; save R6

    add,dc  S3,S13,S13              ; new R7 doubleword

    ldd     SV3,s3                  ; restore s3

    std     S13,56(pR)              ; save R7

    add,dc  S4,S15,S15              ; new R8 doubleword

    ldd     SV4,s4                  ; restore s4

    std     S15,64(pR)              ; save result[8]

    add,dc  zero,zero,v0            ; return carry from R8

    CMPIB,*= 0,v0,$L0               ; if no overflow, exit

    LDO     8(pR),pR

$FINAL1                             ; Final carry propagation

    LDD     64(pR),v0

    LDO     8(pR),pR

    ADDI    1,v0,v0

    CMPIB,*= 0,v0,$FINAL1           ; Keep looping if there is a carry.

    STD     v0,56(pR)

$L0

    bv      zero(rp)                ; -> caller

    ldo     -ST_SZ(sp),sp           ; pop stack

/*  ======================================================================  */

/*      end of module                                                       */

/*  ======================================================================  */

        bve (rp)

        .EXIT

        nop

                .PROCEND

                .SPACE         $TEXT$

                .SUBSPA        $CODE$

                .EXPORT        multacc512,ENTRY

        .end