1; mc88100 __gmpn_mul_1 -- Multiply a limb vector with a single limb and
2; store the product in a second limb vector.
3
4; Copyright 1992, 1994, 1995, 2000 Free Software Foundation, Inc.
5
6;  This file is part of the GNU MP Library.
7;
8;  The GNU MP Library is free software; you can redistribute it and/or modify
9;  it under the terms of either:
10;
11;    * the GNU Lesser General Public License as published by the Free
12;      Software Foundation; either version 3 of the License, or (at your
13;      option) any later version.
14;
15;  or
16;
17;    * the GNU General Public License as published by the Free Software
18;      Foundation; either version 2 of the License, or (at your option) any
19;      later version.
20;
21;  or both in parallel, as here.
22;
23;  The GNU MP Library is distributed in the hope that it will be useful, but
24;  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
25;  or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
26;  for more details.
27;
28;  You should have received copies of the GNU General Public License and the
29;  GNU Lesser General Public License along with the GNU MP Library.  If not,
30;  see https://www.gnu.org/licenses/.
31
32
33; INPUT PARAMETERS
34; res_ptr r2
35; s1_ptr  r3
36; size              r4
37; s2_limb r5
38
39; Common overhead is about 11 cycles/invocation.
40
41; The speed for S2_LIMB >= 0x10000 is approximately 21 cycles/limb.  (The
42; pipeline stalls 2 cycles due to WB contention.)
43
44; The speed for S2_LIMB < 0x10000 is approximately 16 cycles/limb.  (The
45; pipeline stalls 2 cycles due to WB contention and 1 cycle due to latency.)
46
47; To enhance speed:
48; 1. Unroll main loop 4-8 times.
49; 2. Schedule code to avoid WB contention.  It might be tempting to move the
50;    ld instruction in the loops down to save 2 cycles (less WB contention),
51;    but that looses because the ultimate value will be read from outside
52;    the allocated space.  But if we handle the ultimate multiplication in
53;    the tail, we can do this.
54; 3. Make the multiplication with less instructions.  I think the code for
55;    (S2_LIMB >= 0x10000) is not minimal.
56; With these techniques the (S2_LIMB >= 0x10000) case would run in 17 or
57; less cycles/limb; the (S2_LIMB < 0x10000) case would run in 11
58; cycles/limb.  (Assuming infinite unrolling.)
59
60          text
61          align      16
62          global     ___gmpn_mul_1
63___gmpn_mul_1:
64
65          ; Make S1_PTR and RES_PTR point at the end of their blocks
66          ; and negate SIZE.
67          lda        r3,r3[r4]
68          lda        r6,r2[r4]          ; RES_PTR in r6 since r2 is retval
69          subu       r4,r0,r4
70
71          addu.co    r2,r0,r0 ; r2 = cy = 0
72          ld         r9,r3[r4]
73          mask       r7,r5,0xffff       ; r7 = lo(S2_LIMB)
74          extu       r8,r5,16 ; r8 = hi(S2_LIMB)
75          bcnd.n     eq0,r8,Lsmall      ; jump if (hi(S2_LIMB) == 0)
76           subu      r6,r6,4
77
78; General code for any value of S2_LIMB.
79
80          ; Make a stack frame and save r25 and r26
81          subu       r31,r31,16
82          st.d       r25,r31,8
83
84          ; Enter the loop in the middle
85          br.n      L1
86          addu       r4,r4,1
87
88Loop:     ld         r9,r3[r4]
89          st         r26,r6[r4]
90; bcnd    ne0,r0,0            ; bubble
91          addu       r4,r4,1
92L1:       mul        r26,r9,r5          ; low word of product         mul_1     WB ld
93          mask       r12,r9,0xffff      ; r12 = lo(s1_limb) mask_1
94          mul        r11,r12,r7         ; r11 =  prod_0               mul_2     WB mask_1
95          mul        r10,r12,r8         ; r10 = prod_1a               mul_3
96          extu       r13,r9,16          ; r13 = hi(s1_limb) extu_1    WB mul_1
97          mul        r12,r13,r7         ; r12 = prod_1b               mul_4     WB extu_1
98          mul        r25,r13,r8         ; r25  = prod_2               mul_5     WB mul_2
99          extu       r11,r11,16         ; r11 = hi(prod_0)  extu_2    WB mul_3
100          addu       r10,r10,r11        ;                             addu_1    WB extu_2
101; bcnd    ne0,r0,0            ; bubble                      WB addu_1
102          addu.co    r10,r10,r12        ;                                       WB mul_4
103          mask.u     r10,r10,0xffff     ; move the 16 most significant bits...
104          addu.ci    r10,r10,r0         ; ...to the low half of the word...
105          rot        r10,r10,16         ; ...and put carry in pos 16.
106          addu.co    r26,r26,r2         ; add old carry limb
107          bcnd.n     ne0,r4,Loop
108           addu.ci r2,r25,r10 ; compute new carry limb
109
110          st         r26,r6[r4]
111          ld.d       r25,r31,8
112          jmp.n      r1
113           addu      r31,r31,16
114
115; Fast code for S2_LIMB < 0x10000
116Lsmall:
117          ; Enter the loop in the middle
118          br.n      SL1
119          addu       r4,r4,1
120
121SLoop:    ld         r9,r3[r4]          ;
122          st         r8,r6[r4]          ;
123          addu       r4,r4,1  ;
124SL1:      mul        r8,r9,r5 ; low word of product
125          mask       r12,r9,0xffff      ; r12 = lo(s1_limb)
126          extu       r13,r9,16          ; r13 = hi(s1_limb)
127          mul        r11,r12,r7         ; r11 =  prod_0
128          mul        r12,r13,r7         ; r12 = prod_1b
129          addu.cio r8,r8,r2   ; add old carry limb
130          extu       r10,r11,16         ; r11 = hi(prod_0)
131          addu       r10,r10,r12        ;
132          bcnd.n     ne0,r4,SLoop
133          extu       r2,r10,16          ; r2 = new carry limb
134
135          jmp.n      r1
136          st         r8,r6[r4]
137