aboutsummaryrefslogtreecommitdiffstats
path: root/lib/mpi
diff options
context:
space:
mode:
authorDmitry Torokhov <dmitry.torokhov@gmail.com>2012-03-19 20:02:01 -0400
committerDmitry Torokhov <dmitry.torokhov@gmail.com>2012-03-19 20:02:01 -0400
commit10ce3cc919f50c2043b41ca968b43c26a3672600 (patch)
treeea409366a5208aced495bc0516a08b81fd43222e /lib/mpi
parent24e3e5ae1e4c2a3a32f5b1f96b4e3fd721806acd (diff)
parent5c6a7a62c130afef3d61c1dee153012231ff5cd9 (diff)
Merge branch 'next' into for-linus
Diffstat (limited to 'lib/mpi')
-rw-r--r--lib/mpi/Makefile32
-rw-r--r--lib/mpi/generic_mpi-asm-defs.h4
-rw-r--r--lib/mpi/generic_mpih-add1.c61
-rw-r--r--lib/mpi/generic_mpih-lshift.c63
-rw-r--r--lib/mpi/generic_mpih-mul1.c57
-rw-r--r--lib/mpi/generic_mpih-mul2.c60
-rw-r--r--lib/mpi/generic_mpih-mul3.c61
-rw-r--r--lib/mpi/generic_mpih-rshift.c63
-rw-r--r--lib/mpi/generic_mpih-sub1.c60
-rw-r--r--lib/mpi/longlong.h1500
-rw-r--r--lib/mpi/mpi-add.c234
-rw-r--r--lib/mpi/mpi-bit.c217
-rw-r--r--lib/mpi/mpi-cmp.c68
-rw-r--r--lib/mpi/mpi-div.c338
-rw-r--r--lib/mpi/mpi-gcd.c59
-rw-r--r--lib/mpi/mpi-inline.c31
-rw-r--r--lib/mpi/mpi-inline.h122
-rw-r--r--lib/mpi/mpi-internal.h261
-rw-r--r--lib/mpi/mpi-inv.c187
-rw-r--r--lib/mpi/mpi-mpow.c134
-rw-r--r--lib/mpi/mpi-mul.c194
-rw-r--r--lib/mpi/mpi-pow.c323
-rw-r--r--lib/mpi/mpi-scan.c136
-rw-r--r--lib/mpi/mpicoder.c280
-rw-r--r--lib/mpi/mpih-cmp.c56
-rw-r--r--lib/mpi/mpih-div.c545
-rw-r--r--lib/mpi/mpih-mul.c527
-rw-r--r--lib/mpi/mpiutil.c211
28 files changed, 5884 insertions, 0 deletions
diff --git a/lib/mpi/Makefile b/lib/mpi/Makefile
new file mode 100644
index 000000000000..567d52e74d77
--- /dev/null
+++ b/lib/mpi/Makefile
@@ -0,0 +1,32 @@
1#
2# MPI multiprecision maths library (from gpg)
3#
4
5obj-$(CONFIG_MPILIB) = mpi.o
6
7mpi-y = \
8 generic_mpih-lshift.o \
9 generic_mpih-mul1.o \
10 generic_mpih-mul2.o \
11 generic_mpih-mul3.o \
12 generic_mpih-rshift.o \
13 generic_mpih-sub1.o \
14 generic_mpih-add1.o \
15 mpicoder.o \
16 mpi-bit.o \
17 mpih-cmp.o \
18 mpih-div.o \
19 mpih-mul.o \
20 mpi-pow.o \
21 mpiutil.o
22
23mpi-$(CONFIG_MPILIB_EXTRA) += \
24 mpi-add.o \
25 mpi-div.o \
26 mpi-cmp.o \
27 mpi-gcd.o \
28 mpi-inline.o \
29 mpi-inv.o \
30 mpi-mpow.o \
31 mpi-mul.o \
32 mpi-scan.o
diff --git a/lib/mpi/generic_mpi-asm-defs.h b/lib/mpi/generic_mpi-asm-defs.h
new file mode 100644
index 000000000000..047d1f5a7249
--- /dev/null
+++ b/lib/mpi/generic_mpi-asm-defs.h
@@ -0,0 +1,4 @@
1/* This file defines some basic constants for the MPI machinery. We
2 * need to define the types on a per-CPU basis, so it is done with
3 * this file here. */
4#define BYTES_PER_MPI_LIMB (SIZEOF_UNSIGNED_LONG)
diff --git a/lib/mpi/generic_mpih-add1.c b/lib/mpi/generic_mpih-add1.c
new file mode 100644
index 000000000000..c94c7dd344b3
--- /dev/null
+++ b/lib/mpi/generic_mpih-add1.c
@@ -0,0 +1,61 @@
1/* mpihelp-add_1.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1997, 1998,
3 * 2000 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31#include "longlong.h"
32
33mpi_limb_t
34mpihelp_add_n(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
35 mpi_ptr_t s2_ptr, mpi_size_t size)
36{
37 mpi_limb_t x, y, cy;
38 mpi_size_t j;
39
40 /* The loop counter and index J goes from -SIZE to -1. This way
41 the loop becomes faster. */
42 j = -size;
43
44 /* Offset the base pointers to compensate for the negative indices. */
45 s1_ptr -= j;
46 s2_ptr -= j;
47 res_ptr -= j;
48
49 cy = 0;
50 do {
51 y = s2_ptr[j];
52 x = s1_ptr[j];
53 y += cy; /* add previous carry to one addend */
54 cy = y < cy; /* get out carry from that addition */
55 y += x; /* add other addend */
56 cy += y < x; /* get out carry from that add, combine */
57 res_ptr[j] = y;
58 } while (++j);
59
60 return cy;
61}
diff --git a/lib/mpi/generic_mpih-lshift.c b/lib/mpi/generic_mpih-lshift.c
new file mode 100644
index 000000000000..86318927231a
--- /dev/null
+++ b/lib/mpi/generic_mpih-lshift.c
@@ -0,0 +1,63 @@
1/* mpihelp-lshift.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1998, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include "mpi-internal.h"
30
31/* Shift U (pointed to by UP and USIZE digits long) CNT bits to the left
32 * and store the USIZE least significant digits of the result at WP.
33 * Return the bits shifted out from the most significant digit.
34 *
35 * Argument constraints:
36 * 1. 0 < CNT < BITS_PER_MP_LIMB
37 * 2. If the result is to be written over the input, WP must be >= UP.
38 */
39
40mpi_limb_t
41mpihelp_lshift(mpi_ptr_t wp, mpi_ptr_t up, mpi_size_t usize, unsigned int cnt)
42{
43 mpi_limb_t high_limb, low_limb;
44 unsigned sh_1, sh_2;
45 mpi_size_t i;
46 mpi_limb_t retval;
47
48 sh_1 = cnt;
49 wp += 1;
50 sh_2 = BITS_PER_MPI_LIMB - sh_1;
51 i = usize - 1;
52 low_limb = up[i];
53 retval = low_limb >> sh_2;
54 high_limb = low_limb;
55 while (--i >= 0) {
56 low_limb = up[i];
57 wp[i] = (high_limb << sh_1) | (low_limb >> sh_2);
58 high_limb = low_limb;
59 }
60 wp[i] = high_limb << sh_1;
61
62 return retval;
63}
diff --git a/lib/mpi/generic_mpih-mul1.c b/lib/mpi/generic_mpih-mul1.c
new file mode 100644
index 000000000000..1668dfd9092c
--- /dev/null
+++ b/lib/mpi/generic_mpih-mul1.c
@@ -0,0 +1,57 @@
1/* mpihelp-mul_1.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1997, 1998, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include "mpi-internal.h"
30#include "longlong.h"
31
32mpi_limb_t
33mpihelp_mul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr, mpi_size_t s1_size,
34 mpi_limb_t s2_limb)
35{
36 mpi_limb_t cy_limb;
37 mpi_size_t j;
38 mpi_limb_t prod_high, prod_low;
39
40 /* The loop counter and index J goes from -S1_SIZE to -1. This way
41 * the loop becomes faster. */
42 j = -s1_size;
43
44 /* Offset the base pointers to compensate for the negative indices. */
45 s1_ptr -= j;
46 res_ptr -= j;
47
48 cy_limb = 0;
49 do {
50 umul_ppmm(prod_high, prod_low, s1_ptr[j], s2_limb);
51 prod_low += cy_limb;
52 cy_limb = (prod_low < cy_limb ? 1 : 0) + prod_high;
53 res_ptr[j] = prod_low;
54 } while (++j);
55
56 return cy_limb;
57}
diff --git a/lib/mpi/generic_mpih-mul2.c b/lib/mpi/generic_mpih-mul2.c
new file mode 100644
index 000000000000..8a7b29ee1740
--- /dev/null
+++ b/lib/mpi/generic_mpih-mul2.c
@@ -0,0 +1,60 @@
1/* mpihelp-mul_2.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1997, 1998, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include "mpi-internal.h"
30#include "longlong.h"
31
32mpi_limb_t
33mpihelp_addmul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
34 mpi_size_t s1_size, mpi_limb_t s2_limb)
35{
36 mpi_limb_t cy_limb;
37 mpi_size_t j;
38 mpi_limb_t prod_high, prod_low;
39 mpi_limb_t x;
40
41 /* The loop counter and index J goes from -SIZE to -1. This way
42 * the loop becomes faster. */
43 j = -s1_size;
44 res_ptr -= j;
45 s1_ptr -= j;
46
47 cy_limb = 0;
48 do {
49 umul_ppmm(prod_high, prod_low, s1_ptr[j], s2_limb);
50
51 prod_low += cy_limb;
52 cy_limb = (prod_low < cy_limb ? 1 : 0) + prod_high;
53
54 x = res_ptr[j];
55 prod_low = x + prod_low;
56 cy_limb += prod_low < x ? 1 : 0;
57 res_ptr[j] = prod_low;
58 } while (++j);
59 return cy_limb;
60}
diff --git a/lib/mpi/generic_mpih-mul3.c b/lib/mpi/generic_mpih-mul3.c
new file mode 100644
index 000000000000..f96df327be63
--- /dev/null
+++ b/lib/mpi/generic_mpih-mul3.c
@@ -0,0 +1,61 @@
1/* mpihelp-mul_3.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1997, 1998, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include "mpi-internal.h"
30#include "longlong.h"
31
32mpi_limb_t
33mpihelp_submul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
34 mpi_size_t s1_size, mpi_limb_t s2_limb)
35{
36 mpi_limb_t cy_limb;
37 mpi_size_t j;
38 mpi_limb_t prod_high, prod_low;
39 mpi_limb_t x;
40
41 /* The loop counter and index J goes from -SIZE to -1. This way
42 * the loop becomes faster. */
43 j = -s1_size;
44 res_ptr -= j;
45 s1_ptr -= j;
46
47 cy_limb = 0;
48 do {
49 umul_ppmm(prod_high, prod_low, s1_ptr[j], s2_limb);
50
51 prod_low += cy_limb;
52 cy_limb = (prod_low < cy_limb ? 1 : 0) + prod_high;
53
54 x = res_ptr[j];
55 prod_low = x - prod_low;
56 cy_limb += prod_low > x ? 1 : 0;
57 res_ptr[j] = prod_low;
58 } while (++j);
59
60 return cy_limb;
61}
diff --git a/lib/mpi/generic_mpih-rshift.c b/lib/mpi/generic_mpih-rshift.c
new file mode 100644
index 000000000000..ffa328818ca6
--- /dev/null
+++ b/lib/mpi/generic_mpih-rshift.c
@@ -0,0 +1,63 @@
1/* mpih-rshift.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1998, 1999,
3 * 2000, 2001 Free Software Foundation, Inc.
4 *
5 * This file is part of GNUPG
6 *
7 * GNUPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GNUPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31
32/* Shift U (pointed to by UP and USIZE limbs long) CNT bits to the right
33 * and store the USIZE least significant limbs of the result at WP.
34 * The bits shifted out to the right are returned.
35 *
36 * Argument constraints:
37 * 1. 0 < CNT < BITS_PER_MP_LIMB
38 * 2. If the result is to be written over the input, WP must be <= UP.
39 */
40
41mpi_limb_t
42mpihelp_rshift(mpi_ptr_t wp, mpi_ptr_t up, mpi_size_t usize, unsigned cnt)
43{
44 mpi_limb_t high_limb, low_limb;
45 unsigned sh_1, sh_2;
46 mpi_size_t i;
47 mpi_limb_t retval;
48
49 sh_1 = cnt;
50 wp -= 1;
51 sh_2 = BITS_PER_MPI_LIMB - sh_1;
52 high_limb = up[0];
53 retval = high_limb << sh_2;
54 low_limb = high_limb;
55 for (i = 1; i < usize; i++) {
56 high_limb = up[i];
57 wp[i] = (low_limb >> sh_1) | (high_limb << sh_2);
58 low_limb = high_limb;
59 }
60 wp[i] = low_limb >> sh_1;
61
62 return retval;
63}
diff --git a/lib/mpi/generic_mpih-sub1.c b/lib/mpi/generic_mpih-sub1.c
new file mode 100644
index 000000000000..5d98ab7d6853
--- /dev/null
+++ b/lib/mpi/generic_mpih-sub1.c
@@ -0,0 +1,60 @@
1/* mpihelp-add_2.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1997, 1998, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include "mpi-internal.h"
30#include "longlong.h"
31
32mpi_limb_t
33mpihelp_sub_n(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
34 mpi_ptr_t s2_ptr, mpi_size_t size)
35{
36 mpi_limb_t x, y, cy;
37 mpi_size_t j;
38
39 /* The loop counter and index J goes from -SIZE to -1. This way
40 the loop becomes faster. */
41 j = -size;
42
43 /* Offset the base pointers to compensate for the negative indices. */
44 s1_ptr -= j;
45 s2_ptr -= j;
46 res_ptr -= j;
47
48 cy = 0;
49 do {
50 y = s2_ptr[j];
51 x = s1_ptr[j];
52 y += cy; /* add previous carry to subtrahend */
53 cy = y < cy; /* get out carry from that addition */
54 y = x - y; /* main subtract */
55 cy += y > x; /* get out carry from the subtract, combine */
56 res_ptr[j] = y;
57 } while (++j);
58
59 return cy;
60}
diff --git a/lib/mpi/longlong.h b/lib/mpi/longlong.h
new file mode 100644
index 000000000000..29f98624ef93
--- /dev/null
+++ b/lib/mpi/longlong.h
@@ -0,0 +1,1500 @@
1/* longlong.h -- definitions for mixed size 32/64 bit arithmetic.
2 * Note: I added some stuff for use with gnupg
3 *
4 * Copyright (C) 1991, 1992, 1993, 1994, 1996, 1998,
5 * 2000, 2001, 2002, 2003 Free Software Foundation, Inc.
6 *
7 * This file is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU Library General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or (at your
10 * option) any later version.
11 *
12 * This file is distributed in the hope that it will be useful, but
13 * WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
15 * License for more details.
16 *
17 * You should have received a copy of the GNU Library General Public License
18 * along with this file; see the file COPYING.LIB. If not, write to
19 * the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
20 * MA 02111-1307, USA. */
21
22/* You have to define the following before including this file:
23 *
24 * UWtype -- An unsigned type, default type for operations (typically a "word")
25 * UHWtype -- An unsigned type, at least half the size of UWtype.
26 * UDWtype -- An unsigned type, at least twice as large a UWtype
27 * W_TYPE_SIZE -- size in bits of UWtype
28 *
29 * SItype, USItype -- Signed and unsigned 32 bit types.
30 * DItype, UDItype -- Signed and unsigned 64 bit types.
31 *
32 * On a 32 bit machine UWtype should typically be USItype;
33 * on a 64 bit machine, UWtype should typically be UDItype.
34*/
35
36#define __BITS4 (W_TYPE_SIZE / 4)
37#define __ll_B ((UWtype) 1 << (W_TYPE_SIZE / 2))
38#define __ll_lowpart(t) ((UWtype) (t) & (__ll_B - 1))
39#define __ll_highpart(t) ((UWtype) (t) >> (W_TYPE_SIZE / 2))
40
41/* This is used to make sure no undesirable sharing between different libraries
42 that use this file takes place. */
43#ifndef __MPN
44#define __MPN(x) __##x
45#endif
46
47/* Define auxiliary asm macros.
48 *
49 * 1) umul_ppmm(high_prod, low_prod, multipler, multiplicand) multiplies two
50 * UWtype integers MULTIPLER and MULTIPLICAND, and generates a two UWtype
51 * word product in HIGH_PROD and LOW_PROD.
52 *
53 * 2) __umulsidi3(a,b) multiplies two UWtype integers A and B, and returns a
54 * UDWtype product. This is just a variant of umul_ppmm.
55
56 * 3) udiv_qrnnd(quotient, remainder, high_numerator, low_numerator,
57 * denominator) divides a UDWtype, composed by the UWtype integers
58 * HIGH_NUMERATOR and LOW_NUMERATOR, by DENOMINATOR and places the quotient
59 * in QUOTIENT and the remainder in REMAINDER. HIGH_NUMERATOR must be less
60 * than DENOMINATOR for correct operation. If, in addition, the most
61 * significant bit of DENOMINATOR must be 1, then the pre-processor symbol
62 * UDIV_NEEDS_NORMALIZATION is defined to 1.
63 * 4) sdiv_qrnnd(quotient, remainder, high_numerator, low_numerator,
64 * denominator). Like udiv_qrnnd but the numbers are signed. The quotient
65 * is rounded towards 0.
66 *
67 * 5) count_leading_zeros(count, x) counts the number of zero-bits from the
68 * msb to the first non-zero bit in the UWtype X. This is the number of
69 * steps X needs to be shifted left to set the msb. Undefined for X == 0,
70 * unless the symbol COUNT_LEADING_ZEROS_0 is defined to some value.
71 *
72 * 6) count_trailing_zeros(count, x) like count_leading_zeros, but counts
73 * from the least significant end.
74 *
75 * 7) add_ssaaaa(high_sum, low_sum, high_addend_1, low_addend_1,
76 * high_addend_2, low_addend_2) adds two UWtype integers, composed by
77 * HIGH_ADDEND_1 and LOW_ADDEND_1, and HIGH_ADDEND_2 and LOW_ADDEND_2
78 * respectively. The result is placed in HIGH_SUM and LOW_SUM. Overflow
79 * (i.e. carry out) is not stored anywhere, and is lost.
80 *
81 * 8) sub_ddmmss(high_difference, low_difference, high_minuend, low_minuend,
82 * high_subtrahend, low_subtrahend) subtracts two two-word UWtype integers,
83 * composed by HIGH_MINUEND_1 and LOW_MINUEND_1, and HIGH_SUBTRAHEND_2 and
84 * LOW_SUBTRAHEND_2 respectively. The result is placed in HIGH_DIFFERENCE
85 * and LOW_DIFFERENCE. Overflow (i.e. carry out) is not stored anywhere,
86 * and is lost.
87 *
88 * If any of these macros are left undefined for a particular CPU,
89 * C macros are used. */
90
91/* The CPUs come in alphabetical order below.
92 *
93 * Please add support for more CPUs here, or improve the current support
94 * for the CPUs below! */
95
96#if defined(__GNUC__) && !defined(NO_ASM)
97
98/* We sometimes need to clobber "cc" with gcc2, but that would not be
99 understood by gcc1. Use cpp to avoid major code duplication. */
100#if __GNUC__ < 2
101#define __CLOBBER_CC
102#define __AND_CLOBBER_CC
103#else /* __GNUC__ >= 2 */
104#define __CLOBBER_CC : "cc"
105#define __AND_CLOBBER_CC , "cc"
106#endif /* __GNUC__ < 2 */
107
108/***************************************
109 ************** A29K *****************
110 ***************************************/
111#if (defined(__a29k__) || defined(_AM29K)) && W_TYPE_SIZE == 32
112#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
113 __asm__ ("add %1,%4,%5\n" \
114 "addc %0,%2,%3" \
115 : "=r" ((USItype)(sh)), \
116 "=&r" ((USItype)(sl)) \
117 : "%r" ((USItype)(ah)), \
118 "rI" ((USItype)(bh)), \
119 "%r" ((USItype)(al)), \
120 "rI" ((USItype)(bl)))
121#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
122 __asm__ ("sub %1,%4,%5\n" \
123 "subc %0,%2,%3" \
124 : "=r" ((USItype)(sh)), \
125 "=&r" ((USItype)(sl)) \
126 : "r" ((USItype)(ah)), \
127 "rI" ((USItype)(bh)), \
128 "r" ((USItype)(al)), \
129 "rI" ((USItype)(bl)))
130#define umul_ppmm(xh, xl, m0, m1) \
131do { \
132 USItype __m0 = (m0), __m1 = (m1); \
133 __asm__ ("multiplu %0,%1,%2" \
134 : "=r" ((USItype)(xl)) \
135 : "r" (__m0), \
136 "r" (__m1)); \
137 __asm__ ("multmu %0,%1,%2" \
138 : "=r" ((USItype)(xh)) \
139 : "r" (__m0), \
140 "r" (__m1)); \
141} while (0)
142#define udiv_qrnnd(q, r, n1, n0, d) \
143 __asm__ ("dividu %0,%3,%4" \
144 : "=r" ((USItype)(q)), \
145 "=q" ((USItype)(r)) \
146 : "1" ((USItype)(n1)), \
147 "r" ((USItype)(n0)), \
148 "r" ((USItype)(d)))
149
150#define count_leading_zeros(count, x) \
151 __asm__ ("clz %0,%1" \
152 : "=r" ((USItype)(count)) \
153 : "r" ((USItype)(x)))
154#define COUNT_LEADING_ZEROS_0 32
155#endif /* __a29k__ */
156
157#if defined(__alpha) && W_TYPE_SIZE == 64
158#define umul_ppmm(ph, pl, m0, m1) \
159do { \
160 UDItype __m0 = (m0), __m1 = (m1); \
161 __asm__ ("umulh %r1,%2,%0" \
162 : "=r" ((UDItype) ph) \
163 : "%rJ" (__m0), \
164 "rI" (__m1)); \
165 (pl) = __m0 * __m1; \
166 } while (0)
167#define UMUL_TIME 46
168#ifndef LONGLONG_STANDALONE
169#define udiv_qrnnd(q, r, n1, n0, d) \
170do { UDItype __r; \
171 (q) = __udiv_qrnnd(&__r, (n1), (n0), (d)); \
172 (r) = __r; \
173} while (0)
174extern UDItype __udiv_qrnnd();
175#define UDIV_TIME 220
176#endif /* LONGLONG_STANDALONE */
177#endif /* __alpha */
178
179/***************************************
180 ************** ARM ******************
181 ***************************************/
182#if defined(__arm__) && W_TYPE_SIZE == 32
183#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
184 __asm__ ("adds %1, %4, %5\n" \
185 "adc %0, %2, %3" \
186 : "=r" ((USItype)(sh)), \
187 "=&r" ((USItype)(sl)) \
188 : "%r" ((USItype)(ah)), \
189 "rI" ((USItype)(bh)), \
190 "%r" ((USItype)(al)), \
191 "rI" ((USItype)(bl)))
192#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
193 __asm__ ("subs %1, %4, %5\n" \
194 "sbc %0, %2, %3" \
195 : "=r" ((USItype)(sh)), \
196 "=&r" ((USItype)(sl)) \
197 : "r" ((USItype)(ah)), \
198 "rI" ((USItype)(bh)), \
199 "r" ((USItype)(al)), \
200 "rI" ((USItype)(bl)))
201#if defined __ARM_ARCH_2__ || defined __ARM_ARCH_3__
202#define umul_ppmm(xh, xl, a, b) \
203 __asm__ ("%@ Inlined umul_ppmm\n" \
204 "mov %|r0, %2, lsr #16 @ AAAA\n" \
205 "mov %|r2, %3, lsr #16 @ BBBB\n" \
206 "bic %|r1, %2, %|r0, lsl #16 @ aaaa\n" \
207 "bic %0, %3, %|r2, lsl #16 @ bbbb\n" \
208 "mul %1, %|r1, %|r2 @ aaaa * BBBB\n" \
209 "mul %|r2, %|r0, %|r2 @ AAAA * BBBB\n" \
210 "mul %|r1, %0, %|r1 @ aaaa * bbbb\n" \
211 "mul %0, %|r0, %0 @ AAAA * bbbb\n" \
212 "adds %|r0, %1, %0 @ central sum\n" \
213 "addcs %|r2, %|r2, #65536\n" \
214 "adds %1, %|r1, %|r0, lsl #16\n" \
215 "adc %0, %|r2, %|r0, lsr #16" \
216 : "=&r" ((USItype)(xh)), \
217 "=r" ((USItype)(xl)) \
218 : "r" ((USItype)(a)), \
219 "r" ((USItype)(b)) \
220 : "r0", "r1", "r2")
221#else
222#define umul_ppmm(xh, xl, a, b) \
223 __asm__ ("%@ Inlined umul_ppmm\n" \
224 "umull %r1, %r0, %r2, %r3" \
225 : "=&r" ((USItype)(xh)), \
226 "=r" ((USItype)(xl)) \
227 : "r" ((USItype)(a)), \
228 "r" ((USItype)(b)) \
229 : "r0", "r1")
230#endif
231#define UMUL_TIME 20
232#define UDIV_TIME 100
233#endif /* __arm__ */
234
235/***************************************
236 ************** CLIPPER **************
237 ***************************************/
238#if defined(__clipper__) && W_TYPE_SIZE == 32
239#define umul_ppmm(w1, w0, u, v) \
240 ({union {UDItype __ll; \
241 struct {USItype __l, __h; } __i; \
242 } __xx; \
243 __asm__ ("mulwux %2,%0" \
244 : "=r" (__xx.__ll) \
245 : "%0" ((USItype)(u)), \
246 "r" ((USItype)(v))); \
247 (w1) = __xx.__i.__h; (w0) = __xx.__i.__l; })
248#define smul_ppmm(w1, w0, u, v) \
249 ({union {DItype __ll; \
250 struct {SItype __l, __h; } __i; \
251 } __xx; \
252 __asm__ ("mulwx %2,%0" \
253 : "=r" (__xx.__ll) \
254 : "%0" ((SItype)(u)), \
255 "r" ((SItype)(v))); \
256 (w1) = __xx.__i.__h; (w0) = __xx.__i.__l; })
257#define __umulsidi3(u, v) \
258 ({UDItype __w; \
259 __asm__ ("mulwux %2,%0" \
260 : "=r" (__w) \
261 : "%0" ((USItype)(u)), \
262 "r" ((USItype)(v))); \
263 __w; })
264#endif /* __clipper__ */
265
266/***************************************
267 ************** GMICRO ***************
268 ***************************************/
269#if defined(__gmicro__) && W_TYPE_SIZE == 32
270#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
271 __asm__ ("add.w %5,%1\n" \
272 "addx %3,%0" \
273 : "=g" ((USItype)(sh)), \
274 "=&g" ((USItype)(sl)) \
275 : "%0" ((USItype)(ah)), \
276 "g" ((USItype)(bh)), \
277 "%1" ((USItype)(al)), \
278 "g" ((USItype)(bl)))
279#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
280 __asm__ ("sub.w %5,%1\n" \
281 "subx %3,%0" \
282 : "=g" ((USItype)(sh)), \
283 "=&g" ((USItype)(sl)) \
284 : "0" ((USItype)(ah)), \
285 "g" ((USItype)(bh)), \
286 "1" ((USItype)(al)), \
287 "g" ((USItype)(bl)))
288#define umul_ppmm(ph, pl, m0, m1) \
289 __asm__ ("mulx %3,%0,%1" \
290 : "=g" ((USItype)(ph)), \
291 "=r" ((USItype)(pl)) \
292 : "%0" ((USItype)(m0)), \
293 "g" ((USItype)(m1)))
294#define udiv_qrnnd(q, r, nh, nl, d) \
295 __asm__ ("divx %4,%0,%1" \
296 : "=g" ((USItype)(q)), \
297 "=r" ((USItype)(r)) \
298 : "1" ((USItype)(nh)), \
299 "0" ((USItype)(nl)), \
300 "g" ((USItype)(d)))
301#define count_leading_zeros(count, x) \
302 __asm__ ("bsch/1 %1,%0" \
303 : "=g" (count) \
304 : "g" ((USItype)(x)), \
305 "0" ((USItype)0))
306#endif
307
308/***************************************
309 ************** HPPA *****************
310 ***************************************/
311#if defined(__hppa) && W_TYPE_SIZE == 32
312#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
313 __asm__ ("add %4,%5,%1\n" \
314 "addc %2,%3,%0" \
315 : "=r" ((USItype)(sh)), \
316 "=&r" ((USItype)(sl)) \
317 : "%rM" ((USItype)(ah)), \
318 "rM" ((USItype)(bh)), \
319 "%rM" ((USItype)(al)), \
320 "rM" ((USItype)(bl)))
321#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
322 __asm__ ("sub %4,%5,%1\n" \
323 "subb %2,%3,%0" \
324 : "=r" ((USItype)(sh)), \
325 "=&r" ((USItype)(sl)) \
326 : "rM" ((USItype)(ah)), \
327 "rM" ((USItype)(bh)), \
328 "rM" ((USItype)(al)), \
329 "rM" ((USItype)(bl)))
330#if defined(_PA_RISC1_1)
331#define umul_ppmm(wh, wl, u, v) \
332do { \
333 union {UDItype __ll; \
334 struct {USItype __h, __l; } __i; \
335 } __xx; \
336 __asm__ ("xmpyu %1,%2,%0" \
337 : "=*f" (__xx.__ll) \
338 : "*f" ((USItype)(u)), \
339 "*f" ((USItype)(v))); \
340 (wh) = __xx.__i.__h; \
341 (wl) = __xx.__i.__l; \
342} while (0)
343#define UMUL_TIME 8
344#define UDIV_TIME 60
345#else
346#define UMUL_TIME 40
347#define UDIV_TIME 80
348#endif
349#ifndef LONGLONG_STANDALONE
350#define udiv_qrnnd(q, r, n1, n0, d) \
351do { USItype __r; \
352 (q) = __udiv_qrnnd(&__r, (n1), (n0), (d)); \
353 (r) = __r; \
354} while (0)
355extern USItype __udiv_qrnnd();
356#endif /* LONGLONG_STANDALONE */
357#define count_leading_zeros(count, x) \
358do { \
359 USItype __tmp; \
360 __asm__ ( \
361 "ldi 1,%0\n" \
362 "extru,= %1,15,16,%%r0 ; Bits 31..16 zero?\n" \
363 "extru,tr %1,15,16,%1 ; No. Shift down, skip add.\n" \
364 "ldo 16(%0),%0 ; Yes. Perform add.\n" \
365 "extru,= %1,23,8,%%r0 ; Bits 15..8 zero?\n" \
366 "extru,tr %1,23,8,%1 ; No. Shift down, skip add.\n" \
367 "ldo 8(%0),%0 ; Yes. Perform add.\n" \
368 "extru,= %1,27,4,%%r0 ; Bits 7..4 zero?\n" \
369 "extru,tr %1,27,4,%1 ; No. Shift down, skip add.\n" \
370 "ldo 4(%0),%0 ; Yes. Perform add.\n" \
371 "extru,= %1,29,2,%%r0 ; Bits 3..2 zero?\n" \
372 "extru,tr %1,29,2,%1 ; No. Shift down, skip add.\n" \
373 "ldo 2(%0),%0 ; Yes. Perform add.\n" \
374 "extru %1,30,1,%1 ; Extract bit 1.\n" \
375 "sub %0,%1,%0 ; Subtract it. " \
376 : "=r" (count), "=r" (__tmp) : "1" (x)); \
377} while (0)
378#endif /* hppa */
379
380/***************************************
381 ************** I370 *****************
382 ***************************************/
383#if (defined(__i370__) || defined(__mvs__)) && W_TYPE_SIZE == 32
384#define umul_ppmm(xh, xl, m0, m1) \
385do { \
386 union {UDItype __ll; \
387 struct {USItype __h, __l; } __i; \
388 } __xx; \
389 USItype __m0 = (m0), __m1 = (m1); \
390 __asm__ ("mr %0,%3" \
391 : "=r" (__xx.__i.__h), \
392 "=r" (__xx.__i.__l) \
393 : "%1" (__m0), \
394 "r" (__m1)); \
395 (xh) = __xx.__i.__h; (xl) = __xx.__i.__l; \
396 (xh) += ((((SItype) __m0 >> 31) & __m1) \
397 + (((SItype) __m1 >> 31) & __m0)); \
398} while (0)
399#define smul_ppmm(xh, xl, m0, m1) \
400do { \
401 union {DItype __ll; \
402 struct {USItype __h, __l; } __i; \
403 } __xx; \
404 __asm__ ("mr %0,%3" \
405 : "=r" (__xx.__i.__h), \
406 "=r" (__xx.__i.__l) \
407 : "%1" (m0), \
408 "r" (m1)); \
409 (xh) = __xx.__i.__h; (xl) = __xx.__i.__l; \
410} while (0)
411#define sdiv_qrnnd(q, r, n1, n0, d) \
412do { \
413 union {DItype __ll; \
414 struct {USItype __h, __l; } __i; \
415 } __xx; \
416 __xx.__i.__h = n1; __xx.__i.__l = n0; \
417 __asm__ ("dr %0,%2" \
418 : "=r" (__xx.__ll) \
419 : "0" (__xx.__ll), "r" (d)); \
420 (q) = __xx.__i.__l; (r) = __xx.__i.__h; \
421} while (0)
422#endif
423
424/***************************************
425 ************** I386 *****************
426 ***************************************/
427#undef __i386__
428#if (defined(__i386__) || defined(__i486__)) && W_TYPE_SIZE == 32
429#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
430 __asm__ ("addl %5,%1\n" \
431 "adcl %3,%0" \
432 : "=r" ((USItype)(sh)), \
433 "=&r" ((USItype)(sl)) \
434 : "%0" ((USItype)(ah)), \
435 "g" ((USItype)(bh)), \
436 "%1" ((USItype)(al)), \
437 "g" ((USItype)(bl)))
438#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
439 __asm__ ("subl %5,%1\n" \
440 "sbbl %3,%0" \
441 : "=r" ((USItype)(sh)), \
442 "=&r" ((USItype)(sl)) \
443 : "0" ((USItype)(ah)), \
444 "g" ((USItype)(bh)), \
445 "1" ((USItype)(al)), \
446 "g" ((USItype)(bl)))
447#define umul_ppmm(w1, w0, u, v) \
448 __asm__ ("mull %3" \
449 : "=a" ((USItype)(w0)), \
450 "=d" ((USItype)(w1)) \
451 : "%0" ((USItype)(u)), \
452 "rm" ((USItype)(v)))
453#define udiv_qrnnd(q, r, n1, n0, d) \
454 __asm__ ("divl %4" \
455 : "=a" ((USItype)(q)), \
456 "=d" ((USItype)(r)) \
457 : "0" ((USItype)(n0)), \
458 "1" ((USItype)(n1)), \
459 "rm" ((USItype)(d)))
460#define count_leading_zeros(count, x) \
461do { \
462 USItype __cbtmp; \
463 __asm__ ("bsrl %1,%0" \
464 : "=r" (__cbtmp) : "rm" ((USItype)(x))); \
465 (count) = __cbtmp ^ 31; \
466} while (0)
467#define count_trailing_zeros(count, x) \
468 __asm__ ("bsfl %1,%0" : "=r" (count) : "rm" ((USItype)(x)))
469#ifndef UMUL_TIME
470#define UMUL_TIME 40
471#endif
472#ifndef UDIV_TIME
473#define UDIV_TIME 40
474#endif
475#endif /* 80x86 */
476
477/***************************************
478 ************** I860 *****************
479 ***************************************/
480#if defined(__i860__) && W_TYPE_SIZE == 32
481#define rshift_rhlc(r, h, l, c) \
482 __asm__ ("shr %3,r0,r0\n" \
483 "shrd %1,%2,%0" \
484 "=r" (r) : "r" (h), "r" (l), "rn" (c))
485#endif /* i860 */
486
487/***************************************
488 ************** I960 *****************
489 ***************************************/
490#if defined(__i960__) && W_TYPE_SIZE == 32
491#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
492 __asm__ ("cmpo 1,0\n" \
493 "addc %5,%4,%1\n" \
494 "addc %3,%2,%0" \
495 : "=r" ((USItype)(sh)), \
496 "=&r" ((USItype)(sl)) \
497 : "%dI" ((USItype)(ah)), \
498 "dI" ((USItype)(bh)), \
499 "%dI" ((USItype)(al)), \
500 "dI" ((USItype)(bl)))
501#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
502 __asm__ ("cmpo 0,0\n" \
503 "subc %5,%4,%1\n" \
504 "subc %3,%2,%0" \
505 : "=r" ((USItype)(sh)), \
506 "=&r" ((USItype)(sl)) \
507 : "dI" ((USItype)(ah)), \
508 "dI" ((USItype)(bh)), \
509 "dI" ((USItype)(al)), \
510 "dI" ((USItype)(bl)))
511#define umul_ppmm(w1, w0, u, v) \
512 ({union {UDItype __ll; \
513 struct {USItype __l, __h; } __i; \
514 } __xx; \
515 __asm__ ("emul %2,%1,%0" \
516 : "=d" (__xx.__ll) \
517 : "%dI" ((USItype)(u)), \
518 "dI" ((USItype)(v))); \
519 (w1) = __xx.__i.__h; (w0) = __xx.__i.__l; })
520#define __umulsidi3(u, v) \
521 ({UDItype __w; \
522 __asm__ ("emul %2,%1,%0" \
523 : "=d" (__w) \
524 : "%dI" ((USItype)(u)), \
525 "dI" ((USItype)(v))); \
526 __w; })
527#define udiv_qrnnd(q, r, nh, nl, d) \
528do { \
529 union {UDItype __ll; \
530 struct {USItype __l, __h; } __i; \
531 } __nn; \
532 __nn.__i.__h = (nh); __nn.__i.__l = (nl); \
533 __asm__ ("ediv %d,%n,%0" \
534 : "=d" (__rq.__ll) \
535 : "dI" (__nn.__ll), \
536 "dI" ((USItype)(d))); \
537 (r) = __rq.__i.__l; (q) = __rq.__i.__h; \
538} while (0)
539#define count_leading_zeros(count, x) \
540do { \
541 USItype __cbtmp; \
542 __asm__ ("scanbit %1,%0" \
543 : "=r" (__cbtmp) \
544 : "r" ((USItype)(x))); \
545 (count) = __cbtmp ^ 31; \
546} while (0)
547#define COUNT_LEADING_ZEROS_0 (-32) /* sic */
548#if defined(__i960mx) /* what is the proper symbol to test??? */
549#define rshift_rhlc(r, h, l, c) \
550do { \
551 union {UDItype __ll; \
552 struct {USItype __l, __h; } __i; \
553 } __nn; \
554 __nn.__i.__h = (h); __nn.__i.__l = (l); \
555 __asm__ ("shre %2,%1,%0" \
556 : "=d" (r) : "dI" (__nn.__ll), "dI" (c)); \
557}
558#endif /* i960mx */
559#endif /* i960 */
560
561/***************************************
562 ************** 68000 ****************
563 ***************************************/
564#if (defined(__mc68000__) || defined(__mc68020__) || defined(__NeXT__) || defined(mc68020)) && W_TYPE_SIZE == 32
565#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
566 __asm__ ("add%.l %5,%1\n" \
567 "addx%.l %3,%0" \
568 : "=d" ((USItype)(sh)), \
569 "=&d" ((USItype)(sl)) \
570 : "%0" ((USItype)(ah)), \
571 "d" ((USItype)(bh)), \
572 "%1" ((USItype)(al)), \
573 "g" ((USItype)(bl)))
574#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
575 __asm__ ("sub%.l %5,%1\n" \
576 "subx%.l %3,%0" \
577 : "=d" ((USItype)(sh)), \
578 "=&d" ((USItype)(sl)) \
579 : "0" ((USItype)(ah)), \
580 "d" ((USItype)(bh)), \
581 "1" ((USItype)(al)), \
582 "g" ((USItype)(bl)))
583#if (defined(__mc68020__) || defined(__NeXT__) || defined(mc68020))
584#define umul_ppmm(w1, w0, u, v) \
585 __asm__ ("mulu%.l %3,%1:%0" \
586 : "=d" ((USItype)(w0)), \
587 "=d" ((USItype)(w1)) \
588 : "%0" ((USItype)(u)), \
589 "dmi" ((USItype)(v)))
590#define UMUL_TIME 45
591#define udiv_qrnnd(q, r, n1, n0, d) \
592 __asm__ ("divu%.l %4,%1:%0" \
593 : "=d" ((USItype)(q)), \
594 "=d" ((USItype)(r)) \
595 : "0" ((USItype)(n0)), \
596 "1" ((USItype)(n1)), \
597 "dmi" ((USItype)(d)))
598#define UDIV_TIME 90
599#define sdiv_qrnnd(q, r, n1, n0, d) \
600 __asm__ ("divs%.l %4,%1:%0" \
601 : "=d" ((USItype)(q)), \
602 "=d" ((USItype)(r)) \
603 : "0" ((USItype)(n0)), \
604 "1" ((USItype)(n1)), \
605 "dmi" ((USItype)(d)))
606#define count_leading_zeros(count, x) \
607 __asm__ ("bfffo %1{%b2:%b2},%0" \
608 : "=d" ((USItype)(count)) \
609 : "od" ((USItype)(x)), "n" (0))
610#define COUNT_LEADING_ZEROS_0 32
611#else /* not mc68020 */
612#define umul_ppmm(xh, xl, a, b) \
613do { USItype __umul_tmp1, __umul_tmp2; \
614 __asm__ ("| Inlined umul_ppmm\n" \
615 "move%.l %5,%3\n" \
616 "move%.l %2,%0\n" \
617 "move%.w %3,%1\n" \
618 "swap %3\n" \
619 "swap %0\n" \
620 "mulu %2,%1\n" \
621 "mulu %3,%0\n" \
622 "mulu %2,%3\n" \
623 "swap %2\n" \
624 "mulu %5,%2\n" \
625 "add%.l %3,%2\n" \
626 "jcc 1f\n" \
627 "add%.l %#0x10000,%0\n" \
628 "1: move%.l %2,%3\n" \
629 "clr%.w %2\n" \
630 "swap %2\n" \
631 "swap %3\n" \
632 "clr%.w %3\n" \
633 "add%.l %3,%1\n" \
634 "addx%.l %2,%0\n" \
635 "| End inlined umul_ppmm" \
636 : "=&d" ((USItype)(xh)), "=&d" ((USItype)(xl)), \
637 "=d" (__umul_tmp1), "=&d" (__umul_tmp2) \
638 : "%2" ((USItype)(a)), "d" ((USItype)(b))); \
639} while (0)
640#define UMUL_TIME 100
641#define UDIV_TIME 400
642#endif /* not mc68020 */
643#endif /* mc68000 */
644
645/***************************************
646 ************** 88000 ****************
647 ***************************************/
648#if defined(__m88000__) && W_TYPE_SIZE == 32
649#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
650 __asm__ ("addu.co %1,%r4,%r5\n" \
651 "addu.ci %0,%r2,%r3" \
652 : "=r" ((USItype)(sh)), \
653 "=&r" ((USItype)(sl)) \
654 : "%rJ" ((USItype)(ah)), \
655 "rJ" ((USItype)(bh)), \
656 "%rJ" ((USItype)(al)), \
657 "rJ" ((USItype)(bl)))
658#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
659 __asm__ ("subu.co %1,%r4,%r5\n" \
660 "subu.ci %0,%r2,%r3" \
661 : "=r" ((USItype)(sh)), \
662 "=&r" ((USItype)(sl)) \
663 : "rJ" ((USItype)(ah)), \
664 "rJ" ((USItype)(bh)), \
665 "rJ" ((USItype)(al)), \
666 "rJ" ((USItype)(bl)))
667#define count_leading_zeros(count, x) \
668do { \
669 USItype __cbtmp; \
670 __asm__ ("ff1 %0,%1" \
671 : "=r" (__cbtmp) \
672 : "r" ((USItype)(x))); \
673 (count) = __cbtmp ^ 31; \
674} while (0)
675#define COUNT_LEADING_ZEROS_0 63 /* sic */
676#if defined(__m88110__)
677#define umul_ppmm(wh, wl, u, v) \
678do { \
679 union {UDItype __ll; \
680 struct {USItype __h, __l; } __i; \
681 } __x; \
682 __asm__ ("mulu.d %0,%1,%2" : "=r" (__x.__ll) : "r" (u), "r" (v)); \
683 (wh) = __x.__i.__h; \
684 (wl) = __x.__i.__l; \
685} while (0)
686#define udiv_qrnnd(q, r, n1, n0, d) \
687 ({union {UDItype __ll; \
688 struct {USItype __h, __l; } __i; \
689 } __x, __q; \
690 __x.__i.__h = (n1); __x.__i.__l = (n0); \
691 __asm__ ("divu.d %0,%1,%2" \
692 : "=r" (__q.__ll) : "r" (__x.__ll), "r" (d)); \
693 (r) = (n0) - __q.__l * (d); (q) = __q.__l; })
694#define UMUL_TIME 5
695#define UDIV_TIME 25
696#else
697#define UMUL_TIME 17
698#define UDIV_TIME 150
699#endif /* __m88110__ */
700#endif /* __m88000__ */
701
702/***************************************
703 ************** MIPS *****************
704 ***************************************/
705#if defined(__mips__) && W_TYPE_SIZE == 32
706#if __GNUC__ > 2 || __GNUC_MINOR__ >= 7
707#define umul_ppmm(w1, w0, u, v) \
708 __asm__ ("multu %2,%3" \
709 : "=l" ((USItype)(w0)), \
710 "=h" ((USItype)(w1)) \
711 : "d" ((USItype)(u)), \
712 "d" ((USItype)(v)))
713#else
714#define umul_ppmm(w1, w0, u, v) \
715 __asm__ ("multu %2,%3\n" \
716 "mflo %0\n" \
717 "mfhi %1" \
718 : "=d" ((USItype)(w0)), \
719 "=d" ((USItype)(w1)) \
720 : "d" ((USItype)(u)), \
721 "d" ((USItype)(v)))
722#endif
723#define UMUL_TIME 10
724#define UDIV_TIME 100
725#endif /* __mips__ */
726
727/***************************************
728 ************** MIPS/64 **************
729 ***************************************/
730#if (defined(__mips) && __mips >= 3) && W_TYPE_SIZE == 64
731#if __GNUC__ > 2 || __GNUC_MINOR__ >= 7
732#define umul_ppmm(w1, w0, u, v) \
733 __asm__ ("dmultu %2,%3" \
734 : "=l" ((UDItype)(w0)), \
735 "=h" ((UDItype)(w1)) \
736 : "d" ((UDItype)(u)), \
737 "d" ((UDItype)(v)))
738#else
739#define umul_ppmm(w1, w0, u, v) \
740 __asm__ ("dmultu %2,%3\n" \
741 "mflo %0\n" \
742 "mfhi %1" \
743 : "=d" ((UDItype)(w0)), \
744 "=d" ((UDItype)(w1)) \
745 : "d" ((UDItype)(u)), \
746 "d" ((UDItype)(v)))
747#endif
748#define UMUL_TIME 20
749#define UDIV_TIME 140
750#endif /* __mips__ */
751
752/***************************************
753 ************** 32000 ****************
754 ***************************************/
755#if defined(__ns32000__) && W_TYPE_SIZE == 32
756#define umul_ppmm(w1, w0, u, v) \
757 ({union {UDItype __ll; \
758 struct {USItype __l, __h; } __i; \
759 } __xx; \
760 __asm__ ("meid %2,%0" \
761 : "=g" (__xx.__ll) \
762 : "%0" ((USItype)(u)), \
763 "g" ((USItype)(v))); \
764 (w1) = __xx.__i.__h; (w0) = __xx.__i.__l; })
765#define __umulsidi3(u, v) \
766 ({UDItype __w; \
767 __asm__ ("meid %2,%0" \
768 : "=g" (__w) \
769 : "%0" ((USItype)(u)), \
770 "g" ((USItype)(v))); \
771 __w; })
772#define udiv_qrnnd(q, r, n1, n0, d) \
773 ({union {UDItype __ll; \
774 struct {USItype __l, __h; } __i; \
775 } __xx; \
776 __xx.__i.__h = (n1); __xx.__i.__l = (n0); \
777 __asm__ ("deid %2,%0" \
778 : "=g" (__xx.__ll) \
779 : "0" (__xx.__ll), \
780 "g" ((USItype)(d))); \
781 (r) = __xx.__i.__l; (q) = __xx.__i.__h; })
782#define count_trailing_zeros(count, x) \
783do { \
784 __asm__("ffsd %2,%0" \
785 : "=r"((USItype) (count)) \
786 : "0"((USItype) 0), "r"((USItype) (x))); \
787 } while (0)
788#endif /* __ns32000__ */
789
790/***************************************
791 ************** PPC ******************
792 ***************************************/
793#if (defined(_ARCH_PPC) || defined(_IBMR2)) && W_TYPE_SIZE == 32
794#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
795do { \
796 if (__builtin_constant_p(bh) && (bh) == 0) \
797 __asm__ ("{a%I4|add%I4c} %1,%3,%4\n\t{aze|addze} %0,%2" \
798 : "=r" ((USItype)(sh)), \
799 "=&r" ((USItype)(sl)) \
800 : "%r" ((USItype)(ah)), \
801 "%r" ((USItype)(al)), \
802 "rI" ((USItype)(bl))); \
803 else if (__builtin_constant_p(bh) && (bh) == ~(USItype) 0) \
804 __asm__ ("{a%I4|add%I4c} %1,%3,%4\n\t{ame|addme} %0,%2" \
805 : "=r" ((USItype)(sh)), \
806 "=&r" ((USItype)(sl)) \
807 : "%r" ((USItype)(ah)), \
808 "%r" ((USItype)(al)), \
809 "rI" ((USItype)(bl))); \
810 else \
811 __asm__ ("{a%I5|add%I5c} %1,%4,%5\n\t{ae|adde} %0,%2,%3" \
812 : "=r" ((USItype)(sh)), \
813 "=&r" ((USItype)(sl)) \
814 : "%r" ((USItype)(ah)), \
815 "r" ((USItype)(bh)), \
816 "%r" ((USItype)(al)), \
817 "rI" ((USItype)(bl))); \
818} while (0)
819#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
820do { \
821 if (__builtin_constant_p(ah) && (ah) == 0) \
822 __asm__ ("{sf%I3|subf%I3c} %1,%4,%3\n\t{sfze|subfze} %0,%2" \
823 : "=r" ((USItype)(sh)), \
824 "=&r" ((USItype)(sl)) \
825 : "r" ((USItype)(bh)), \
826 "rI" ((USItype)(al)), \
827 "r" ((USItype)(bl))); \
828 else if (__builtin_constant_p(ah) && (ah) == ~(USItype) 0) \
829 __asm__ ("{sf%I3|subf%I3c} %1,%4,%3\n\t{sfme|subfme} %0,%2" \
830 : "=r" ((USItype)(sh)), \
831 "=&r" ((USItype)(sl)) \
832 : "r" ((USItype)(bh)), \
833 "rI" ((USItype)(al)), \
834 "r" ((USItype)(bl))); \
835 else if (__builtin_constant_p(bh) && (bh) == 0) \
836 __asm__ ("{sf%I3|subf%I3c} %1,%4,%3\n\t{ame|addme} %0,%2" \
837 : "=r" ((USItype)(sh)), \
838 "=&r" ((USItype)(sl)) \
839 : "r" ((USItype)(ah)), \
840 "rI" ((USItype)(al)), \
841 "r" ((USItype)(bl))); \
842 else if (__builtin_constant_p(bh) && (bh) == ~(USItype) 0) \
843 __asm__ ("{sf%I3|subf%I3c} %1,%4,%3\n\t{aze|addze} %0,%2" \
844 : "=r" ((USItype)(sh)), \
845 "=&r" ((USItype)(sl)) \
846 : "r" ((USItype)(ah)), \
847 "rI" ((USItype)(al)), \
848 "r" ((USItype)(bl))); \
849 else \
850 __asm__ ("{sf%I4|subf%I4c} %1,%5,%4\n\t{sfe|subfe} %0,%3,%2" \
851 : "=r" ((USItype)(sh)), \
852 "=&r" ((USItype)(sl)) \
853 : "r" ((USItype)(ah)), \
854 "r" ((USItype)(bh)), \
855 "rI" ((USItype)(al)), \
856 "r" ((USItype)(bl))); \
857} while (0)
858#define count_leading_zeros(count, x) \
859 __asm__ ("{cntlz|cntlzw} %0,%1" \
860 : "=r" ((USItype)(count)) \
861 : "r" ((USItype)(x)))
862#define COUNT_LEADING_ZEROS_0 32
863#if defined(_ARCH_PPC)
864#define umul_ppmm(ph, pl, m0, m1) \
865do { \
866 USItype __m0 = (m0), __m1 = (m1); \
867 __asm__ ("mulhwu %0,%1,%2" \
868 : "=r" ((USItype) ph) \
869 : "%r" (__m0), \
870 "r" (__m1)); \
871 (pl) = __m0 * __m1; \
872} while (0)
873#define UMUL_TIME 15
874#define smul_ppmm(ph, pl, m0, m1) \
875do { \
876 SItype __m0 = (m0), __m1 = (m1); \
877 __asm__ ("mulhw %0,%1,%2" \
878 : "=r" ((SItype) ph) \
879 : "%r" (__m0), \
880 "r" (__m1)); \
881 (pl) = __m0 * __m1; \
882} while (0)
883#define SMUL_TIME 14
884#define UDIV_TIME 120
885#else
886#define umul_ppmm(xh, xl, m0, m1) \
887do { \
888 USItype __m0 = (m0), __m1 = (m1); \
889 __asm__ ("mul %0,%2,%3" \
890 : "=r" ((USItype)(xh)), \
891 "=q" ((USItype)(xl)) \
892 : "r" (__m0), \
893 "r" (__m1)); \
894 (xh) += ((((SItype) __m0 >> 31) & __m1) \
895 + (((SItype) __m1 >> 31) & __m0)); \
896} while (0)
897#define UMUL_TIME 8
898#define smul_ppmm(xh, xl, m0, m1) \
899 __asm__ ("mul %0,%2,%3" \
900 : "=r" ((SItype)(xh)), \
901 "=q" ((SItype)(xl)) \
902 : "r" (m0), \
903 "r" (m1))
904#define SMUL_TIME 4
905#define sdiv_qrnnd(q, r, nh, nl, d) \
906 __asm__ ("div %0,%2,%4" \
907 : "=r" ((SItype)(q)), "=q" ((SItype)(r)) \
908 : "r" ((SItype)(nh)), "1" ((SItype)(nl)), "r" ((SItype)(d)))
909#define UDIV_TIME 100
910#endif
911#endif /* Power architecture variants. */
912
913/***************************************
914 ************** PYR ******************
915 ***************************************/
916#if defined(__pyr__) && W_TYPE_SIZE == 32
917#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
918 __asm__ ("addw %5,%1\n" \
919 "addwc %3,%0" \
920 : "=r" ((USItype)(sh)), \
921 "=&r" ((USItype)(sl)) \
922 : "%0" ((USItype)(ah)), \
923 "g" ((USItype)(bh)), \
924 "%1" ((USItype)(al)), \
925 "g" ((USItype)(bl)))
926#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
927 __asm__ ("subw %5,%1\n" \
928 "subwb %3,%0" \
929 : "=r" ((USItype)(sh)), \
930 "=&r" ((USItype)(sl)) \
931 : "0" ((USItype)(ah)), \
932 "g" ((USItype)(bh)), \
933 "1" ((USItype)(al)), \
934 "g" ((USItype)(bl)))
935 /* This insn works on Pyramids with AP, XP, or MI CPUs, but not with SP. */
936#define umul_ppmm(w1, w0, u, v) \
937 ({union {UDItype __ll; \
938 struct {USItype __h, __l; } __i; \
939 } __xx; \
940 __asm__ ("movw %1,%R0\n" \
941 "uemul %2,%0" \
942 : "=&r" (__xx.__ll) \
943 : "g" ((USItype) (u)), \
944 "g" ((USItype)(v))); \
945 (w1) = __xx.__i.__h; (w0) = __xx.__i.__l; })
946#endif /* __pyr__ */
947
948/***************************************
949 ************** RT/ROMP **************
950 ***************************************/
951#if defined(__ibm032__) /* RT/ROMP */ && W_TYPE_SIZE == 32
952#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
953 __asm__ ("a %1,%5\n" \
954 "ae %0,%3" \
955 : "=r" ((USItype)(sh)), \
956 "=&r" ((USItype)(sl)) \
957 : "%0" ((USItype)(ah)), \
958 "r" ((USItype)(bh)), \
959 "%1" ((USItype)(al)), \
960 "r" ((USItype)(bl)))
961#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
962 __asm__ ("s %1,%5\n" \
963 "se %0,%3" \
964 : "=r" ((USItype)(sh)), \
965 "=&r" ((USItype)(sl)) \
966 : "0" ((USItype)(ah)), \
967 "r" ((USItype)(bh)), \
968 "1" ((USItype)(al)), \
969 "r" ((USItype)(bl)))
970#define umul_ppmm(ph, pl, m0, m1) \
971do { \
972 USItype __m0 = (m0), __m1 = (m1); \
973 __asm__ ( \
974 "s r2,r2\n" \
975 "mts r10,%2\n" \
976 "m r2,%3\n" \
977 "m r2,%3\n" \
978 "m r2,%3\n" \
979 "m r2,%3\n" \
980 "m r2,%3\n" \
981 "m r2,%3\n" \
982 "m r2,%3\n" \
983 "m r2,%3\n" \
984 "m r2,%3\n" \
985 "m r2,%3\n" \
986 "m r2,%3\n" \
987 "m r2,%3\n" \
988 "m r2,%3\n" \
989 "m r2,%3\n" \
990 "m r2,%3\n" \
991 "m r2,%3\n" \
992 "cas %0,r2,r0\n" \
993 "mfs r10,%1" \
994 : "=r" ((USItype)(ph)), \
995 "=r" ((USItype)(pl)) \
996 : "%r" (__m0), \
997 "r" (__m1) \
998 : "r2"); \
999 (ph) += ((((SItype) __m0 >> 31) & __m1) \
1000 + (((SItype) __m1 >> 31) & __m0)); \
1001} while (0)
1002#define UMUL_TIME 20
1003#define UDIV_TIME 200
1004#define count_leading_zeros(count, x) \
1005do { \
1006 if ((x) >= 0x10000) \
1007 __asm__ ("clz %0,%1" \
1008 : "=r" ((USItype)(count)) \
1009 : "r" ((USItype)(x) >> 16)); \
1010 else { \
1011 __asm__ ("clz %0,%1" \
1012 : "=r" ((USItype)(count)) \
1013 : "r" ((USItype)(x))); \
1014 (count) += 16; \
1015 } \
1016} while (0)
1017#endif /* RT/ROMP */
1018
1019/***************************************
1020 ************** SH2 ******************
1021 ***************************************/
1022#if (defined(__sh2__) || defined(__sh3__) || defined(__SH4__)) \
1023 && W_TYPE_SIZE == 32
1024#define umul_ppmm(w1, w0, u, v) \
1025 __asm__ ( \
1026 "dmulu.l %2,%3\n" \
1027 "sts macl,%1\n" \
1028 "sts mach,%0" \
1029 : "=r" ((USItype)(w1)), \
1030 "=r" ((USItype)(w0)) \
1031 : "r" ((USItype)(u)), \
1032 "r" ((USItype)(v)) \
1033 : "macl", "mach")
1034#define UMUL_TIME 5
1035#endif
1036
1037/***************************************
1038 ************** SPARC ****************
1039 ***************************************/
1040#if defined(__sparc__) && W_TYPE_SIZE == 32
1041#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
1042 __asm__ ("addcc %r4,%5,%1\n" \
1043 "addx %r2,%3,%0" \
1044 : "=r" ((USItype)(sh)), \
1045 "=&r" ((USItype)(sl)) \
1046 : "%rJ" ((USItype)(ah)), \
1047 "rI" ((USItype)(bh)), \
1048 "%rJ" ((USItype)(al)), \
1049 "rI" ((USItype)(bl)) \
1050 __CLOBBER_CC)
1051#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
1052 __asm__ ("subcc %r4,%5,%1\n" \
1053 "subx %r2,%3,%0" \
1054 : "=r" ((USItype)(sh)), \
1055 "=&r" ((USItype)(sl)) \
1056 : "rJ" ((USItype)(ah)), \
1057 "rI" ((USItype)(bh)), \
1058 "rJ" ((USItype)(al)), \
1059 "rI" ((USItype)(bl)) \
1060 __CLOBBER_CC)
1061#if defined(__sparc_v8__)
1062/* Don't match immediate range because, 1) it is not often useful,
1063 2) the 'I' flag thinks of the range as a 13 bit signed interval,
1064 while we want to match a 13 bit interval, sign extended to 32 bits,
1065 but INTERPRETED AS UNSIGNED. */
1066#define umul_ppmm(w1, w0, u, v) \
1067 __asm__ ("umul %2,%3,%1;rd %%y,%0" \
1068 : "=r" ((USItype)(w1)), \
1069 "=r" ((USItype)(w0)) \
1070 : "r" ((USItype)(u)), \
1071 "r" ((USItype)(v)))
1072#define UMUL_TIME 5
1073#ifndef SUPERSPARC /* SuperSPARC's udiv only handles 53 bit dividends */
1074#define udiv_qrnnd(q, r, n1, n0, d) \
1075do { \
1076 USItype __q; \
1077 __asm__ ("mov %1,%%y;nop;nop;nop;udiv %2,%3,%0" \
1078 : "=r" ((USItype)(__q)) \
1079 : "r" ((USItype)(n1)), \
1080 "r" ((USItype)(n0)), \
1081 "r" ((USItype)(d))); \
1082 (r) = (n0) - __q * (d); \
1083 (q) = __q; \
1084} while (0)
1085#define UDIV_TIME 25
1086#endif /* SUPERSPARC */
1087#else /* ! __sparc_v8__ */
1088#if defined(__sparclite__)
1089/* This has hardware multiply but not divide. It also has two additional
1090 instructions scan (ffs from high bit) and divscc. */
1091#define umul_ppmm(w1, w0, u, v) \
1092 __asm__ ("umul %2,%3,%1;rd %%y,%0" \
1093 : "=r" ((USItype)(w1)), \
1094 "=r" ((USItype)(w0)) \
1095 : "r" ((USItype)(u)), \
1096 "r" ((USItype)(v)))
1097#define UMUL_TIME 5
1098#define udiv_qrnnd(q, r, n1, n0, d) \
1099 __asm__ ("! Inlined udiv_qrnnd\n" \
1100 "wr %%g0,%2,%%y ! Not a delayed write for sparclite\n" \
1101 "tst %%g0\n" \
1102 "divscc %3,%4,%%g1\n" \
1103 "divscc %%g1,%4,%%g1\n" \
1104 "divscc %%g1,%4,%%g1\n" \
1105 "divscc %%g1,%4,%%g1\n" \
1106 "divscc %%g1,%4,%%g1\n" \
1107 "divscc %%g1,%4,%%g1\n" \
1108 "divscc %%g1,%4,%%g1\n" \
1109 "divscc %%g1,%4,%%g1\n" \
1110 "divscc %%g1,%4,%%g1\n" \
1111 "divscc %%g1,%4,%%g1\n" \
1112 "divscc %%g1,%4,%%g1\n" \
1113 "divscc %%g1,%4,%%g1\n" \
1114 "divscc %%g1,%4,%%g1\n" \
1115 "divscc %%g1,%4,%%g1\n" \
1116 "divscc %%g1,%4,%%g1\n" \
1117 "divscc %%g1,%4,%%g1\n" \
1118 "divscc %%g1,%4,%%g1\n" \
1119 "divscc %%g1,%4,%%g1\n" \
1120 "divscc %%g1,%4,%%g1\n" \
1121 "divscc %%g1,%4,%%g1\n" \
1122 "divscc %%g1,%4,%%g1\n" \
1123 "divscc %%g1,%4,%%g1\n" \
1124 "divscc %%g1,%4,%%g1\n" \
1125 "divscc %%g1,%4,%%g1\n" \
1126 "divscc %%g1,%4,%%g1\n" \
1127 "divscc %%g1,%4,%%g1\n" \
1128 "divscc %%g1,%4,%%g1\n" \
1129 "divscc %%g1,%4,%%g1\n" \
1130 "divscc %%g1,%4,%%g1\n" \
1131 "divscc %%g1,%4,%%g1\n" \
1132 "divscc %%g1,%4,%%g1\n" \
1133 "divscc %%g1,%4,%0\n" \
1134 "rd %%y,%1\n" \
1135 "bl,a 1f\n" \
1136 "add %1,%4,%1\n" \
1137 "1: ! End of inline udiv_qrnnd" \
1138 : "=r" ((USItype)(q)), \
1139 "=r" ((USItype)(r)) \
1140 : "r" ((USItype)(n1)), \
1141 "r" ((USItype)(n0)), \
1142 "rI" ((USItype)(d)) \
1143 : "%g1" __AND_CLOBBER_CC)
1144#define UDIV_TIME 37
1145#define count_leading_zeros(count, x) \
1146 __asm__ ("scan %1,0,%0" \
1147 : "=r" ((USItype)(x)) \
1148 : "r" ((USItype)(count)))
1149/* Early sparclites return 63 for an argument of 0, but they warn that future
1150 implementations might change this. Therefore, leave COUNT_LEADING_ZEROS_0
1151 undefined. */
1152#endif /* __sparclite__ */
1153#endif /* __sparc_v8__ */
1154 /* Default to sparc v7 versions of umul_ppmm and udiv_qrnnd. */
1155#ifndef umul_ppmm
1156#define umul_ppmm(w1, w0, u, v) \
1157 __asm__ ("! Inlined umul_ppmm\n" \
1158 "wr %%g0,%2,%%y ! SPARC has 0-3 delay insn after a wr\n" \
1159 "sra %3,31,%%g2 ! Don't move this insn\n" \
1160 "and %2,%%g2,%%g2 ! Don't move this insn\n" \
1161 "andcc %%g0,0,%%g1 ! Don't move this insn\n" \
1162 "mulscc %%g1,%3,%%g1\n" \
1163 "mulscc %%g1,%3,%%g1\n" \
1164 "mulscc %%g1,%3,%%g1\n" \
1165 "mulscc %%g1,%3,%%g1\n" \
1166 "mulscc %%g1,%3,%%g1\n" \
1167 "mulscc %%g1,%3,%%g1\n" \
1168 "mulscc %%g1,%3,%%g1\n" \
1169 "mulscc %%g1,%3,%%g1\n" \
1170 "mulscc %%g1,%3,%%g1\n" \
1171 "mulscc %%g1,%3,%%g1\n" \
1172 "mulscc %%g1,%3,%%g1\n" \
1173 "mulscc %%g1,%3,%%g1\n" \
1174 "mulscc %%g1,%3,%%g1\n" \
1175 "mulscc %%g1,%3,%%g1\n" \
1176 "mulscc %%g1,%3,%%g1\n" \
1177 "mulscc %%g1,%3,%%g1\n" \
1178 "mulscc %%g1,%3,%%g1\n" \
1179 "mulscc %%g1,%3,%%g1\n" \
1180 "mulscc %%g1,%3,%%g1\n" \
1181 "mulscc %%g1,%3,%%g1\n" \
1182 "mulscc %%g1,%3,%%g1\n" \
1183 "mulscc %%g1,%3,%%g1\n" \
1184 "mulscc %%g1,%3,%%g1\n" \
1185 "mulscc %%g1,%3,%%g1\n" \
1186 "mulscc %%g1,%3,%%g1\n" \
1187 "mulscc %%g1,%3,%%g1\n" \
1188 "mulscc %%g1,%3,%%g1\n" \
1189 "mulscc %%g1,%3,%%g1\n" \
1190 "mulscc %%g1,%3,%%g1\n" \
1191 "mulscc %%g1,%3,%%g1\n" \
1192 "mulscc %%g1,%3,%%g1\n" \
1193 "mulscc %%g1,%3,%%g1\n" \
1194 "mulscc %%g1,0,%%g1\n" \
1195 "add %%g1,%%g2,%0\n" \
1196 "rd %%y,%1" \
1197 : "=r" ((USItype)(w1)), \
1198 "=r" ((USItype)(w0)) \
1199 : "%rI" ((USItype)(u)), \
1200 "r" ((USItype)(v)) \
1201 : "%g1", "%g2" __AND_CLOBBER_CC)
1202#define UMUL_TIME 39 /* 39 instructions */
1203/* It's quite necessary to add this much assembler for the sparc.
1204 The default udiv_qrnnd (in C) is more than 10 times slower! */
1205#define udiv_qrnnd(q, r, n1, n0, d) \
1206 __asm__ ("! Inlined udiv_qrnnd\n\t" \
1207 "mov 32,%%g1\n\t" \
1208 "subcc %1,%2,%%g0\n\t" \
1209 "1: bcs 5f\n\t" \
1210 "addxcc %0,%0,%0 ! shift n1n0 and a q-bit in lsb\n\t" \
1211 "sub %1,%2,%1 ! this kills msb of n\n\t" \
1212 "addx %1,%1,%1 ! so this can't give carry\n\t" \
1213 "subcc %%g1,1,%%g1\n\t" \
1214 "2: bne 1b\n\t" \
1215 "subcc %1,%2,%%g0\n\t" \
1216 "bcs 3f\n\t" \
1217 "addxcc %0,%0,%0 ! shift n1n0 and a q-bit in lsb\n\t" \
1218 "b 3f\n\t" \
1219 "sub %1,%2,%1 ! this kills msb of n\n\t" \
1220 "4: sub %1,%2,%1\n\t" \
1221 "5: addxcc %1,%1,%1\n\t" \
1222 "bcc 2b\n\t" \
1223 "subcc %%g1,1,%%g1\n\t" \
1224 "! Got carry from n. Subtract next step to cancel this carry.\n\t" \
1225 "bne 4b\n\t" \
1226 "addcc %0,%0,%0 ! shift n1n0 and a 0-bit in lsb\n\t" \
1227 "sub %1,%2,%1\n\t" \
1228 "3: xnor %0,0,%0\n\t" \
1229 "! End of inline udiv_qrnnd\n" \
1230 : "=&r" ((USItype)(q)), \
1231 "=&r" ((USItype)(r)) \
1232 : "r" ((USItype)(d)), \
1233 "1" ((USItype)(n1)), \
1234 "0" ((USItype)(n0)) : "%g1", "cc")
1235#define UDIV_TIME (3+7*32) /* 7 instructions/iteration. 32 iterations. */
1236#endif
1237#endif /* __sparc__ */
1238
1239/***************************************
1240 ************** VAX ******************
1241 ***************************************/
1242#if defined(__vax__) && W_TYPE_SIZE == 32
1243#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
1244 __asm__ ("addl2 %5,%1\n" \
1245 "adwc %3,%0" \
1246 : "=g" ((USItype)(sh)), \
1247 "=&g" ((USItype)(sl)) \
1248 : "%0" ((USItype)(ah)), \
1249 "g" ((USItype)(bh)), \
1250 "%1" ((USItype)(al)), \
1251 "g" ((USItype)(bl)))
1252#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
1253 __asm__ ("subl2 %5,%1\n" \
1254 "sbwc %3,%0" \
1255 : "=g" ((USItype)(sh)), \
1256 "=&g" ((USItype)(sl)) \
1257 : "0" ((USItype)(ah)), \
1258 "g" ((USItype)(bh)), \
1259 "1" ((USItype)(al)), \
1260 "g" ((USItype)(bl)))
1261#define umul_ppmm(xh, xl, m0, m1) \
1262do { \
1263 union {UDItype __ll; \
1264 struct {USItype __l, __h; } __i; \
1265 } __xx; \
1266 USItype __m0 = (m0), __m1 = (m1); \
1267 __asm__ ("emul %1,%2,$0,%0" \
1268 : "=g" (__xx.__ll) \
1269 : "g" (__m0), \
1270 "g" (__m1)); \
1271 (xh) = __xx.__i.__h; (xl) = __xx.__i.__l; \
1272 (xh) += ((((SItype) __m0 >> 31) & __m1) \
1273 + (((SItype) __m1 >> 31) & __m0)); \
1274} while (0)
1275#define sdiv_qrnnd(q, r, n1, n0, d) \
1276do { \
1277 union {DItype __ll; \
1278 struct {SItype __l, __h; } __i; \
1279 } __xx; \
1280 __xx.__i.__h = n1; __xx.__i.__l = n0; \
1281 __asm__ ("ediv %3,%2,%0,%1" \
1282 : "=g" (q), "=g" (r) \
1283 : "g" (__xx.__ll), "g" (d)); \
1284} while (0)
1285#endif /* __vax__ */
1286
1287/***************************************
1288 ************** Z8000 ****************
1289 ***************************************/
1290#if defined(__z8000__) && W_TYPE_SIZE == 16
1291#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
1292 __asm__ ("add %H1,%H5\n\tadc %H0,%H3" \
1293 : "=r" ((unsigned int)(sh)), \
1294 "=&r" ((unsigned int)(sl)) \
1295 : "%0" ((unsigned int)(ah)), \
1296 "r" ((unsigned int)(bh)), \
1297 "%1" ((unsigned int)(al)), \
1298 "rQR" ((unsigned int)(bl)))
1299#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
1300 __asm__ ("sub %H1,%H5\n\tsbc %H0,%H3" \
1301 : "=r" ((unsigned int)(sh)), \
1302 "=&r" ((unsigned int)(sl)) \
1303 : "0" ((unsigned int)(ah)), \
1304 "r" ((unsigned int)(bh)), \
1305 "1" ((unsigned int)(al)), \
1306 "rQR" ((unsigned int)(bl)))
1307#define umul_ppmm(xh, xl, m0, m1) \
1308do { \
1309 union {long int __ll; \
1310 struct {unsigned int __h, __l; } __i; \
1311 } __xx; \
1312 unsigned int __m0 = (m0), __m1 = (m1); \
1313 __asm__ ("mult %S0,%H3" \
1314 : "=r" (__xx.__i.__h), \
1315 "=r" (__xx.__i.__l) \
1316 : "%1" (__m0), \
1317 "rQR" (__m1)); \
1318 (xh) = __xx.__i.__h; (xl) = __xx.__i.__l; \
1319 (xh) += ((((signed int) __m0 >> 15) & __m1) \
1320 + (((signed int) __m1 >> 15) & __m0)); \
1321} while (0)
1322#endif /* __z8000__ */
1323
1324#endif /* __GNUC__ */
1325
1326/***************************************
1327 *********** Generic Versions ********
1328 ***************************************/
1329#if !defined(umul_ppmm) && defined(__umulsidi3)
1330#define umul_ppmm(ph, pl, m0, m1) \
1331{ \
1332 UDWtype __ll = __umulsidi3(m0, m1); \
1333 ph = (UWtype) (__ll >> W_TYPE_SIZE); \
1334 pl = (UWtype) __ll; \
1335}
1336#endif
1337
1338#if !defined(__umulsidi3)
1339#define __umulsidi3(u, v) \
1340 ({UWtype __hi, __lo; \
1341 umul_ppmm(__hi, __lo, u, v); \
1342 ((UDWtype) __hi << W_TYPE_SIZE) | __lo; })
1343#endif
1344
1345 /* If this machine has no inline assembler, use C macros. */
1346
1347#if !defined(add_ssaaaa)
1348#define add_ssaaaa(sh, sl, ah, al, bh, bl) \
1349do { \
1350 UWtype __x; \
1351 __x = (al) + (bl); \
1352 (sh) = (ah) + (bh) + (__x < (al)); \
1353 (sl) = __x; \
1354} while (0)
1355#endif
1356
1357#if !defined(sub_ddmmss)
1358#define sub_ddmmss(sh, sl, ah, al, bh, bl) \
1359do { \
1360 UWtype __x; \
1361 __x = (al) - (bl); \
1362 (sh) = (ah) - (bh) - (__x > (al)); \
1363 (sl) = __x; \
1364} while (0)
1365#endif
1366
1367#if !defined(umul_ppmm)
1368#define umul_ppmm(w1, w0, u, v) \
1369do { \
1370 UWtype __x0, __x1, __x2, __x3; \
1371 UHWtype __ul, __vl, __uh, __vh; \
1372 UWtype __u = (u), __v = (v); \
1373 \
1374 __ul = __ll_lowpart(__u); \
1375 __uh = __ll_highpart(__u); \
1376 __vl = __ll_lowpart(__v); \
1377 __vh = __ll_highpart(__v); \
1378 \
1379 __x0 = (UWtype) __ul * __vl; \
1380 __x1 = (UWtype) __ul * __vh; \
1381 __x2 = (UWtype) __uh * __vl; \
1382 __x3 = (UWtype) __uh * __vh; \
1383 \
1384 __x1 += __ll_highpart(__x0);/* this can't give carry */ \
1385 __x1 += __x2; /* but this indeed can */ \
1386 if (__x1 < __x2) /* did we get it? */ \
1387 __x3 += __ll_B; /* yes, add it in the proper pos. */ \
1388 \
1389 (w1) = __x3 + __ll_highpart(__x1); \
1390 (w0) = (__ll_lowpart(__x1) << W_TYPE_SIZE/2) + __ll_lowpart(__x0); \
1391} while (0)
1392#endif
1393
1394#if !defined(umul_ppmm)
1395#define smul_ppmm(w1, w0, u, v) \
1396do { \
1397 UWtype __w1; \
1398 UWtype __m0 = (u), __m1 = (v); \
1399 umul_ppmm(__w1, w0, __m0, __m1); \
1400 (w1) = __w1 - (-(__m0 >> (W_TYPE_SIZE - 1)) & __m1) \
1401 - (-(__m1 >> (W_TYPE_SIZE - 1)) & __m0); \
1402} while (0)
1403#endif
1404
1405 /* Define this unconditionally, so it can be used for debugging. */
1406#define __udiv_qrnnd_c(q, r, n1, n0, d) \
1407do { \
1408 UWtype __d1, __d0, __q1, __q0, __r1, __r0, __m; \
1409 __d1 = __ll_highpart(d); \
1410 __d0 = __ll_lowpart(d); \
1411 \
1412 __r1 = (n1) % __d1; \
1413 __q1 = (n1) / __d1; \
1414 __m = (UWtype) __q1 * __d0; \
1415 __r1 = __r1 * __ll_B | __ll_highpart(n0); \
1416 if (__r1 < __m) { \
1417 __q1--, __r1 += (d); \
1418 if (__r1 >= (d)) /* i.e. we didn't get carry when adding to __r1 */ \
1419 if (__r1 < __m) \
1420 __q1--, __r1 += (d); \
1421 } \
1422 __r1 -= __m; \
1423 \
1424 __r0 = __r1 % __d1; \
1425 __q0 = __r1 / __d1; \
1426 __m = (UWtype) __q0 * __d0; \
1427 __r0 = __r0 * __ll_B | __ll_lowpart(n0); \
1428 if (__r0 < __m) { \
1429 __q0--, __r0 += (d); \
1430 if (__r0 >= (d)) \
1431 if (__r0 < __m) \
1432 __q0--, __r0 += (d); \
1433 } \
1434 __r0 -= __m; \
1435 \
1436 (q) = (UWtype) __q1 * __ll_B | __q0; \
1437 (r) = __r0; \
1438} while (0)
1439
1440/* If the processor has no udiv_qrnnd but sdiv_qrnnd, go through
1441 __udiv_w_sdiv (defined in libgcc or elsewhere). */
1442#if !defined(udiv_qrnnd) && defined(sdiv_qrnnd)
1443#define udiv_qrnnd(q, r, nh, nl, d) \
1444do { \
1445 UWtype __r; \
1446 (q) = __MPN(udiv_w_sdiv) (&__r, nh, nl, d); \
1447 (r) = __r; \
1448} while (0)
1449#endif
1450
1451 /* If udiv_qrnnd was not defined for this processor, use __udiv_qrnnd_c. */
1452#if !defined(udiv_qrnnd)
1453#define UDIV_NEEDS_NORMALIZATION 1
1454#define udiv_qrnnd __udiv_qrnnd_c
1455#endif
1456
1457#undef count_leading_zeros
1458#if !defined(count_leading_zeros)
1459 extern
1460#ifdef __STDC__
1461 const
1462#endif
1463 unsigned char __clz_tab[];
1464#define count_leading_zeros(count, x) \
1465do { \
1466 UWtype __xr = (x); \
1467 UWtype __a; \
1468 \
1469 if (W_TYPE_SIZE <= 32) { \
1470 __a = __xr < ((UWtype) 1 << 2*__BITS4) \
1471 ? (__xr < ((UWtype) 1 << __BITS4) ? 0 : __BITS4) \
1472 : (__xr < ((UWtype) 1 << 3*__BITS4) ? 2*__BITS4 : 3*__BITS4); \
1473 } \
1474 else { \
1475 for (__a = W_TYPE_SIZE - 8; __a > 0; __a -= 8) \
1476 if (((__xr >> __a) & 0xff) != 0) \
1477 break; \
1478 } \
1479 \
1480 (count) = W_TYPE_SIZE - (__clz_tab[__xr >> __a] + __a); \
1481} while (0)
1482 /* This version gives a well-defined value for zero. */
1483#define COUNT_LEADING_ZEROS_0 W_TYPE_SIZE
1484#endif
1485
1486#if !defined(count_trailing_zeros)
1487/* Define count_trailing_zeros using count_leading_zeros. The latter might be
1488 defined in asm, but if it is not, the C version above is good enough. */
1489#define count_trailing_zeros(count, x) \
1490do { \
1491 UWtype __ctz_x = (x); \
1492 UWtype __ctz_c; \
1493 count_leading_zeros(__ctz_c, __ctz_x & -__ctz_x); \
1494 (count) = W_TYPE_SIZE - 1 - __ctz_c; \
1495} while (0)
1496#endif
1497
1498#ifndef UDIV_NEEDS_NORMALIZATION
1499#define UDIV_NEEDS_NORMALIZATION 0
1500#endif
diff --git a/lib/mpi/mpi-add.c b/lib/mpi/mpi-add.c
new file mode 100644
index 000000000000..f56b9ba295e6
--- /dev/null
+++ b/lib/mpi/mpi-add.c
@@ -0,0 +1,234 @@
1/* mpi-add.c - MPI functions
2 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
3 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31
32/****************
33 * Add the unsigned integer V to the mpi-integer U and store the
34 * result in W. U and V may be the same.
35 */
36int mpi_add_ui(MPI w, const MPI u, unsigned long v)
37{
38 mpi_ptr_t wp, up;
39 mpi_size_t usize, wsize;
40 int usign, wsign;
41
42 usize = u->nlimbs;
43 usign = u->sign;
44 wsign = 0;
45
46 /* If not space for W (and possible carry), increase space. */
47 wsize = usize + 1;
48 if (w->alloced < wsize)
49 if (mpi_resize(w, wsize) < 0)
50 return -ENOMEM;
51
52 /* These must be after realloc (U may be the same as W). */
53 up = u->d;
54 wp = w->d;
55
56 if (!usize) { /* simple */
57 wp[0] = v;
58 wsize = v ? 1 : 0;
59 } else if (!usign) { /* mpi is not negative */
60 mpi_limb_t cy;
61 cy = mpihelp_add_1(wp, up, usize, v);
62 wp[usize] = cy;
63 wsize = usize + cy;
64 } else { /* The signs are different. Need exact comparison to determine
65 * which operand to subtract from which. */
66 if (usize == 1 && up[0] < v) {
67 wp[0] = v - up[0];
68 wsize = 1;
69 } else {
70 mpihelp_sub_1(wp, up, usize, v);
71 /* Size can decrease with at most one limb. */
72 wsize = usize - (wp[usize - 1] == 0);
73 wsign = 1;
74 }
75 }
76
77 w->nlimbs = wsize;
78 w->sign = wsign;
79 return 0;
80}
81
82int mpi_add(MPI w, MPI u, MPI v)
83{
84 mpi_ptr_t wp, up, vp;
85 mpi_size_t usize, vsize, wsize;
86 int usign, vsign, wsign;
87
88 if (u->nlimbs < v->nlimbs) { /* Swap U and V. */
89 usize = v->nlimbs;
90 usign = v->sign;
91 vsize = u->nlimbs;
92 vsign = u->sign;
93 wsize = usize + 1;
94 if (RESIZE_IF_NEEDED(w, wsize) < 0)
95 return -ENOMEM;
96 /* These must be after realloc (u or v may be the same as w). */
97 up = v->d;
98 vp = u->d;
99 } else {
100 usize = u->nlimbs;
101 usign = u->sign;
102 vsize = v->nlimbs;
103 vsign = v->sign;
104 wsize = usize + 1;
105 if (RESIZE_IF_NEEDED(w, wsize) < 0)
106 return -ENOMEM;
107 /* These must be after realloc (u or v may be the same as w). */
108 up = u->d;
109 vp = v->d;
110 }
111 wp = w->d;
112 wsign = 0;
113
114 if (!vsize) { /* simple */
115 MPN_COPY(wp, up, usize);
116 wsize = usize;
117 wsign = usign;
118 } else if (usign != vsign) { /* different sign */
119 /* This test is right since USIZE >= VSIZE */
120 if (usize != vsize) {
121 mpihelp_sub(wp, up, usize, vp, vsize);
122 wsize = usize;
123 MPN_NORMALIZE(wp, wsize);
124 wsign = usign;
125 } else if (mpihelp_cmp(up, vp, usize) < 0) {
126 mpihelp_sub_n(wp, vp, up, usize);
127 wsize = usize;
128 MPN_NORMALIZE(wp, wsize);
129 if (!usign)
130 wsign = 1;
131 } else {
132 mpihelp_sub_n(wp, up, vp, usize);
133 wsize = usize;
134 MPN_NORMALIZE(wp, wsize);
135 if (usign)
136 wsign = 1;
137 }
138 } else { /* U and V have same sign. Add them. */
139 mpi_limb_t cy = mpihelp_add(wp, up, usize, vp, vsize);
140 wp[usize] = cy;
141 wsize = usize + cy;
142 if (usign)
143 wsign = 1;
144 }
145
146 w->nlimbs = wsize;
147 w->sign = wsign;
148 return 0;
149}
150
151/****************
152 * Subtract the unsigned integer V from the mpi-integer U and store the
153 * result in W.
154 */
155int mpi_sub_ui(MPI w, MPI u, unsigned long v)
156{
157 mpi_ptr_t wp, up;
158 mpi_size_t usize, wsize;
159 int usign, wsign;
160
161 usize = u->nlimbs;
162 usign = u->sign;
163 wsign = 0;
164
165 /* If not space for W (and possible carry), increase space. */
166 wsize = usize + 1;
167 if (w->alloced < wsize)
168 if (mpi_resize(w, wsize) < 0)
169 return -ENOMEM;
170
171 /* These must be after realloc (U may be the same as W). */
172 up = u->d;
173 wp = w->d;
174
175 if (!usize) { /* simple */
176 wp[0] = v;
177 wsize = v ? 1 : 0;
178 wsign = 1;
179 } else if (usign) { /* mpi and v are negative */
180 mpi_limb_t cy;
181 cy = mpihelp_add_1(wp, up, usize, v);
182 wp[usize] = cy;
183 wsize = usize + cy;
184 } else { /* The signs are different. Need exact comparison to determine
185 * which operand to subtract from which. */
186 if (usize == 1 && up[0] < v) {
187 wp[0] = v - up[0];
188 wsize = 1;
189 wsign = 1;
190 } else {
191 mpihelp_sub_1(wp, up, usize, v);
192 /* Size can decrease with at most one limb. */
193 wsize = usize - (wp[usize - 1] == 0);
194 }
195 }
196
197 w->nlimbs = wsize;
198 w->sign = wsign;
199 return 0;
200}
201
202int mpi_sub(MPI w, MPI u, MPI v)
203{
204 int rc;
205
206 if (w == v) {
207 MPI vv;
208 if (mpi_copy(&vv, v) < 0)
209 return -ENOMEM;
210 vv->sign = !vv->sign;
211 rc = mpi_add(w, u, vv);
212 mpi_free(vv);
213 } else {
214 /* fixme: this is not thread-save (we temp. modify v) */
215 v->sign = !v->sign;
216 rc = mpi_add(w, u, v);
217 v->sign = !v->sign;
218 }
219 return rc;
220}
221
222int mpi_addm(MPI w, MPI u, MPI v, MPI m)
223{
224 if (mpi_add(w, u, v) < 0 || mpi_fdiv_r(w, w, m) < 0)
225 return -ENOMEM;
226 return 0;
227}
228
229int mpi_subm(MPI w, MPI u, MPI v, MPI m)
230{
231 if (mpi_sub(w, u, v) < 0 || mpi_fdiv_r(w, w, m) < 0)
232 return -ENOMEM;
233 return 0;
234}
diff --git a/lib/mpi/mpi-bit.c b/lib/mpi/mpi-bit.c
new file mode 100644
index 000000000000..2f526627e4f5
--- /dev/null
+++ b/lib/mpi/mpi-bit.c
@@ -0,0 +1,217 @@
1/* mpi-bit.c - MPI bit level fucntions
2 * Copyright (C) 1998, 1999 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22#include "longlong.h"
23
24#define A_LIMB_1 ((mpi_limb_t) 1)
25
26/****************
27 * Sometimes we have MSL (most significant limbs) which are 0;
28 * this is for some reasons not good, so this function removes them.
29 */
30void mpi_normalize(MPI a)
31{
32 for (; a->nlimbs && !a->d[a->nlimbs - 1]; a->nlimbs--)
33 ;
34}
35
36/****************
37 * Return the number of bits in A.
38 */
39unsigned mpi_get_nbits(MPI a)
40{
41 unsigned n;
42
43 mpi_normalize(a);
44
45 if (a->nlimbs) {
46 mpi_limb_t alimb = a->d[a->nlimbs - 1];
47 if (alimb)
48 count_leading_zeros(n, alimb);
49 else
50 n = BITS_PER_MPI_LIMB;
51 n = BITS_PER_MPI_LIMB - n + (a->nlimbs - 1) * BITS_PER_MPI_LIMB;
52 } else
53 n = 0;
54 return n;
55}
56EXPORT_SYMBOL_GPL(mpi_get_nbits);
57
58/****************
59 * Test whether bit N is set.
60 */
61int mpi_test_bit(MPI a, unsigned n)
62{
63 unsigned limbno, bitno;
64 mpi_limb_t limb;
65
66 limbno = n / BITS_PER_MPI_LIMB;
67 bitno = n % BITS_PER_MPI_LIMB;
68
69 if (limbno >= a->nlimbs)
70 return 0; /* too far left: this is a 0 */
71 limb = a->d[limbno];
72 return (limb & (A_LIMB_1 << bitno)) ? 1 : 0;
73}
74
75/****************
76 * Set bit N of A.
77 */
78int mpi_set_bit(MPI a, unsigned n)
79{
80 unsigned limbno, bitno;
81
82 limbno = n / BITS_PER_MPI_LIMB;
83 bitno = n % BITS_PER_MPI_LIMB;
84
85 if (limbno >= a->nlimbs) { /* resize */
86 if (a->alloced >= limbno)
87 if (mpi_resize(a, limbno + 1) < 0)
88 return -ENOMEM;
89 a->nlimbs = limbno + 1;
90 }
91 a->d[limbno] |= (A_LIMB_1 << bitno);
92 return 0;
93}
94
95/****************
96 * Set bit N of A. and clear all bits above
97 */
98int mpi_set_highbit(MPI a, unsigned n)
99{
100 unsigned limbno, bitno;
101
102 limbno = n / BITS_PER_MPI_LIMB;
103 bitno = n % BITS_PER_MPI_LIMB;
104
105 if (limbno >= a->nlimbs) { /* resize */
106 if (a->alloced >= limbno)
107 if (mpi_resize(a, limbno + 1) < 0)
108 return -ENOMEM;
109 a->nlimbs = limbno + 1;
110 }
111 a->d[limbno] |= (A_LIMB_1 << bitno);
112 for (bitno++; bitno < BITS_PER_MPI_LIMB; bitno++)
113 a->d[limbno] &= ~(A_LIMB_1 << bitno);
114 a->nlimbs = limbno + 1;
115 return 0;
116}
117
118/****************
119 * clear bit N of A and all bits above
120 */
121void mpi_clear_highbit(MPI a, unsigned n)
122{
123 unsigned limbno, bitno;
124
125 limbno = n / BITS_PER_MPI_LIMB;
126 bitno = n % BITS_PER_MPI_LIMB;
127
128 if (limbno >= a->nlimbs)
129 return; /* not allocated, so need to clear bits :-) */
130
131 for (; bitno < BITS_PER_MPI_LIMB; bitno++)
132 a->d[limbno] &= ~(A_LIMB_1 << bitno);
133 a->nlimbs = limbno + 1;
134}
135
136/****************
137 * Clear bit N of A.
138 */
139void mpi_clear_bit(MPI a, unsigned n)
140{
141 unsigned limbno, bitno;
142
143 limbno = n / BITS_PER_MPI_LIMB;
144 bitno = n % BITS_PER_MPI_LIMB;
145
146 if (limbno >= a->nlimbs)
147 return; /* don't need to clear this bit, it's to far to left */
148 a->d[limbno] &= ~(A_LIMB_1 << bitno);
149}
150
151/****************
152 * Shift A by N bits to the right
153 * FIXME: should use alloc_limb if X and A are same.
154 */
155int mpi_rshift(MPI x, MPI a, unsigned n)
156{
157 mpi_ptr_t xp;
158 mpi_size_t xsize;
159
160 xsize = a->nlimbs;
161 x->sign = a->sign;
162 if (RESIZE_IF_NEEDED(x, (size_t) xsize) < 0)
163 return -ENOMEM;
164 xp = x->d;
165
166 if (xsize) {
167 mpihelp_rshift(xp, a->d, xsize, n);
168 MPN_NORMALIZE(xp, xsize);
169 }
170 x->nlimbs = xsize;
171 return 0;
172}
173
174/****************
175 * Shift A by COUNT limbs to the left
176 * This is used only within the MPI library
177 */
178int mpi_lshift_limbs(MPI a, unsigned int count)
179{
180 mpi_ptr_t ap = a->d;
181 int n = a->nlimbs;
182 int i;
183
184 if (!count || !n)
185 return 0;
186
187 if (RESIZE_IF_NEEDED(a, n + count) < 0)
188 return -ENOMEM;
189
190 for (i = n - 1; i >= 0; i--)
191 ap[i + count] = ap[i];
192 for (i = 0; i < count; i++)
193 ap[i] = 0;
194 a->nlimbs += count;
195 return 0;
196}
197
198/****************
199 * Shift A by COUNT limbs to the right
200 * This is used only within the MPI library
201 */
202void mpi_rshift_limbs(MPI a, unsigned int count)
203{
204 mpi_ptr_t ap = a->d;
205 mpi_size_t n = a->nlimbs;
206 unsigned int i;
207
208 if (count >= n) {
209 a->nlimbs = 0;
210 return;
211 }
212
213 for (i = 0; i < n - count; i++)
214 ap[i] = ap[i + count];
215 ap[i] = 0;
216 a->nlimbs -= count;
217}
diff --git a/lib/mpi/mpi-cmp.c b/lib/mpi/mpi-cmp.c
new file mode 100644
index 000000000000..914bc42a8a80
--- /dev/null
+++ b/lib/mpi/mpi-cmp.c
@@ -0,0 +1,68 @@
1/* mpi-cmp.c - MPI functions
2 * Copyright (C) 1998, 1999 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22
23int mpi_cmp_ui(MPI u, unsigned long v)
24{
25 mpi_limb_t limb = v;
26
27 mpi_normalize(u);
28 if (!u->nlimbs && !limb)
29 return 0;
30 if (u->sign)
31 return -1;
32 if (u->nlimbs > 1)
33 return 1;
34
35 if (u->d[0] == limb)
36 return 0;
37 else if (u->d[0] > limb)
38 return 1;
39 else
40 return -1;
41}
42
43int mpi_cmp(MPI u, MPI v)
44{
45 mpi_size_t usize, vsize;
46 int cmp;
47
48 mpi_normalize(u);
49 mpi_normalize(v);
50 usize = u->nlimbs;
51 vsize = v->nlimbs;
52 if (!u->sign && v->sign)
53 return 1;
54 if (u->sign && !v->sign)
55 return -1;
56 if (usize != vsize && !u->sign && !v->sign)
57 return usize - vsize;
58 if (usize != vsize && u->sign && v->sign)
59 return vsize + usize;
60 if (!usize)
61 return 0;
62 cmp = mpihelp_cmp(u->d, v->d, usize);
63 if (!cmp)
64 return 0;
65 if ((cmp < 0 ? 1 : 0) == (u->sign ? 1 : 0))
66 return 1;
67 return -1;
68}
diff --git a/lib/mpi/mpi-div.c b/lib/mpi/mpi-div.c
new file mode 100644
index 000000000000..f68cbbb4d4a4
--- /dev/null
+++ b/lib/mpi/mpi-div.c
@@ -0,0 +1,338 @@
1/* mpi-div.c - MPI functions
2 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
3 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include <linux/string.h>
31#include "mpi-internal.h"
32#include "longlong.h"
33
34int mpi_fdiv_r(MPI rem, MPI dividend, MPI divisor)
35{
36 int rc = -ENOMEM;
37 int divisor_sign = divisor->sign;
38 MPI temp_divisor = NULL;
39
40 /* We need the original value of the divisor after the remainder has been
41 * preliminary calculated. We have to copy it to temporary space if it's
42 * the same variable as REM. */
43 if (rem == divisor) {
44 if (mpi_copy(&temp_divisor, divisor) < 0)
45 goto nomem;
46 divisor = temp_divisor;
47 }
48
49 if (mpi_tdiv_qr(NULL, rem, dividend, divisor) < 0)
50 goto nomem;
51 if (((divisor_sign ? 1 : 0) ^ (dividend->sign ? 1 : 0)) && rem->nlimbs)
52 if (mpi_add(rem, rem, divisor) < 0)
53 goto nomem;
54
55 rc = 0;
56
57nomem:
58 if (temp_divisor)
59 mpi_free(temp_divisor);
60 return rc;
61}
62
63/****************
64 * Division rounding the quotient towards -infinity.
65 * The remainder gets the same sign as the denominator.
66 * rem is optional
67 */
68
69ulong mpi_fdiv_r_ui(MPI rem, MPI dividend, ulong divisor)
70{
71 mpi_limb_t rlimb;
72
73 rlimb = mpihelp_mod_1(dividend->d, dividend->nlimbs, divisor);
74 if (rlimb && dividend->sign)
75 rlimb = divisor - rlimb;
76
77 if (rem) {
78 rem->d[0] = rlimb;
79 rem->nlimbs = rlimb ? 1 : 0;
80 }
81 return rlimb;
82}
83
84int mpi_fdiv_q(MPI quot, MPI dividend, MPI divisor)
85{
86 MPI tmp = mpi_alloc(mpi_get_nlimbs(quot));
87 if (!tmp)
88 return -ENOMEM;
89 mpi_fdiv_qr(quot, tmp, dividend, divisor);
90 mpi_free(tmp);
91 return 0;
92}
93
94int mpi_fdiv_qr(MPI quot, MPI rem, MPI dividend, MPI divisor)
95{
96 int divisor_sign = divisor->sign;
97 MPI temp_divisor = NULL;
98
99 if (quot == divisor || rem == divisor) {
100 if (mpi_copy(&temp_divisor, divisor) < 0)
101 return -ENOMEM;
102 divisor = temp_divisor;
103 }
104
105 if (mpi_tdiv_qr(quot, rem, dividend, divisor) < 0)
106 goto nomem;
107
108 if ((divisor_sign ^ dividend->sign) && rem->nlimbs) {
109 if (mpi_sub_ui(quot, quot, 1) < 0)
110 goto nomem;
111 if (mpi_add(rem, rem, divisor) < 0)
112 goto nomem;
113 }
114
115 if (temp_divisor)
116 mpi_free(temp_divisor);
117
118 return 0;
119
120nomem:
121 mpi_free(temp_divisor);
122 return -ENOMEM;
123}
124
125/* If den == quot, den needs temporary storage.
126 * If den == rem, den needs temporary storage.
127 * If num == quot, num needs temporary storage.
128 * If den has temporary storage, it can be normalized while being copied,
129 * i.e no extra storage should be allocated.
130 */
131
132int mpi_tdiv_r(MPI rem, MPI num, MPI den)
133{
134 return mpi_tdiv_qr(NULL, rem, num, den);
135}
136
137int mpi_tdiv_qr(MPI quot, MPI rem, MPI num, MPI den)
138{
139 int rc = -ENOMEM;
140 mpi_ptr_t np, dp;
141 mpi_ptr_t qp, rp;
142 mpi_size_t nsize = num->nlimbs;
143 mpi_size_t dsize = den->nlimbs;
144 mpi_size_t qsize, rsize;
145 mpi_size_t sign_remainder = num->sign;
146 mpi_size_t sign_quotient = num->sign ^ den->sign;
147 unsigned normalization_steps;
148 mpi_limb_t q_limb;
149 mpi_ptr_t marker[5];
150 int markidx = 0;
151
152 if (!dsize)
153 return -EINVAL;
154
155 memset(marker, 0, sizeof(marker));
156
157 /* Ensure space is enough for quotient and remainder.
158 * We need space for an extra limb in the remainder, because it's
159 * up-shifted (normalized) below. */
160 rsize = nsize + 1;
161 if (mpi_resize(rem, rsize) < 0)
162 goto nomem;
163
164 qsize = rsize - dsize; /* qsize cannot be bigger than this. */
165 if (qsize <= 0) {
166 if (num != rem) {
167 rem->nlimbs = num->nlimbs;
168 rem->sign = num->sign;
169 MPN_COPY(rem->d, num->d, nsize);
170 }
171 if (quot) {
172 /* This needs to follow the assignment to rem, in case the
173 * numerator and quotient are the same. */
174 quot->nlimbs = 0;
175 quot->sign = 0;
176 }
177 return 0;
178 }
179
180 if (quot)
181 if (mpi_resize(quot, qsize) < 0)
182 goto nomem;
183
184 /* Read pointers here, when reallocation is finished. */
185 np = num->d;
186 dp = den->d;
187 rp = rem->d;
188
189 /* Optimize division by a single-limb divisor. */
190 if (dsize == 1) {
191 mpi_limb_t rlimb;
192 if (quot) {
193 qp = quot->d;
194 rlimb = mpihelp_divmod_1(qp, np, nsize, dp[0]);
195 qsize -= qp[qsize - 1] == 0;
196 quot->nlimbs = qsize;
197 quot->sign = sign_quotient;
198 } else
199 rlimb = mpihelp_mod_1(np, nsize, dp[0]);
200 rp[0] = rlimb;
201 rsize = rlimb != 0 ? 1 : 0;
202 rem->nlimbs = rsize;
203 rem->sign = sign_remainder;
204 return 0;
205 }
206
207 if (quot) {
208 qp = quot->d;
209 /* Make sure QP and NP point to different objects. Otherwise the
210 * numerator would be gradually overwritten by the quotient limbs. */
211 if (qp == np) { /* Copy NP object to temporary space. */
212 np = marker[markidx++] = mpi_alloc_limb_space(nsize);
213 if (!np)
214 goto nomem;
215 MPN_COPY(np, qp, nsize);
216 }
217 } else /* Put quotient at top of remainder. */
218 qp = rp + dsize;
219
220 count_leading_zeros(normalization_steps, dp[dsize - 1]);
221
222 /* Normalize the denominator, i.e. make its most significant bit set by
223 * shifting it NORMALIZATION_STEPS bits to the left. Also shift the
224 * numerator the same number of steps (to keep the quotient the same!).
225 */
226 if (normalization_steps) {
227 mpi_ptr_t tp;
228 mpi_limb_t nlimb;
229
230 /* Shift up the denominator setting the most significant bit of
231 * the most significant word. Use temporary storage not to clobber
232 * the original contents of the denominator. */
233 tp = marker[markidx++] = mpi_alloc_limb_space(dsize);
234 if (!tp)
235 goto nomem;
236 mpihelp_lshift(tp, dp, dsize, normalization_steps);
237 dp = tp;
238
239 /* Shift up the numerator, possibly introducing a new most
240 * significant word. Move the shifted numerator in the remainder
241 * meanwhile. */
242 nlimb = mpihelp_lshift(rp, np, nsize, normalization_steps);
243 if (nlimb) {
244 rp[nsize] = nlimb;
245 rsize = nsize + 1;
246 } else
247 rsize = nsize;
248 } else {
249 /* The denominator is already normalized, as required. Copy it to
250 * temporary space if it overlaps with the quotient or remainder. */
251 if (dp == rp || (quot && (dp == qp))) {
252 mpi_ptr_t tp;
253
254 tp = marker[markidx++] = mpi_alloc_limb_space(dsize);
255 if (!tp)
256 goto nomem;
257 MPN_COPY(tp, dp, dsize);
258 dp = tp;
259 }
260
261 /* Move the numerator to the remainder. */
262 if (rp != np)
263 MPN_COPY(rp, np, nsize);
264
265 rsize = nsize;
266 }
267
268 q_limb = mpihelp_divrem(qp, 0, rp, rsize, dp, dsize);
269
270 if (quot) {
271 qsize = rsize - dsize;
272 if (q_limb) {
273 qp[qsize] = q_limb;
274 qsize += 1;
275 }
276
277 quot->nlimbs = qsize;
278 quot->sign = sign_quotient;
279 }
280
281 rsize = dsize;
282 MPN_NORMALIZE(rp, rsize);
283
284 if (normalization_steps && rsize) {
285 mpihelp_rshift(rp, rp, rsize, normalization_steps);
286 rsize -= rp[rsize - 1] == 0 ? 1 : 0;
287 }
288
289 rem->nlimbs = rsize;
290 rem->sign = sign_remainder;
291
292 rc = 0;
293nomem:
294 while (markidx)
295 mpi_free_limb_space(marker[--markidx]);
296 return rc;
297}
298
299int mpi_tdiv_q_2exp(MPI w, MPI u, unsigned count)
300{
301 mpi_size_t usize, wsize;
302 mpi_size_t limb_cnt;
303
304 usize = u->nlimbs;
305 limb_cnt = count / BITS_PER_MPI_LIMB;
306 wsize = usize - limb_cnt;
307 if (limb_cnt >= usize)
308 w->nlimbs = 0;
309 else {
310 mpi_ptr_t wp;
311 mpi_ptr_t up;
312
313 if (RESIZE_IF_NEEDED(w, wsize) < 0)
314 return -ENOMEM;
315 wp = w->d;
316 up = u->d;
317
318 count %= BITS_PER_MPI_LIMB;
319 if (count) {
320 mpihelp_rshift(wp, up + limb_cnt, wsize, count);
321 wsize -= !wp[wsize - 1];
322 } else {
323 MPN_COPY_INCR(wp, up + limb_cnt, wsize);
324 }
325
326 w->nlimbs = wsize;
327 }
328 return 0;
329}
330
331/****************
332 * Check whether dividend is divisible by divisor
333 * (note: divisor must fit into a limb)
334 */
335int mpi_divisible_ui(MPI dividend, ulong divisor)
336{
337 return !mpihelp_mod_1(dividend->d, dividend->nlimbs, divisor);
338}
diff --git a/lib/mpi/mpi-gcd.c b/lib/mpi/mpi-gcd.c
new file mode 100644
index 000000000000..13c48aef9c4e
--- /dev/null
+++ b/lib/mpi/mpi-gcd.c
@@ -0,0 +1,59 @@
1/* mpi-gcd.c - MPI functions
2 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22
23/****************
24 * Find the greatest common divisor G of A and B.
25 * Return: true if this 1, false in all other cases
26 */
27int mpi_gcd(MPI g, const MPI xa, const MPI xb)
28{
29 MPI a = NULL, b = NULL;
30
31 if (mpi_copy(&a, xa) < 0)
32 goto nomem;
33
34 if (mpi_copy(&b, xb) < 0)
35 goto nomem;
36
37 /* TAOCP Vol II, 4.5.2, Algorithm A */
38 a->sign = 0;
39 b->sign = 0;
40 while (mpi_cmp_ui(b, 0)) {
41 if (mpi_fdiv_r(g, a, b) < 0) /* g used as temorary variable */
42 goto nomem;
43 if (mpi_set(a, b) < 0)
44 goto nomem;
45 if (mpi_set(b, g) < 0)
46 goto nomem;
47 }
48 if (mpi_set(g, a) < 0)
49 goto nomem;
50
51 mpi_free(a);
52 mpi_free(b);
53 return !mpi_cmp_ui(g, 1);
54
55nomem:
56 mpi_free(a);
57 mpi_free(b);
58 return -ENOMEM;
59}
diff --git a/lib/mpi/mpi-inline.c b/lib/mpi/mpi-inline.c
new file mode 100644
index 000000000000..654f68aeed8b
--- /dev/null
+++ b/lib/mpi/mpi-inline.c
@@ -0,0 +1,31 @@
1/* mpi-inline.c
2 * Copyright (C) 1999, 2000, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21/* put the inline functions as real functions into the lib */
22#define G10_MPI_INLINE_DECL
23
24#include "mpi-internal.h"
25
26/* always include the header becuase it is only
27 * included by mpi-internal if __GCC__ is defined but we
28 * need it here in all cases and the above definition of
29 * of the macro allows us to do so
30 */
31#include "mpi-inline.h"
diff --git a/lib/mpi/mpi-inline.h b/lib/mpi/mpi-inline.h
new file mode 100644
index 000000000000..e2b39852b30a
--- /dev/null
+++ b/lib/mpi/mpi-inline.h
@@ -0,0 +1,122 @@
1/* mpi-inline.h - Internal to the Multi Precision Integers
2 * Copyright (C) 1994, 1996, 1998, 1999 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#ifndef G10_MPI_INLINE_H
30#define G10_MPI_INLINE_H
31
32#ifndef G10_MPI_INLINE_DECL
33#define G10_MPI_INLINE_DECL extern inline
34#endif
35
36G10_MPI_INLINE_DECL mpi_limb_t
37mpihelp_add_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
38 mpi_size_t s1_size, mpi_limb_t s2_limb)
39{
40 mpi_limb_t x;
41
42 x = *s1_ptr++;
43 s2_limb += x;
44 *res_ptr++ = s2_limb;
45 if (s2_limb < x) { /* sum is less than the left operand: handle carry */
46 while (--s1_size) {
47 x = *s1_ptr++ + 1; /* add carry */
48 *res_ptr++ = x; /* and store */
49 if (x) /* not 0 (no overflow): we can stop */
50 goto leave;
51 }
52 return 1; /* return carry (size of s1 to small) */
53 }
54
55leave:
56 if (res_ptr != s1_ptr) { /* not the same variable */
57 mpi_size_t i; /* copy the rest */
58 for (i = 0; i < s1_size - 1; i++)
59 res_ptr[i] = s1_ptr[i];
60 }
61 return 0; /* no carry */
62}
63
64G10_MPI_INLINE_DECL mpi_limb_t
65mpihelp_add(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr, mpi_size_t s1_size,
66 mpi_ptr_t s2_ptr, mpi_size_t s2_size)
67{
68 mpi_limb_t cy = 0;
69
70 if (s2_size)
71 cy = mpihelp_add_n(res_ptr, s1_ptr, s2_ptr, s2_size);
72
73 if (s1_size - s2_size)
74 cy = mpihelp_add_1(res_ptr + s2_size, s1_ptr + s2_size,
75 s1_size - s2_size, cy);
76 return cy;
77}
78
79G10_MPI_INLINE_DECL mpi_limb_t
80mpihelp_sub_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
81 mpi_size_t s1_size, mpi_limb_t s2_limb)
82{
83 mpi_limb_t x;
84
85 x = *s1_ptr++;
86 s2_limb = x - s2_limb;
87 *res_ptr++ = s2_limb;
88 if (s2_limb > x) {
89 while (--s1_size) {
90 x = *s1_ptr++;
91 *res_ptr++ = x - 1;
92 if (x)
93 goto leave;
94 }
95 return 1;
96 }
97
98leave:
99 if (res_ptr != s1_ptr) {
100 mpi_size_t i;
101 for (i = 0; i < s1_size - 1; i++)
102 res_ptr[i] = s1_ptr[i];
103 }
104 return 0;
105}
106
107G10_MPI_INLINE_DECL mpi_limb_t
108mpihelp_sub(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr, mpi_size_t s1_size,
109 mpi_ptr_t s2_ptr, mpi_size_t s2_size)
110{
111 mpi_limb_t cy = 0;
112
113 if (s2_size)
114 cy = mpihelp_sub_n(res_ptr, s1_ptr, s2_ptr, s2_size);
115
116 if (s1_size - s2_size)
117 cy = mpihelp_sub_1(res_ptr + s2_size, s1_ptr + s2_size,
118 s1_size - s2_size, cy);
119 return cy;
120}
121
122#endif /*G10_MPI_INLINE_H */
diff --git a/lib/mpi/mpi-internal.h b/lib/mpi/mpi-internal.h
new file mode 100644
index 000000000000..77adcf6bc257
--- /dev/null
+++ b/lib/mpi/mpi-internal.h
@@ -0,0 +1,261 @@
1/* mpi-internal.h - Internal to the Multi Precision Integers
2 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
3 * Copyright (C) 1998, 2000 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#ifndef G10_MPI_INTERNAL_H
31#define G10_MPI_INTERNAL_H
32
33#include <linux/module.h>
34#include <linux/kernel.h>
35#include <linux/slab.h>
36#include <linux/string.h>
37#include <linux/mpi.h>
38#include <linux/errno.h>
39
40#define log_debug printk
41#define log_bug printk
42
43#define assert(x) \
44 do { \
45 if (!x) \
46 log_bug("failed assertion\n"); \
47 } while (0);
48
49/* If KARATSUBA_THRESHOLD is not already defined, define it to a
50 * value which is good on most machines. */
51
52/* tested 4, 16, 32 and 64, where 16 gave the best performance when
53 * checking a 768 and a 1024 bit ElGamal signature.
54 * (wk 22.12.97) */
55#ifndef KARATSUBA_THRESHOLD
56#define KARATSUBA_THRESHOLD 16
57#endif
58
59/* The code can't handle KARATSUBA_THRESHOLD smaller than 2. */
60#if KARATSUBA_THRESHOLD < 2
61#undef KARATSUBA_THRESHOLD
62#define KARATSUBA_THRESHOLD 2
63#endif
64
65typedef mpi_limb_t *mpi_ptr_t; /* pointer to a limb */
66typedef int mpi_size_t; /* (must be a signed type) */
67
68#define ABS(x) (x >= 0 ? x : -x)
69#define MIN(l, o) ((l) < (o) ? (l) : (o))
70#define MAX(h, i) ((h) > (i) ? (h) : (i))
71
72static inline int RESIZE_IF_NEEDED(MPI a, unsigned b)
73{
74 if (a->alloced < b)
75 return mpi_resize(a, b);
76 return 0;
77}
78
79/* Copy N limbs from S to D. */
80#define MPN_COPY(d, s, n) \
81 do { \
82 mpi_size_t _i; \
83 for (_i = 0; _i < (n); _i++) \
84 (d)[_i] = (s)[_i]; \
85 } while (0)
86
87#define MPN_COPY_INCR(d, s, n) \
88 do { \
89 mpi_size_t _i; \
90 for (_i = 0; _i < (n); _i++) \
91 (d)[_i] = (d)[_i]; \
92 } while (0)
93
94#define MPN_COPY_DECR(d, s, n) \
95 do { \
96 mpi_size_t _i; \
97 for (_i = (n)-1; _i >= 0; _i--) \
98 (d)[_i] = (s)[_i]; \
99 } while (0)
100
101/* Zero N limbs at D */
102#define MPN_ZERO(d, n) \
103 do { \
104 int _i; \
105 for (_i = 0; _i < (n); _i++) \
106 (d)[_i] = 0; \
107 } while (0)
108
109#define MPN_NORMALIZE(d, n) \
110 do { \
111 while ((n) > 0) { \
112 if ((d)[(n)-1]) \
113 break; \
114 (n)--; \
115 } \
116 } while (0)
117
118#define MPN_NORMALIZE_NOT_ZERO(d, n) \
119 do { \
120 for (;;) { \
121 if ((d)[(n)-1]) \
122 break; \
123 (n)--; \
124 } \
125 } while (0)
126
127#define MPN_MUL_N_RECURSE(prodp, up, vp, size, tspace) \
128 do { \
129 if ((size) < KARATSUBA_THRESHOLD) \
130 mul_n_basecase(prodp, up, vp, size); \
131 else \
132 mul_n(prodp, up, vp, size, tspace); \
133 } while (0);
134
135/* Divide the two-limb number in (NH,,NL) by D, with DI being the largest
136 * limb not larger than (2**(2*BITS_PER_MP_LIMB))/D - (2**BITS_PER_MP_LIMB).
137 * If this would yield overflow, DI should be the largest possible number
138 * (i.e., only ones). For correct operation, the most significant bit of D
139 * has to be set. Put the quotient in Q and the remainder in R.
140 */
141#define UDIV_QRNND_PREINV(q, r, nh, nl, d, di) \
142 do { \
143 mpi_limb_t _q, _ql, _r; \
144 mpi_limb_t _xh, _xl; \
145 umul_ppmm(_q, _ql, (nh), (di)); \
146 _q += (nh); /* DI is 2**BITS_PER_MPI_LIMB too small */ \
147 umul_ppmm(_xh, _xl, _q, (d)); \
148 sub_ddmmss(_xh, _r, (nh), (nl), _xh, _xl); \
149 if (_xh) { \
150 sub_ddmmss(_xh, _r, _xh, _r, 0, (d)); \
151 _q++; \
152 if (_xh) { \
153 sub_ddmmss(_xh, _r, _xh, _r, 0, (d)); \
154 _q++; \
155 } \
156 } \
157 if (_r >= (d)) { \
158 _r -= (d); \
159 _q++; \
160 } \
161 (r) = _r; \
162 (q) = _q; \
163 } while (0)
164
165/*-- mpiutil.c --*/
166mpi_ptr_t mpi_alloc_limb_space(unsigned nlimbs);
167void mpi_free_limb_space(mpi_ptr_t a);
168void mpi_assign_limb_space(MPI a, mpi_ptr_t ap, unsigned nlimbs);
169
170/*-- mpi-bit.c --*/
171void mpi_rshift_limbs(MPI a, unsigned int count);
172int mpi_lshift_limbs(MPI a, unsigned int count);
173
174/*-- mpihelp-add.c --*/
175mpi_limb_t mpihelp_add_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
176 mpi_size_t s1_size, mpi_limb_t s2_limb);
177mpi_limb_t mpihelp_add_n(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
178 mpi_ptr_t s2_ptr, mpi_size_t size);
179mpi_limb_t mpihelp_add(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr, mpi_size_t s1_size,
180 mpi_ptr_t s2_ptr, mpi_size_t s2_size);
181
182/*-- mpihelp-sub.c --*/
183mpi_limb_t mpihelp_sub_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
184 mpi_size_t s1_size, mpi_limb_t s2_limb);
185mpi_limb_t mpihelp_sub_n(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
186 mpi_ptr_t s2_ptr, mpi_size_t size);
187mpi_limb_t mpihelp_sub(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr, mpi_size_t s1_size,
188 mpi_ptr_t s2_ptr, mpi_size_t s2_size);
189
190/*-- mpihelp-cmp.c --*/
191int mpihelp_cmp(mpi_ptr_t op1_ptr, mpi_ptr_t op2_ptr, mpi_size_t size);
192
193/*-- mpihelp-mul.c --*/
194
195struct karatsuba_ctx {
196 struct karatsuba_ctx *next;
197 mpi_ptr_t tspace;
198 mpi_size_t tspace_size;
199 mpi_ptr_t tp;
200 mpi_size_t tp_size;
201};
202
203void mpihelp_release_karatsuba_ctx(struct karatsuba_ctx *ctx);
204
205mpi_limb_t mpihelp_addmul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
206 mpi_size_t s1_size, mpi_limb_t s2_limb);
207mpi_limb_t mpihelp_submul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
208 mpi_size_t s1_size, mpi_limb_t s2_limb);
209int mpihelp_mul_n(mpi_ptr_t prodp, mpi_ptr_t up, mpi_ptr_t vp, mpi_size_t size);
210int mpihelp_mul(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t usize,
211 mpi_ptr_t vp, mpi_size_t vsize, mpi_limb_t *_result);
212void mpih_sqr_n_basecase(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t size);
213void mpih_sqr_n(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t size,
214 mpi_ptr_t tspace);
215
216int mpihelp_mul_karatsuba_case(mpi_ptr_t prodp,
217 mpi_ptr_t up, mpi_size_t usize,
218 mpi_ptr_t vp, mpi_size_t vsize,
219 struct karatsuba_ctx *ctx);
220
221/*-- mpihelp-mul_1.c (or xxx/cpu/ *.S) --*/
222mpi_limb_t mpihelp_mul_1(mpi_ptr_t res_ptr, mpi_ptr_t s1_ptr,
223 mpi_size_t s1_size, mpi_limb_t s2_limb);
224
225/*-- mpihelp-div.c --*/
226mpi_limb_t mpihelp_mod_1(mpi_ptr_t dividend_ptr, mpi_size_t dividend_size,
227 mpi_limb_t divisor_limb);
228mpi_limb_t mpihelp_divrem(mpi_ptr_t qp, mpi_size_t qextra_limbs,
229 mpi_ptr_t np, mpi_size_t nsize,
230 mpi_ptr_t dp, mpi_size_t dsize);
231mpi_limb_t mpihelp_divmod_1(mpi_ptr_t quot_ptr,
232 mpi_ptr_t dividend_ptr, mpi_size_t dividend_size,
233 mpi_limb_t divisor_limb);
234
235/*-- mpihelp-shift.c --*/
236mpi_limb_t mpihelp_lshift(mpi_ptr_t wp, mpi_ptr_t up, mpi_size_t usize,
237 unsigned cnt);
238mpi_limb_t mpihelp_rshift(mpi_ptr_t wp, mpi_ptr_t up, mpi_size_t usize,
239 unsigned cnt);
240
241/* Define stuff for longlong.h. */
242#define W_TYPE_SIZE BITS_PER_MPI_LIMB
243typedef mpi_limb_t UWtype;
244typedef unsigned int UHWtype;
245#if defined(__GNUC__)
246typedef unsigned int UQItype __attribute__ ((mode(QI)));
247typedef int SItype __attribute__ ((mode(SI)));
248typedef unsigned int USItype __attribute__ ((mode(SI)));
249typedef int DItype __attribute__ ((mode(DI)));
250typedef unsigned int UDItype __attribute__ ((mode(DI)));
251#else
252typedef unsigned char UQItype;
253typedef long SItype;
254typedef unsigned long USItype;
255#endif
256
257#ifdef __GNUC__
258#include "mpi-inline.h"
259#endif
260
261#endif /*G10_MPI_INTERNAL_H */
diff --git a/lib/mpi/mpi-inv.c b/lib/mpi/mpi-inv.c
new file mode 100644
index 000000000000..0951f9847745
--- /dev/null
+++ b/lib/mpi/mpi-inv.c
@@ -0,0 +1,187 @@
1/* mpi-inv.c - MPI functions
2 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22
23/****************
24 * Calculate the multiplicative inverse X of A mod N
25 * That is: Find the solution x for
26 * 1 = (a*x) mod n
27 */
28int mpi_invm(MPI x, const MPI a, const MPI n)
29{
30 /* Extended Euclid's algorithm (See TAOPC Vol II, 4.5.2, Alg X)
31 * modified according to Michael Penk's solution for Exercice 35
32 * with further enhancement */
33 MPI u = NULL, v = NULL;
34 MPI u1 = NULL, u2 = NULL, u3 = NULL;
35 MPI v1 = NULL, v2 = NULL, v3 = NULL;
36 MPI t1 = NULL, t2 = NULL, t3 = NULL;
37 unsigned k;
38 int sign;
39 int odd = 0;
40 int rc = -ENOMEM;
41
42 if (mpi_copy(&u, a) < 0)
43 goto cleanup;
44 if (mpi_copy(&v, n) < 0)
45 goto cleanup;
46
47 for (k = 0; !mpi_test_bit(u, 0) && !mpi_test_bit(v, 0); k++) {
48 if (mpi_rshift(u, u, 1) < 0)
49 goto cleanup;
50 if (mpi_rshift(v, v, 1) < 0)
51 goto cleanup;
52 }
53 odd = mpi_test_bit(v, 0);
54
55 u1 = mpi_alloc_set_ui(1);
56 if (!u1)
57 goto cleanup;
58 if (!odd) {
59 u2 = mpi_alloc_set_ui(0);
60 if (!u2)
61 goto cleanup;
62 }
63 if (mpi_copy(&u3, u) < 0)
64 goto cleanup;
65 if (mpi_copy(&v1, v) < 0)
66 goto cleanup;
67 if (!odd) {
68 v2 = mpi_alloc(mpi_get_nlimbs(u));
69 if (!v2)
70 goto cleanup;
71 if (mpi_sub(v2, u1, u) < 0)
72 goto cleanup; /* U is used as const 1 */
73 }
74 if (mpi_copy(&v3, v) < 0)
75 goto cleanup;
76 if (mpi_test_bit(u, 0)) { /* u is odd */
77 t1 = mpi_alloc_set_ui(0);
78 if (!t1)
79 goto cleanup;
80 if (!odd) {
81 t2 = mpi_alloc_set_ui(1);
82 if (!t2)
83 goto cleanup;
84 t2->sign = 1;
85 }
86 if (mpi_copy(&t3, v) < 0)
87 goto cleanup;
88 t3->sign = !t3->sign;
89 goto Y4;
90 } else {
91 t1 = mpi_alloc_set_ui(1);
92 if (!t1)
93 goto cleanup;
94 if (!odd) {
95 t2 = mpi_alloc_set_ui(0);
96 if (!t2)
97 goto cleanup;
98 }
99 if (mpi_copy(&t3, u) < 0)
100 goto cleanup;
101 }
102 do {
103 do {
104 if (!odd) {
105 if (mpi_test_bit(t1, 0) || mpi_test_bit(t2, 0)) { /* one is odd */
106 if (mpi_add(t1, t1, v) < 0)
107 goto cleanup;
108 if (mpi_sub(t2, t2, u) < 0)
109 goto cleanup;
110 }
111 if (mpi_rshift(t1, t1, 1) < 0)
112 goto cleanup;
113 if (mpi_rshift(t2, t2, 1) < 0)
114 goto cleanup;
115 if (mpi_rshift(t3, t3, 1) < 0)
116 goto cleanup;
117 } else {
118 if (mpi_test_bit(t1, 0))
119 if (mpi_add(t1, t1, v) < 0)
120 goto cleanup;
121 if (mpi_rshift(t1, t1, 1) < 0)
122 goto cleanup;
123 if (mpi_rshift(t3, t3, 1) < 0)
124 goto cleanup;
125 }
126Y4:
127 ;
128 } while (!mpi_test_bit(t3, 0)); /* while t3 is even */
129
130 if (!t3->sign) {
131 if (mpi_set(u1, t1) < 0)
132 goto cleanup;
133 if (!odd)
134 if (mpi_set(u2, t2) < 0)
135 goto cleanup;
136 if (mpi_set(u3, t3) < 0)
137 goto cleanup;
138 } else {
139 if (mpi_sub(v1, v, t1) < 0)
140 goto cleanup;
141 sign = u->sign;
142 u->sign = !u->sign;
143 if (!odd)
144 if (mpi_sub(v2, u, t2) < 0)
145 goto cleanup;
146 u->sign = sign;
147 sign = t3->sign;
148 t3->sign = !t3->sign;
149 if (mpi_set(v3, t3) < 0)
150 goto cleanup;
151 t3->sign = sign;
152 }
153 if (mpi_sub(t1, u1, v1) < 0)
154 goto cleanup;
155 if (!odd)
156 if (mpi_sub(t2, u2, v2) < 0)
157 goto cleanup;
158 if (mpi_sub(t3, u3, v3) < 0)
159 goto cleanup;
160 if (t1->sign) {
161 if (mpi_add(t1, t1, v) < 0)
162 goto cleanup;
163 if (!odd)
164 if (mpi_sub(t2, t2, u) < 0)
165 goto cleanup;
166 }
167 } while (mpi_cmp_ui(t3, 0)); /* while t3 != 0 */
168 /* mpi_lshift( u3, k ); */
169 rc = mpi_set(x, u1);
170
171cleanup:
172 mpi_free(u1);
173 mpi_free(v1);
174 mpi_free(t1);
175 if (!odd) {
176 mpi_free(u2);
177 mpi_free(v2);
178 mpi_free(t2);
179 }
180 mpi_free(u3);
181 mpi_free(v3);
182 mpi_free(t3);
183
184 mpi_free(u);
185 mpi_free(v);
186 return rc;
187}
diff --git a/lib/mpi/mpi-mpow.c b/lib/mpi/mpi-mpow.c
new file mode 100644
index 000000000000..7328d0d6c748
--- /dev/null
+++ b/lib/mpi/mpi-mpow.c
@@ -0,0 +1,134 @@
1/* mpi-mpow.c - MPI functions
2 * Copyright (C) 1998, 1999, 2000 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22#include "longlong.h"
23
24static int build_index(const MPI *exparray, int k, int i, int t)
25{
26 int j, bitno;
27 int index = 0;
28
29 bitno = t - i;
30 for (j = k - 1; j >= 0; j--) {
31 index <<= 1;
32 if (mpi_test_bit(exparray[j], bitno))
33 index |= 1;
34 }
35 return index;
36}
37
38/****************
39 * RES = (BASE[0] ^ EXP[0]) * (BASE[1] ^ EXP[1]) * ... * mod M
40 */
41int mpi_mulpowm(MPI res, MPI *basearray, MPI *exparray, MPI m)
42{
43 int rc = -ENOMEM;
44 int k; /* number of elements */
45 int t; /* bit size of largest exponent */
46 int i, j, idx;
47 MPI *G = NULL; /* table with precomputed values of size 2^k */
48 MPI tmp = NULL;
49
50 for (k = 0; basearray[k]; k++)
51 ;
52 if (!k) {
53 pr_emerg("mpi_mulpowm: assert(k) failed\n");
54 BUG();
55 }
56 for (t = 0, i = 0; (tmp = exparray[i]); i++) {
57 j = mpi_get_nbits(tmp);
58 if (j > t)
59 t = j;
60 }
61 if (i != k) {
62 pr_emerg("mpi_mulpowm: assert(i==k) failed\n");
63 BUG();
64 }
65 if (!t) {
66 pr_emerg("mpi_mulpowm: assert(t) failed\n");
67 BUG();
68 }
69 if (k >= 10) {
70 pr_emerg("mpi_mulpowm: assert(k<10) failed\n");
71 BUG();
72 }
73
74 G = kzalloc((1 << k) * sizeof *G, GFP_KERNEL);
75 if (!G)
76 goto err_out;
77
78 /* and calculate */
79 tmp = mpi_alloc(mpi_get_nlimbs(m) + 1);
80 if (!tmp)
81 goto nomem;
82 if (mpi_set_ui(res, 1) < 0)
83 goto nomem;
84 for (i = 1; i <= t; i++) {
85 if (mpi_mulm(tmp, res, res, m) < 0)
86 goto nomem;
87 idx = build_index(exparray, k, i, t);
88 if (!(idx >= 0 && idx < (1 << k))) {
89 pr_emerg("mpi_mulpowm: assert(idx >= 0 && idx < (1<<k)) failed\n");
90 BUG();
91 }
92 if (!G[idx]) {
93 if (!idx) {
94 G[0] = mpi_alloc_set_ui(1);
95 if (!G[0])
96 goto nomem;
97 } else {
98 for (j = 0; j < k; j++) {
99 if ((idx & (1 << j))) {
100 if (!G[idx]) {
101 if (mpi_copy
102 (&G[idx],
103 basearray[j]) < 0)
104 goto nomem;
105 } else {
106 if (mpi_mulm
107 (G[idx], G[idx],
108 basearray[j],
109 m) < 0)
110 goto nomem;
111 }
112 }
113 }
114 if (!G[idx]) {
115 G[idx] = mpi_alloc(0);
116 if (!G[idx])
117 goto nomem;
118 }
119 }
120 }
121 if (mpi_mulm(res, tmp, G[idx], m) < 0)
122 goto nomem;
123 }
124
125 rc = 0;
126nomem:
127 /* cleanup */
128 mpi_free(tmp);
129 for (i = 0; i < (1 << k); i++)
130 mpi_free(G[i]);
131 kfree(G);
132err_out:
133 return rc;
134}
diff --git a/lib/mpi/mpi-mul.c b/lib/mpi/mpi-mul.c
new file mode 100644
index 000000000000..1f3219e27292
--- /dev/null
+++ b/lib/mpi/mpi-mul.c
@@ -0,0 +1,194 @@
1/* mpi-mul.c - MPI functions
2 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
3 * Copyright (C) 1998, 2001 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31
32int mpi_mul_ui(MPI prod, MPI mult, unsigned long small_mult)
33{
34 mpi_size_t size, prod_size;
35 mpi_ptr_t prod_ptr;
36 mpi_limb_t cy;
37 int sign;
38
39 size = mult->nlimbs;
40 sign = mult->sign;
41
42 if (!size || !small_mult) {
43 prod->nlimbs = 0;
44 prod->sign = 0;
45 return 0;
46 }
47
48 prod_size = size + 1;
49 if (prod->alloced < prod_size)
50 if (mpi_resize(prod, prod_size) < 0)
51 return -ENOMEM;
52 prod_ptr = prod->d;
53
54 cy = mpihelp_mul_1(prod_ptr, mult->d, size, (mpi_limb_t) small_mult);
55 if (cy)
56 prod_ptr[size++] = cy;
57 prod->nlimbs = size;
58 prod->sign = sign;
59 return 0;
60}
61
62int mpi_mul_2exp(MPI w, MPI u, unsigned long cnt)
63{
64 mpi_size_t usize, wsize, limb_cnt;
65 mpi_ptr_t wp;
66 mpi_limb_t wlimb;
67 int usign, wsign;
68
69 usize = u->nlimbs;
70 usign = u->sign;
71
72 if (!usize) {
73 w->nlimbs = 0;
74 w->sign = 0;
75 return 0;
76 }
77
78 limb_cnt = cnt / BITS_PER_MPI_LIMB;
79 wsize = usize + limb_cnt + 1;
80 if (w->alloced < wsize)
81 if (mpi_resize(w, wsize) < 0)
82 return -ENOMEM;
83 wp = w->d;
84 wsize = usize + limb_cnt;
85 wsign = usign;
86
87 cnt %= BITS_PER_MPI_LIMB;
88 if (cnt) {
89 wlimb = mpihelp_lshift(wp + limb_cnt, u->d, usize, cnt);
90 if (wlimb) {
91 wp[wsize] = wlimb;
92 wsize++;
93 }
94 } else {
95 MPN_COPY_DECR(wp + limb_cnt, u->d, usize);
96 }
97
98 /* Zero all whole limbs at low end. Do it here and not before calling
99 * mpn_lshift, not to lose for U == W. */
100 MPN_ZERO(wp, limb_cnt);
101
102 w->nlimbs = wsize;
103 w->sign = wsign;
104 return 0;
105}
106
107int mpi_mul(MPI w, MPI u, MPI v)
108{
109 int rc = -ENOMEM;
110 mpi_size_t usize, vsize, wsize;
111 mpi_ptr_t up, vp, wp;
112 mpi_limb_t cy;
113 int usign, vsign, sign_product;
114 int assign_wp = 0;
115 mpi_ptr_t tmp_limb = NULL;
116
117 if (u->nlimbs < v->nlimbs) { /* Swap U and V. */
118 usize = v->nlimbs;
119 usign = v->sign;
120 up = v->d;
121 vsize = u->nlimbs;
122 vsign = u->sign;
123 vp = u->d;
124 } else {
125 usize = u->nlimbs;
126 usign = u->sign;
127 up = u->d;
128 vsize = v->nlimbs;
129 vsign = v->sign;
130 vp = v->d;
131 }
132 sign_product = usign ^ vsign;
133 wp = w->d;
134
135 /* Ensure W has space enough to store the result. */
136 wsize = usize + vsize;
137 if (w->alloced < (size_t) wsize) {
138 if (wp == up || wp == vp) {
139 wp = mpi_alloc_limb_space(wsize);
140 if (!wp)
141 goto nomem;
142 assign_wp = 1;
143 } else {
144 if (mpi_resize(w, wsize) < 0)
145 goto nomem;
146 wp = w->d;
147 }
148 } else { /* Make U and V not overlap with W. */
149 if (wp == up) {
150 /* W and U are identical. Allocate temporary space for U. */
151 up = tmp_limb = mpi_alloc_limb_space(usize);
152 if (!up)
153 goto nomem;
154 /* Is V identical too? Keep it identical with U. */
155 if (wp == vp)
156 vp = up;
157 /* Copy to the temporary space. */
158 MPN_COPY(up, wp, usize);
159 } else if (wp == vp) {
160 /* W and V are identical. Allocate temporary space for V. */
161 vp = tmp_limb = mpi_alloc_limb_space(vsize);
162 if (!vp)
163 goto nomem;
164 /* Copy to the temporary space. */
165 MPN_COPY(vp, wp, vsize);
166 }
167 }
168
169 if (!vsize)
170 wsize = 0;
171 else {
172 if (mpihelp_mul(wp, up, usize, vp, vsize, &cy) < 0)
173 goto nomem;
174 wsize -= cy ? 0 : 1;
175 }
176
177 if (assign_wp)
178 mpi_assign_limb_space(w, wp, wsize);
179
180 w->nlimbs = wsize;
181 w->sign = sign_product;
182 rc = 0;
183nomem:
184 if (tmp_limb)
185 mpi_free_limb_space(tmp_limb);
186 return rc;
187}
188
189int mpi_mulm(MPI w, MPI u, MPI v, MPI m)
190{
191 if (mpi_mul(w, u, v) < 0)
192 return -ENOMEM;
193 return mpi_fdiv_r(w, w, m);
194}
diff --git a/lib/mpi/mpi-pow.c b/lib/mpi/mpi-pow.c
new file mode 100644
index 000000000000..67f3e79af914
--- /dev/null
+++ b/lib/mpi/mpi-pow.c
@@ -0,0 +1,323 @@
1/* mpi-pow.c - MPI functions
2 * Copyright (C) 1994, 1996, 1998, 2000 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 *
20 * Note: This code is heavily based on the GNU MP Library.
21 * Actually it's the same code with only minor changes in the
22 * way the data is stored; this is to support the abstraction
23 * of an optional secure memory allocation which may be used
24 * to avoid revealing of sensitive data due to paging etc.
25 * The GNU MP Library itself is published under the LGPL;
26 * however I decided to publish this code under the plain GPL.
27 */
28
29#include <linux/string.h>
30#include "mpi-internal.h"
31#include "longlong.h"
32
33/****************
34 * RES = BASE ^ EXP mod MOD
35 */
36int mpi_powm(MPI res, MPI base, MPI exp, MPI mod)
37{
38 mpi_ptr_t mp_marker = NULL, bp_marker = NULL, ep_marker = NULL;
39 mpi_ptr_t xp_marker = NULL;
40 mpi_ptr_t tspace = NULL;
41 mpi_ptr_t rp, ep, mp, bp;
42 mpi_size_t esize, msize, bsize, rsize;
43 int esign, msign, bsign, rsign;
44 mpi_size_t size;
45 int mod_shift_cnt;
46 int negative_result;
47 int assign_rp = 0;
48 mpi_size_t tsize = 0; /* to avoid compiler warning */
49 /* fixme: we should check that the warning is void */
50 int rc = -ENOMEM;
51
52 esize = exp->nlimbs;
53 msize = mod->nlimbs;
54 size = 2 * msize;
55 esign = exp->sign;
56 msign = mod->sign;
57
58 rp = res->d;
59 ep = exp->d;
60
61 if (!msize)
62 return -EINVAL;
63
64 if (!esize) {
65 /* Exponent is zero, result is 1 mod MOD, i.e., 1 or 0
66 * depending on if MOD equals 1. */
67 rp[0] = 1;
68 res->nlimbs = (msize == 1 && mod->d[0] == 1) ? 0 : 1;
69 res->sign = 0;
70 goto leave;
71 }
72
73 /* Normalize MOD (i.e. make its most significant bit set) as required by
74 * mpn_divrem. This will make the intermediate values in the calculation
75 * slightly larger, but the correct result is obtained after a final
76 * reduction using the original MOD value. */
77 mp = mp_marker = mpi_alloc_limb_space(msize);
78 if (!mp)
79 goto enomem;
80 count_leading_zeros(mod_shift_cnt, mod->d[msize - 1]);
81 if (mod_shift_cnt)
82 mpihelp_lshift(mp, mod->d, msize, mod_shift_cnt);
83 else
84 MPN_COPY(mp, mod->d, msize);
85
86 bsize = base->nlimbs;
87 bsign = base->sign;
88 if (bsize > msize) { /* The base is larger than the module. Reduce it. */
89 /* Allocate (BSIZE + 1) with space for remainder and quotient.
90 * (The quotient is (bsize - msize + 1) limbs.) */
91 bp = bp_marker = mpi_alloc_limb_space(bsize + 1);
92 if (!bp)
93 goto enomem;
94 MPN_COPY(bp, base->d, bsize);
95 /* We don't care about the quotient, store it above the remainder,
96 * at BP + MSIZE. */
97 mpihelp_divrem(bp + msize, 0, bp, bsize, mp, msize);
98 bsize = msize;
99 /* Canonicalize the base, since we are going to multiply with it
100 * quite a few times. */
101 MPN_NORMALIZE(bp, bsize);
102 } else
103 bp = base->d;
104
105 if (!bsize) {
106 res->nlimbs = 0;
107 res->sign = 0;
108 goto leave;
109 }
110
111 if (res->alloced < size) {
112 /* We have to allocate more space for RES. If any of the input
113 * parameters are identical to RES, defer deallocation of the old
114 * space. */
115 if (rp == ep || rp == mp || rp == bp) {
116 rp = mpi_alloc_limb_space(size);
117 if (!rp)
118 goto enomem;
119 assign_rp = 1;
120 } else {
121 if (mpi_resize(res, size) < 0)
122 goto enomem;
123 rp = res->d;
124 }
125 } else { /* Make BASE, EXP and MOD not overlap with RES. */
126 if (rp == bp) {
127 /* RES and BASE are identical. Allocate temp. space for BASE. */
128 BUG_ON(bp_marker);
129 bp = bp_marker = mpi_alloc_limb_space(bsize);
130 if (!bp)
131 goto enomem;
132 MPN_COPY(bp, rp, bsize);
133 }
134 if (rp == ep) {
135 /* RES and EXP are identical. Allocate temp. space for EXP. */
136 ep = ep_marker = mpi_alloc_limb_space(esize);
137 if (!ep)
138 goto enomem;
139 MPN_COPY(ep, rp, esize);
140 }
141 if (rp == mp) {
142 /* RES and MOD are identical. Allocate temporary space for MOD. */
143 BUG_ON(mp_marker);
144 mp = mp_marker = mpi_alloc_limb_space(msize);
145 if (!mp)
146 goto enomem;
147 MPN_COPY(mp, rp, msize);
148 }
149 }
150
151 MPN_COPY(rp, bp, bsize);
152 rsize = bsize;
153 rsign = bsign;
154
155 {
156 mpi_size_t i;
157 mpi_ptr_t xp;
158 int c;
159 mpi_limb_t e;
160 mpi_limb_t carry_limb;
161 struct karatsuba_ctx karactx;
162
163 xp = xp_marker = mpi_alloc_limb_space(2 * (msize + 1));
164 if (!xp)
165 goto enomem;
166
167 memset(&karactx, 0, sizeof karactx);
168 negative_result = (ep[0] & 1) && base->sign;
169
170 i = esize - 1;
171 e = ep[i];
172 count_leading_zeros(c, e);
173 e = (e << c) << 1; /* shift the exp bits to the left, lose msb */
174 c = BITS_PER_MPI_LIMB - 1 - c;
175
176 /* Main loop.
177 *
178 * Make the result be pointed to alternately by XP and RP. This
179 * helps us avoid block copying, which would otherwise be necessary
180 * with the overlap restrictions of mpihelp_divmod. With 50% probability
181 * the result after this loop will be in the area originally pointed
182 * by RP (==RES->d), and with 50% probability in the area originally
183 * pointed to by XP.
184 */
185
186 for (;;) {
187 while (c) {
188 mpi_ptr_t tp;
189 mpi_size_t xsize;
190
191 /*if (mpihelp_mul_n(xp, rp, rp, rsize) < 0) goto enomem */
192 if (rsize < KARATSUBA_THRESHOLD)
193 mpih_sqr_n_basecase(xp, rp, rsize);
194 else {
195 if (!tspace) {
196 tsize = 2 * rsize;
197 tspace =
198 mpi_alloc_limb_space(tsize);
199 if (!tspace)
200 goto enomem;
201 } else if (tsize < (2 * rsize)) {
202 mpi_free_limb_space(tspace);
203 tsize = 2 * rsize;
204 tspace =
205 mpi_alloc_limb_space(tsize);
206 if (!tspace)
207 goto enomem;
208 }
209 mpih_sqr_n(xp, rp, rsize, tspace);
210 }
211
212 xsize = 2 * rsize;
213 if (xsize > msize) {
214 mpihelp_divrem(xp + msize, 0, xp, xsize,
215 mp, msize);
216 xsize = msize;
217 }
218
219 tp = rp;
220 rp = xp;
221 xp = tp;
222 rsize = xsize;
223
224 if ((mpi_limb_signed_t) e < 0) {
225 /*mpihelp_mul( xp, rp, rsize, bp, bsize ); */
226 if (bsize < KARATSUBA_THRESHOLD) {
227 mpi_limb_t tmp;
228 if (mpihelp_mul
229 (xp, rp, rsize, bp, bsize,
230 &tmp) < 0)
231 goto enomem;
232 } else {
233 if (mpihelp_mul_karatsuba_case
234 (xp, rp, rsize, bp, bsize,
235 &karactx) < 0)
236 goto enomem;
237 }
238
239 xsize = rsize + bsize;
240 if (xsize > msize) {
241 mpihelp_divrem(xp + msize, 0,
242 xp, xsize, mp,
243 msize);
244 xsize = msize;
245 }
246
247 tp = rp;
248 rp = xp;
249 xp = tp;
250 rsize = xsize;
251 }
252 e <<= 1;
253 c--;
254 }
255
256 i--;
257 if (i < 0)
258 break;
259 e = ep[i];
260 c = BITS_PER_MPI_LIMB;
261 }
262
263 /* We shifted MOD, the modulo reduction argument, left MOD_SHIFT_CNT
264 * steps. Adjust the result by reducing it with the original MOD.
265 *
266 * Also make sure the result is put in RES->d (where it already
267 * might be, see above).
268 */
269 if (mod_shift_cnt) {
270 carry_limb =
271 mpihelp_lshift(res->d, rp, rsize, mod_shift_cnt);
272 rp = res->d;
273 if (carry_limb) {
274 rp[rsize] = carry_limb;
275 rsize++;
276 }
277 } else {
278 MPN_COPY(res->d, rp, rsize);
279 rp = res->d;
280 }
281
282 if (rsize >= msize) {
283 mpihelp_divrem(rp + msize, 0, rp, rsize, mp, msize);
284 rsize = msize;
285 }
286
287 /* Remove any leading zero words from the result. */
288 if (mod_shift_cnt)
289 mpihelp_rshift(rp, rp, rsize, mod_shift_cnt);
290 MPN_NORMALIZE(rp, rsize);
291
292 mpihelp_release_karatsuba_ctx(&karactx);
293 }
294
295 if (negative_result && rsize) {
296 if (mod_shift_cnt)
297 mpihelp_rshift(mp, mp, msize, mod_shift_cnt);
298 mpihelp_sub(rp, mp, msize, rp, rsize);
299 rsize = msize;
300 rsign = msign;
301 MPN_NORMALIZE(rp, rsize);
302 }
303 res->nlimbs = rsize;
304 res->sign = rsign;
305
306leave:
307 rc = 0;
308enomem:
309 if (assign_rp)
310 mpi_assign_limb_space(res, rp, size);
311 if (mp_marker)
312 mpi_free_limb_space(mp_marker);
313 if (bp_marker)
314 mpi_free_limb_space(bp_marker);
315 if (ep_marker)
316 mpi_free_limb_space(ep_marker);
317 if (xp_marker)
318 mpi_free_limb_space(xp_marker);
319 if (tspace)
320 mpi_free_limb_space(tspace);
321 return rc;
322}
323EXPORT_SYMBOL_GPL(mpi_powm);
diff --git a/lib/mpi/mpi-scan.c b/lib/mpi/mpi-scan.c
new file mode 100644
index 000000000000..b2da5ad96199
--- /dev/null
+++ b/lib/mpi/mpi-scan.c
@@ -0,0 +1,136 @@
1/* mpi-scan.c - MPI functions
2 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22#include "longlong.h"
23
24/****************
25 * Scan through an mpi and return byte for byte. a -1 is returned to indicate
26 * the end of the mpi. Scanning is done from the lsb to the msb, returned
27 * values are in the range of 0 .. 255.
28 *
29 * FIXME: This code is VERY ugly!
30 */
31int mpi_getbyte(const MPI a, unsigned idx)
32{
33 int i, j;
34 unsigned n;
35 mpi_ptr_t ap;
36 mpi_limb_t limb;
37
38 ap = a->d;
39 for (n = 0, i = 0; i < a->nlimbs; i++) {
40 limb = ap[i];
41 for (j = 0; j < BYTES_PER_MPI_LIMB; j++, n++)
42 if (n == idx)
43 return (limb >> j * 8) & 0xff;
44 }
45 return -1;
46}
47
48/****************
49 * Put a value at position IDX into A. idx counts from lsb to msb
50 */
51void mpi_putbyte(MPI a, unsigned idx, int xc)
52{
53 int i, j;
54 unsigned n;
55 mpi_ptr_t ap;
56 mpi_limb_t limb, c;
57
58 c = xc & 0xff;
59 ap = a->d;
60 for (n = 0, i = 0; i < a->alloced; i++) {
61 limb = ap[i];
62 for (j = 0; j < BYTES_PER_MPI_LIMB; j++, n++)
63 if (n == idx) {
64#if BYTES_PER_MPI_LIMB == 4
65 if (j == 0)
66 limb = (limb & 0xffffff00) | c;
67 else if (j == 1)
68 limb = (limb & 0xffff00ff) | (c << 8);
69 else if (j == 2)
70 limb = (limb & 0xff00ffff) | (c << 16);
71 else
72 limb = (limb & 0x00ffffff) | (c << 24);
73#elif BYTES_PER_MPI_LIMB == 8
74 if (j == 0)
75 limb = (limb & 0xffffffffffffff00) | c;
76 else if (j == 1)
77 limb =
78 (limb & 0xffffffffffff00ff) | (c <<
79 8);
80 else if (j == 2)
81 limb =
82 (limb & 0xffffffffff00ffff) | (c <<
83 16);
84 else if (j == 3)
85 limb =
86 (limb & 0xffffffff00ffffff) | (c <<
87 24);
88 else if (j == 4)
89 limb =
90 (limb & 0xffffff00ffffffff) | (c <<
91 32);
92 else if (j == 5)
93 limb =
94 (limb & 0xffff00ffffffffff) | (c <<
95 40);
96 else if (j == 6)
97 limb =
98 (limb & 0xff00ffffffffffff) | (c <<
99 48);
100 else
101 limb =
102 (limb & 0x00ffffffffffffff) | (c <<
103 56);
104#else
105#error please enhance this function, its ugly - i know.
106#endif
107 if (a->nlimbs <= i)
108 a->nlimbs = i + 1;
109 ap[i] = limb;
110 return;
111 }
112 }
113 log_bug("index out of range\n");
114}
115
116/****************
117 * Count the number of zerobits at the low end of A
118 */
119unsigned mpi_trailing_zeros(const MPI a)
120{
121 unsigned n, count = 0;
122
123 for (n = 0; n < a->nlimbs; n++) {
124 if (a->d[n]) {
125 unsigned nn;
126 mpi_limb_t alimb = a->d[n];
127
128 count_trailing_zeros(nn, alimb);
129 count += nn;
130 break;
131 }
132 count += BITS_PER_MPI_LIMB;
133 }
134 return count;
135
136}
diff --git a/lib/mpi/mpicoder.c b/lib/mpi/mpicoder.c
new file mode 100644
index 000000000000..f26b41fcb48c
--- /dev/null
+++ b/lib/mpi/mpicoder.c
@@ -0,0 +1,280 @@
1/* mpicoder.c - Coder for the external representation of MPIs
2 * Copyright (C) 1998, 1999 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22
23#define MAX_EXTERN_MPI_BITS 16384
24
25MPI mpi_read_from_buffer(const void *xbuffer, unsigned *ret_nread)
26{
27 const uint8_t *buffer = xbuffer;
28 int i, j;
29 unsigned nbits, nbytes, nlimbs, nread = 0;
30 mpi_limb_t a;
31 MPI val = NULL;
32
33 if (*ret_nread < 2)
34 goto leave;
35 nbits = buffer[0] << 8 | buffer[1];
36
37 if (nbits > MAX_EXTERN_MPI_BITS) {
38 pr_info("MPI: mpi too large (%u bits)\n", nbits);
39 goto leave;
40 }
41 buffer += 2;
42 nread = 2;
43
44 nbytes = (nbits + 7) / 8;
45 nlimbs = (nbytes + BYTES_PER_MPI_LIMB - 1) / BYTES_PER_MPI_LIMB;
46 val = mpi_alloc(nlimbs);
47 if (!val)
48 return NULL;
49 i = BYTES_PER_MPI_LIMB - nbytes % BYTES_PER_MPI_LIMB;
50 i %= BYTES_PER_MPI_LIMB;
51 val->nbits = nbits;
52 j = val->nlimbs = nlimbs;
53 val->sign = 0;
54 for (; j > 0; j--) {
55 a = 0;
56 for (; i < BYTES_PER_MPI_LIMB; i++) {
57 if (++nread > *ret_nread) {
58 printk
59 ("MPI: mpi larger than buffer nread=%d ret_nread=%d\n",
60 nread, *ret_nread);
61 goto leave;
62 }
63 a <<= 8;
64 a |= *buffer++;
65 }
66 i = 0;
67 val->d[j - 1] = a;
68 }
69
70leave:
71 *ret_nread = nread;
72 return val;
73}
74EXPORT_SYMBOL_GPL(mpi_read_from_buffer);
75
76/****************
77 * Make an mpi from a character string.
78 */
79int mpi_fromstr(MPI val, const char *str)
80{
81 int hexmode = 0, sign = 0, prepend_zero = 0, i, j, c, c1, c2;
82 unsigned nbits, nbytes, nlimbs;
83 mpi_limb_t a;
84
85 if (*str == '-') {
86 sign = 1;
87 str++;
88 }
89 if (*str == '0' && str[1] == 'x')
90 hexmode = 1;
91 else
92 return -EINVAL; /* other bases are not yet supported */
93 str += 2;
94
95 nbits = strlen(str) * 4;
96 if (nbits % 8)
97 prepend_zero = 1;
98 nbytes = (nbits + 7) / 8;
99 nlimbs = (nbytes + BYTES_PER_MPI_LIMB - 1) / BYTES_PER_MPI_LIMB;
100 if (val->alloced < nlimbs)
101 if (!mpi_resize(val, nlimbs))
102 return -ENOMEM;
103 i = BYTES_PER_MPI_LIMB - nbytes % BYTES_PER_MPI_LIMB;
104 i %= BYTES_PER_MPI_LIMB;
105 j = val->nlimbs = nlimbs;
106 val->sign = sign;
107 for (; j > 0; j--) {
108 a = 0;
109 for (; i < BYTES_PER_MPI_LIMB; i++) {
110 if (prepend_zero) {
111 c1 = '0';
112 prepend_zero = 0;
113 } else
114 c1 = *str++;
115 assert(c1);
116 c2 = *str++;
117 assert(c2);
118 if (c1 >= '0' && c1 <= '9')
119 c = c1 - '0';
120 else if (c1 >= 'a' && c1 <= 'f')
121 c = c1 - 'a' + 10;
122 else if (c1 >= 'A' && c1 <= 'F')
123 c = c1 - 'A' + 10;
124 else {
125 mpi_clear(val);
126 return 1;
127 }
128 c <<= 4;
129 if (c2 >= '0' && c2 <= '9')
130 c |= c2 - '0';
131 else if (c2 >= 'a' && c2 <= 'f')
132 c |= c2 - 'a' + 10;
133 else if (c2 >= 'A' && c2 <= 'F')
134 c |= c2 - 'A' + 10;
135 else {
136 mpi_clear(val);
137 return 1;
138 }
139 a <<= 8;
140 a |= c;
141 }
142 i = 0;
143
144 val->d[j - 1] = a;
145 }
146
147 return 0;
148}
149EXPORT_SYMBOL_GPL(mpi_fromstr);
150
151/****************
152 * Return an allocated buffer with the MPI (msb first).
153 * NBYTES receives the length of this buffer. Caller must free the
154 * return string (This function does return a 0 byte buffer with NBYTES
155 * set to zero if the value of A is zero. If sign is not NULL, it will
156 * be set to the sign of the A.
157 */
158void *mpi_get_buffer(MPI a, unsigned *nbytes, int *sign)
159{
160 uint8_t *p, *buffer;
161 mpi_limb_t alimb;
162 int i;
163 unsigned int n;
164
165 if (sign)
166 *sign = a->sign;
167 *nbytes = n = a->nlimbs * BYTES_PER_MPI_LIMB;
168 if (!n)
169 n++; /* avoid zero length allocation */
170 p = buffer = kmalloc(n, GFP_KERNEL);
171 if (!p)
172 return NULL;
173
174 for (i = a->nlimbs - 1; i >= 0; i--) {
175 alimb = a->d[i];
176#if BYTES_PER_MPI_LIMB == 4
177 *p++ = alimb >> 24;
178 *p++ = alimb >> 16;
179 *p++ = alimb >> 8;
180 *p++ = alimb;
181#elif BYTES_PER_MPI_LIMB == 8
182 *p++ = alimb >> 56;
183 *p++ = alimb >> 48;
184 *p++ = alimb >> 40;
185 *p++ = alimb >> 32;
186 *p++ = alimb >> 24;
187 *p++ = alimb >> 16;
188 *p++ = alimb >> 8;
189 *p++ = alimb;
190#else
191#error please implement for this limb size.
192#endif
193 }
194
195 /* this is sub-optimal but we need to do the shift operation
196 * because the caller has to free the returned buffer */
197 for (p = buffer; !*p && *nbytes; p++, --*nbytes)
198 ;
199 if (p != buffer)
200 memmove(buffer, p, *nbytes);
201
202 return buffer;
203}
204EXPORT_SYMBOL_GPL(mpi_get_buffer);
205
206/****************
207 * Use BUFFER to update MPI.
208 */
209int mpi_set_buffer(MPI a, const void *xbuffer, unsigned nbytes, int sign)
210{
211 const uint8_t *buffer = xbuffer, *p;
212 mpi_limb_t alimb;
213 int nlimbs;
214 int i;
215
216 nlimbs = (nbytes + BYTES_PER_MPI_LIMB - 1) / BYTES_PER_MPI_LIMB;
217 if (RESIZE_IF_NEEDED(a, nlimbs) < 0)
218 return -ENOMEM;
219 a->sign = sign;
220
221 for (i = 0, p = buffer + nbytes - 1; p >= buffer + BYTES_PER_MPI_LIMB;) {
222#if BYTES_PER_MPI_LIMB == 4
223 alimb = (mpi_limb_t) *p--;
224 alimb |= (mpi_limb_t) *p-- << 8;
225 alimb |= (mpi_limb_t) *p-- << 16;
226 alimb |= (mpi_limb_t) *p-- << 24;
227#elif BYTES_PER_MPI_LIMB == 8
228 alimb = (mpi_limb_t) *p--;
229 alimb |= (mpi_limb_t) *p-- << 8;
230 alimb |= (mpi_limb_t) *p-- << 16;
231 alimb |= (mpi_limb_t) *p-- << 24;
232 alimb |= (mpi_limb_t) *p-- << 32;
233 alimb |= (mpi_limb_t) *p-- << 40;
234 alimb |= (mpi_limb_t) *p-- << 48;
235 alimb |= (mpi_limb_t) *p-- << 56;
236#else
237#error please implement for this limb size.
238#endif
239 a->d[i++] = alimb;
240 }
241 if (p >= buffer) {
242#if BYTES_PER_MPI_LIMB == 4
243 alimb = *p--;
244 if (p >= buffer)
245 alimb |= (mpi_limb_t) *p-- << 8;
246 if (p >= buffer)
247 alimb |= (mpi_limb_t) *p-- << 16;
248 if (p >= buffer)
249 alimb |= (mpi_limb_t) *p-- << 24;
250#elif BYTES_PER_MPI_LIMB == 8
251 alimb = (mpi_limb_t) *p--;
252 if (p >= buffer)
253 alimb |= (mpi_limb_t) *p-- << 8;
254 if (p >= buffer)
255 alimb |= (mpi_limb_t) *p-- << 16;
256 if (p >= buffer)
257 alimb |= (mpi_limb_t) *p-- << 24;
258 if (p >= buffer)
259 alimb |= (mpi_limb_t) *p-- << 32;
260 if (p >= buffer)
261 alimb |= (mpi_limb_t) *p-- << 40;
262 if (p >= buffer)
263 alimb |= (mpi_limb_t) *p-- << 48;
264 if (p >= buffer)
265 alimb |= (mpi_limb_t) *p-- << 56;
266#else
267#error please implement for this limb size.
268#endif
269 a->d[i++] = alimb;
270 }
271 a->nlimbs = i;
272
273 if (i != nlimbs) {
274 pr_emerg("MPI: mpi_set_buffer: Assertion failed (%d != %d)", i,
275 nlimbs);
276 BUG();
277 }
278 return 0;
279}
280EXPORT_SYMBOL_GPL(mpi_set_buffer);
diff --git a/lib/mpi/mpih-cmp.c b/lib/mpi/mpih-cmp.c
new file mode 100644
index 000000000000..b2fd39677f1b
--- /dev/null
+++ b/lib/mpi/mpih-cmp.c
@@ -0,0 +1,56 @@
1/* mpihelp-sub.c - MPI helper functions
2 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
3 * Copyright (C) 1998, 1999, 2000, 2001 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31
32/****************
33 * Compare OP1_PTR/OP1_SIZE with OP2_PTR/OP2_SIZE.
34 * There are no restrictions on the relative sizes of
35 * the two arguments.
36 * Return 1 if OP1 > OP2, 0 if they are equal, and -1 if OP1 < OP2.
37 */
38int mpihelp_cmp(mpi_ptr_t op1_ptr, mpi_ptr_t op2_ptr, mpi_size_t size)
39{
40 mpi_size_t i;
41 mpi_limb_t op1_word, op2_word;
42
43 for (i = size - 1; i >= 0; i--) {
44 op1_word = op1_ptr[i];
45 op2_word = op2_ptr[i];
46 if (op1_word != op2_word)
47 goto diff;
48 }
49 return 0;
50
51diff:
52 /* This can *not* be simplified to
53 * op2_word - op2_word
54 * since that expression might give signed overflow. */
55 return (op1_word > op2_word) ? 1 : -1;
56}
diff --git a/lib/mpi/mpih-div.c b/lib/mpi/mpih-div.c
new file mode 100644
index 000000000000..cde1aaec18da
--- /dev/null
+++ b/lib/mpi/mpih-div.c
@@ -0,0 +1,545 @@
1/* mpihelp-div.c - MPI helper functions
2 * Copyright (C) 1994, 1996 Free Software Foundation, Inc.
3 * Copyright (C) 1998, 1999 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include "mpi-internal.h"
31#include "longlong.h"
32
33#ifndef UMUL_TIME
34#define UMUL_TIME 1
35#endif
36#ifndef UDIV_TIME
37#define UDIV_TIME UMUL_TIME
38#endif
39
40/* FIXME: We should be using invert_limb (or invert_normalized_limb)
41 * here (not udiv_qrnnd).
42 */
43
44mpi_limb_t
45mpihelp_mod_1(mpi_ptr_t dividend_ptr, mpi_size_t dividend_size,
46 mpi_limb_t divisor_limb)
47{
48 mpi_size_t i;
49 mpi_limb_t n1, n0, r;
50 int dummy;
51
52 /* Botch: Should this be handled at all? Rely on callers? */
53 if (!dividend_size)
54 return 0;
55
56 /* If multiplication is much faster than division, and the
57 * dividend is large, pre-invert the divisor, and use
58 * only multiplications in the inner loop.
59 *
60 * This test should be read:
61 * Does it ever help to use udiv_qrnnd_preinv?
62 * && Does what we save compensate for the inversion overhead?
63 */
64 if (UDIV_TIME > (2 * UMUL_TIME + 6)
65 && (UDIV_TIME - (2 * UMUL_TIME + 6)) * dividend_size > UDIV_TIME) {
66 int normalization_steps;
67
68 count_leading_zeros(normalization_steps, divisor_limb);
69 if (normalization_steps) {
70 mpi_limb_t divisor_limb_inverted;
71
72 divisor_limb <<= normalization_steps;
73
74 /* Compute (2**2N - 2**N * DIVISOR_LIMB) / DIVISOR_LIMB. The
75 * result is a (N+1)-bit approximation to 1/DIVISOR_LIMB, with the
76 * most significant bit (with weight 2**N) implicit.
77 *
78 * Special case for DIVISOR_LIMB == 100...000.
79 */
80 if (!(divisor_limb << 1))
81 divisor_limb_inverted = ~(mpi_limb_t) 0;
82 else
83 udiv_qrnnd(divisor_limb_inverted, dummy,
84 -divisor_limb, 0, divisor_limb);
85
86 n1 = dividend_ptr[dividend_size - 1];
87 r = n1 >> (BITS_PER_MPI_LIMB - normalization_steps);
88
89 /* Possible optimization:
90 * if (r == 0
91 * && divisor_limb > ((n1 << normalization_steps)
92 * | (dividend_ptr[dividend_size - 2] >> ...)))
93 * ...one division less...
94 */
95 for (i = dividend_size - 2; i >= 0; i--) {
96 n0 = dividend_ptr[i];
97 UDIV_QRNND_PREINV(dummy, r, r,
98 ((n1 << normalization_steps)
99 | (n0 >>
100 (BITS_PER_MPI_LIMB -
101 normalization_steps))),
102 divisor_limb,
103 divisor_limb_inverted);
104 n1 = n0;
105 }
106 UDIV_QRNND_PREINV(dummy, r, r,
107 n1 << normalization_steps,
108 divisor_limb, divisor_limb_inverted);
109 return r >> normalization_steps;
110 } else {
111 mpi_limb_t divisor_limb_inverted;
112
113 /* Compute (2**2N - 2**N * DIVISOR_LIMB) / DIVISOR_LIMB. The
114 * result is a (N+1)-bit approximation to 1/DIVISOR_LIMB, with the
115 * most significant bit (with weight 2**N) implicit.
116 *
117 * Special case for DIVISOR_LIMB == 100...000.
118 */
119 if (!(divisor_limb << 1))
120 divisor_limb_inverted = ~(mpi_limb_t) 0;
121 else
122 udiv_qrnnd(divisor_limb_inverted, dummy,
123 -divisor_limb, 0, divisor_limb);
124
125 i = dividend_size - 1;
126 r = dividend_ptr[i];
127
128 if (r >= divisor_limb)
129 r = 0;
130 else
131 i--;
132
133 for (; i >= 0; i--) {
134 n0 = dividend_ptr[i];
135 UDIV_QRNND_PREINV(dummy, r, r,
136 n0, divisor_limb,
137 divisor_limb_inverted);
138 }
139 return r;
140 }
141 } else {
142 if (UDIV_NEEDS_NORMALIZATION) {
143 int normalization_steps;
144
145 count_leading_zeros(normalization_steps, divisor_limb);
146 if (normalization_steps) {
147 divisor_limb <<= normalization_steps;
148
149 n1 = dividend_ptr[dividend_size - 1];
150 r = n1 >> (BITS_PER_MPI_LIMB -
151 normalization_steps);
152
153 /* Possible optimization:
154 * if (r == 0
155 * && divisor_limb > ((n1 << normalization_steps)
156 * | (dividend_ptr[dividend_size - 2] >> ...)))
157 * ...one division less...
158 */
159 for (i = dividend_size - 2; i >= 0; i--) {
160 n0 = dividend_ptr[i];
161 udiv_qrnnd(dummy, r, r,
162 ((n1 << normalization_steps)
163 | (n0 >>
164 (BITS_PER_MPI_LIMB -
165 normalization_steps))),
166 divisor_limb);
167 n1 = n0;
168 }
169 udiv_qrnnd(dummy, r, r,
170 n1 << normalization_steps,
171 divisor_limb);
172 return r >> normalization_steps;
173 }
174 }
175 /* No normalization needed, either because udiv_qrnnd doesn't require
176 * it, or because DIVISOR_LIMB is already normalized. */
177 i = dividend_size - 1;
178 r = dividend_ptr[i];
179
180 if (r >= divisor_limb)
181 r = 0;
182 else
183 i--;
184
185 for (; i >= 0; i--) {
186 n0 = dividend_ptr[i];
187 udiv_qrnnd(dummy, r, r, n0, divisor_limb);
188 }
189 return r;
190 }
191}
192
193/* Divide num (NP/NSIZE) by den (DP/DSIZE) and write
194 * the NSIZE-DSIZE least significant quotient limbs at QP
195 * and the DSIZE long remainder at NP. If QEXTRA_LIMBS is
196 * non-zero, generate that many fraction bits and append them after the
197 * other quotient limbs.
198 * Return the most significant limb of the quotient, this is always 0 or 1.
199 *
200 * Preconditions:
201 * 0. NSIZE >= DSIZE.
202 * 1. The most significant bit of the divisor must be set.
203 * 2. QP must either not overlap with the input operands at all, or
204 * QP + DSIZE >= NP must hold true. (This means that it's
205 * possible to put the quotient in the high part of NUM, right after the
206 * remainder in NUM.
207 * 3. NSIZE >= DSIZE, even if QEXTRA_LIMBS is non-zero.
208 */
209
210mpi_limb_t
211mpihelp_divrem(mpi_ptr_t qp, mpi_size_t qextra_limbs,
212 mpi_ptr_t np, mpi_size_t nsize, mpi_ptr_t dp, mpi_size_t dsize)
213{
214 mpi_limb_t most_significant_q_limb = 0;
215
216 switch (dsize) {
217 case 0:
218 /* We are asked to divide by zero, so go ahead and do it! (To make
219 the compiler not remove this statement, return the value.) */
220 /*
221 * existing clients of this function have been modified
222 * not to call it with dsize == 0, so this should not happen
223 */
224 return 1 / dsize;
225
226 case 1:
227 {
228 mpi_size_t i;
229 mpi_limb_t n1;
230 mpi_limb_t d;
231
232 d = dp[0];
233 n1 = np[nsize - 1];
234
235 if (n1 >= d) {
236 n1 -= d;
237 most_significant_q_limb = 1;
238 }
239
240 qp += qextra_limbs;
241 for (i = nsize - 2; i >= 0; i--)
242 udiv_qrnnd(qp[i], n1, n1, np[i], d);
243 qp -= qextra_limbs;
244
245 for (i = qextra_limbs - 1; i >= 0; i--)
246 udiv_qrnnd(qp[i], n1, n1, 0, d);
247
248 np[0] = n1;
249 }
250 break;
251
252 case 2:
253 {
254 mpi_size_t i;
255 mpi_limb_t n1, n0, n2;
256 mpi_limb_t d1, d0;
257
258 np += nsize - 2;
259 d1 = dp[1];
260 d0 = dp[0];
261 n1 = np[1];
262 n0 = np[0];
263
264 if (n1 >= d1 && (n1 > d1 || n0 >= d0)) {
265 sub_ddmmss(n1, n0, n1, n0, d1, d0);
266 most_significant_q_limb = 1;
267 }
268
269 for (i = qextra_limbs + nsize - 2 - 1; i >= 0; i--) {
270 mpi_limb_t q;
271 mpi_limb_t r;
272
273 if (i >= qextra_limbs)
274 np--;
275 else
276 np[0] = 0;
277
278 if (n1 == d1) {
279 /* Q should be either 111..111 or 111..110. Need special
280 * treatment of this rare case as normal division would
281 * give overflow. */
282 q = ~(mpi_limb_t) 0;
283
284 r = n0 + d1;
285 if (r < d1) { /* Carry in the addition? */
286 add_ssaaaa(n1, n0, r - d0,
287 np[0], 0, d0);
288 qp[i] = q;
289 continue;
290 }
291 n1 = d0 - (d0 != 0 ? 1 : 0);
292 n0 = -d0;
293 } else {
294 udiv_qrnnd(q, r, n1, n0, d1);
295 umul_ppmm(n1, n0, d0, q);
296 }
297
298 n2 = np[0];
299q_test:
300 if (n1 > r || (n1 == r && n0 > n2)) {
301 /* The estimated Q was too large. */
302 q--;
303 sub_ddmmss(n1, n0, n1, n0, 0, d0);
304 r += d1;
305 if (r >= d1) /* If not carry, test Q again. */
306 goto q_test;
307 }
308
309 qp[i] = q;
310 sub_ddmmss(n1, n0, r, n2, n1, n0);
311 }
312 np[1] = n1;
313 np[0] = n0;
314 }
315 break;
316
317 default:
318 {
319 mpi_size_t i;
320 mpi_limb_t dX, d1, n0;
321
322 np += nsize - dsize;
323 dX = dp[dsize - 1];
324 d1 = dp[dsize - 2];
325 n0 = np[dsize - 1];
326
327 if (n0 >= dX) {
328 if (n0 > dX
329 || mpihelp_cmp(np, dp, dsize - 1) >= 0) {
330 mpihelp_sub_n(np, np, dp, dsize);
331 n0 = np[dsize - 1];
332 most_significant_q_limb = 1;
333 }
334 }
335
336 for (i = qextra_limbs + nsize - dsize - 1; i >= 0; i--) {
337 mpi_limb_t q;
338 mpi_limb_t n1, n2;
339 mpi_limb_t cy_limb;
340
341 if (i >= qextra_limbs) {
342 np--;
343 n2 = np[dsize];
344 } else {
345 n2 = np[dsize - 1];
346 MPN_COPY_DECR(np + 1, np, dsize - 1);
347 np[0] = 0;
348 }
349
350 if (n0 == dX) {
351 /* This might over-estimate q, but it's probably not worth
352 * the extra code here to find out. */
353 q = ~(mpi_limb_t) 0;
354 } else {
355 mpi_limb_t r;
356
357 udiv_qrnnd(q, r, n0, np[dsize - 1], dX);
358 umul_ppmm(n1, n0, d1, q);
359
360 while (n1 > r
361 || (n1 == r
362 && n0 > np[dsize - 2])) {
363 q--;
364 r += dX;
365 if (r < dX) /* I.e. "carry in previous addition?" */
366 break;
367 n1 -= n0 < d1;
368 n0 -= d1;
369 }
370 }
371
372 /* Possible optimization: We already have (q * n0) and (1 * n1)
373 * after the calculation of q. Taking advantage of that, we
374 * could make this loop make two iterations less. */
375 cy_limb = mpihelp_submul_1(np, dp, dsize, q);
376
377 if (n2 != cy_limb) {
378 mpihelp_add_n(np, np, dp, dsize);
379 q--;
380 }
381
382 qp[i] = q;
383 n0 = np[dsize - 1];
384 }
385 }
386 }
387
388 return most_significant_q_limb;
389}
390
391/****************
392 * Divide (DIVIDEND_PTR,,DIVIDEND_SIZE) by DIVISOR_LIMB.
393 * Write DIVIDEND_SIZE limbs of quotient at QUOT_PTR.
394 * Return the single-limb remainder.
395 * There are no constraints on the value of the divisor.
396 *
397 * QUOT_PTR and DIVIDEND_PTR might point to the same limb.
398 */
399
400mpi_limb_t
401mpihelp_divmod_1(mpi_ptr_t quot_ptr,
402 mpi_ptr_t dividend_ptr, mpi_size_t dividend_size,
403 mpi_limb_t divisor_limb)
404{
405 mpi_size_t i;
406 mpi_limb_t n1, n0, r;
407 int dummy;
408
409 if (!dividend_size)
410 return 0;
411
412 /* If multiplication is much faster than division, and the
413 * dividend is large, pre-invert the divisor, and use
414 * only multiplications in the inner loop.
415 *
416 * This test should be read:
417 * Does it ever help to use udiv_qrnnd_preinv?
418 * && Does what we save compensate for the inversion overhead?
419 */
420 if (UDIV_TIME > (2 * UMUL_TIME + 6)
421 && (UDIV_TIME - (2 * UMUL_TIME + 6)) * dividend_size > UDIV_TIME) {
422 int normalization_steps;
423
424 count_leading_zeros(normalization_steps, divisor_limb);
425 if (normalization_steps) {
426 mpi_limb_t divisor_limb_inverted;
427
428 divisor_limb <<= normalization_steps;
429
430 /* Compute (2**2N - 2**N * DIVISOR_LIMB) / DIVISOR_LIMB. The
431 * result is a (N+1)-bit approximation to 1/DIVISOR_LIMB, with the
432 * most significant bit (with weight 2**N) implicit.
433 */
434 /* Special case for DIVISOR_LIMB == 100...000. */
435 if (!(divisor_limb << 1))
436 divisor_limb_inverted = ~(mpi_limb_t) 0;
437 else
438 udiv_qrnnd(divisor_limb_inverted, dummy,
439 -divisor_limb, 0, divisor_limb);
440
441 n1 = dividend_ptr[dividend_size - 1];
442 r = n1 >> (BITS_PER_MPI_LIMB - normalization_steps);
443
444 /* Possible optimization:
445 * if (r == 0
446 * && divisor_limb > ((n1 << normalization_steps)
447 * | (dividend_ptr[dividend_size - 2] >> ...)))
448 * ...one division less...
449 */
450 for (i = dividend_size - 2; i >= 0; i--) {
451 n0 = dividend_ptr[i];
452 UDIV_QRNND_PREINV(quot_ptr[i + 1], r, r,
453 ((n1 << normalization_steps)
454 | (n0 >>
455 (BITS_PER_MPI_LIMB -
456 normalization_steps))),
457 divisor_limb,
458 divisor_limb_inverted);
459 n1 = n0;
460 }
461 UDIV_QRNND_PREINV(quot_ptr[0], r, r,
462 n1 << normalization_steps,
463 divisor_limb, divisor_limb_inverted);
464 return r >> normalization_steps;
465 } else {
466 mpi_limb_t divisor_limb_inverted;
467
468 /* Compute (2**2N - 2**N * DIVISOR_LIMB) / DIVISOR_LIMB. The
469 * result is a (N+1)-bit approximation to 1/DIVISOR_LIMB, with the
470 * most significant bit (with weight 2**N) implicit.
471 */
472 /* Special case for DIVISOR_LIMB == 100...000. */
473 if (!(divisor_limb << 1))
474 divisor_limb_inverted = ~(mpi_limb_t) 0;
475 else
476 udiv_qrnnd(divisor_limb_inverted, dummy,
477 -divisor_limb, 0, divisor_limb);
478
479 i = dividend_size - 1;
480 r = dividend_ptr[i];
481
482 if (r >= divisor_limb)
483 r = 0;
484 else
485 quot_ptr[i--] = 0;
486
487 for (; i >= 0; i--) {
488 n0 = dividend_ptr[i];
489 UDIV_QRNND_PREINV(quot_ptr[i], r, r,
490 n0, divisor_limb,
491 divisor_limb_inverted);
492 }
493 return r;
494 }
495 } else {
496 if (UDIV_NEEDS_NORMALIZATION) {
497 int normalization_steps;
498
499 count_leading_zeros(normalization_steps, divisor_limb);
500 if (normalization_steps) {
501 divisor_limb <<= normalization_steps;
502
503 n1 = dividend_ptr[dividend_size - 1];
504 r = n1 >> (BITS_PER_MPI_LIMB -
505 normalization_steps);
506
507 /* Possible optimization:
508 * if (r == 0
509 * && divisor_limb > ((n1 << normalization_steps)
510 * | (dividend_ptr[dividend_size - 2] >> ...)))
511 * ...one division less...
512 */
513 for (i = dividend_size - 2; i >= 0; i--) {
514 n0 = dividend_ptr[i];
515 udiv_qrnnd(quot_ptr[i + 1], r, r,
516 ((n1 << normalization_steps)
517 | (n0 >>
518 (BITS_PER_MPI_LIMB -
519 normalization_steps))),
520 divisor_limb);
521 n1 = n0;
522 }
523 udiv_qrnnd(quot_ptr[0], r, r,
524 n1 << normalization_steps,
525 divisor_limb);
526 return r >> normalization_steps;
527 }
528 }
529 /* No normalization needed, either because udiv_qrnnd doesn't require
530 * it, or because DIVISOR_LIMB is already normalized. */
531 i = dividend_size - 1;
532 r = dividend_ptr[i];
533
534 if (r >= divisor_limb)
535 r = 0;
536 else
537 quot_ptr[i--] = 0;
538
539 for (; i >= 0; i--) {
540 n0 = dividend_ptr[i];
541 udiv_qrnnd(quot_ptr[i], r, r, n0, divisor_limb);
542 }
543 return r;
544 }
545}
diff --git a/lib/mpi/mpih-mul.c b/lib/mpi/mpih-mul.c
new file mode 100644
index 000000000000..c69c5eef233b
--- /dev/null
+++ b/lib/mpi/mpih-mul.c
@@ -0,0 +1,527 @@
1/* mpihelp-mul.c - MPI helper functions
2 * Copyright (C) 1994, 1996, 1998, 1999,
3 * 2000 Free Software Foundation, Inc.
4 *
5 * This file is part of GnuPG.
6 *
7 * GnuPG is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * GnuPG is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, write to the Free Software
19 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
20 *
21 * Note: This code is heavily based on the GNU MP Library.
22 * Actually it's the same code with only minor changes in the
23 * way the data is stored; this is to support the abstraction
24 * of an optional secure memory allocation which may be used
25 * to avoid revealing of sensitive data due to paging etc.
26 * The GNU MP Library itself is published under the LGPL;
27 * however I decided to publish this code under the plain GPL.
28 */
29
30#include <linux/string.h>
31#include "mpi-internal.h"
32#include "longlong.h"
33
34#define MPN_MUL_N_RECURSE(prodp, up, vp, size, tspace) \
35 do { \
36 if ((size) < KARATSUBA_THRESHOLD) \
37 mul_n_basecase(prodp, up, vp, size); \
38 else \
39 mul_n(prodp, up, vp, size, tspace); \
40 } while (0);
41
42#define MPN_SQR_N_RECURSE(prodp, up, size, tspace) \
43 do { \
44 if ((size) < KARATSUBA_THRESHOLD) \
45 mpih_sqr_n_basecase(prodp, up, size); \
46 else \
47 mpih_sqr_n(prodp, up, size, tspace); \
48 } while (0);
49
50/* Multiply the natural numbers u (pointed to by UP) and v (pointed to by VP),
51 * both with SIZE limbs, and store the result at PRODP. 2 * SIZE limbs are
52 * always stored. Return the most significant limb.
53 *
54 * Argument constraints:
55 * 1. PRODP != UP and PRODP != VP, i.e. the destination
56 * must be distinct from the multiplier and the multiplicand.
57 *
58 *
59 * Handle simple cases with traditional multiplication.
60 *
61 * This is the most critical code of multiplication. All multiplies rely
62 * on this, both small and huge. Small ones arrive here immediately. Huge
63 * ones arrive here as this is the base case for Karatsuba's recursive
64 * algorithm below.
65 */
66
67static mpi_limb_t
68mul_n_basecase(mpi_ptr_t prodp, mpi_ptr_t up, mpi_ptr_t vp, mpi_size_t size)
69{
70 mpi_size_t i;
71 mpi_limb_t cy;
72 mpi_limb_t v_limb;
73
74 /* Multiply by the first limb in V separately, as the result can be
75 * stored (not added) to PROD. We also avoid a loop for zeroing. */
76 v_limb = vp[0];
77 if (v_limb <= 1) {
78 if (v_limb == 1)
79 MPN_COPY(prodp, up, size);
80 else
81 MPN_ZERO(prodp, size);
82 cy = 0;
83 } else
84 cy = mpihelp_mul_1(prodp, up, size, v_limb);
85
86 prodp[size] = cy;
87 prodp++;
88
89 /* For each iteration in the outer loop, multiply one limb from
90 * U with one limb from V, and add it to PROD. */
91 for (i = 1; i < size; i++) {
92 v_limb = vp[i];
93 if (v_limb <= 1) {
94 cy = 0;
95 if (v_limb == 1)
96 cy = mpihelp_add_n(prodp, prodp, up, size);
97 } else
98 cy = mpihelp_addmul_1(prodp, up, size, v_limb);
99
100 prodp[size] = cy;
101 prodp++;
102 }
103
104 return cy;
105}
106
107static void
108mul_n(mpi_ptr_t prodp, mpi_ptr_t up, mpi_ptr_t vp,
109 mpi_size_t size, mpi_ptr_t tspace)
110{
111 if (size & 1) {
112 /* The size is odd, and the code below doesn't handle that.
113 * Multiply the least significant (size - 1) limbs with a recursive
114 * call, and handle the most significant limb of S1 and S2
115 * separately.
116 * A slightly faster way to do this would be to make the Karatsuba
117 * code below behave as if the size were even, and let it check for
118 * odd size in the end. I.e., in essence move this code to the end.
119 * Doing so would save us a recursive call, and potentially make the
120 * stack grow a lot less.
121 */
122 mpi_size_t esize = size - 1; /* even size */
123 mpi_limb_t cy_limb;
124
125 MPN_MUL_N_RECURSE(prodp, up, vp, esize, tspace);
126 cy_limb = mpihelp_addmul_1(prodp + esize, up, esize, vp[esize]);
127 prodp[esize + esize] = cy_limb;
128 cy_limb = mpihelp_addmul_1(prodp + esize, vp, size, up[esize]);
129 prodp[esize + size] = cy_limb;
130 } else {
131 /* Anatolij Alekseevich Karatsuba's divide-and-conquer algorithm.
132 *
133 * Split U in two pieces, U1 and U0, such that
134 * U = U0 + U1*(B**n),
135 * and V in V1 and V0, such that
136 * V = V0 + V1*(B**n).
137 *
138 * UV is then computed recursively using the identity
139 *
140 * 2n n n n
141 * UV = (B + B )U V + B (U -U )(V -V ) + (B + 1)U V
142 * 1 1 1 0 0 1 0 0
143 *
144 * Where B = 2**BITS_PER_MP_LIMB.
145 */
146 mpi_size_t hsize = size >> 1;
147 mpi_limb_t cy;
148 int negflg;
149
150 /* Product H. ________________ ________________
151 * |_____U1 x V1____||____U0 x V0_____|
152 * Put result in upper part of PROD and pass low part of TSPACE
153 * as new TSPACE.
154 */
155 MPN_MUL_N_RECURSE(prodp + size, up + hsize, vp + hsize, hsize,
156 tspace);
157
158 /* Product M. ________________
159 * |_(U1-U0)(V0-V1)_|
160 */
161 if (mpihelp_cmp(up + hsize, up, hsize) >= 0) {
162 mpihelp_sub_n(prodp, up + hsize, up, hsize);
163 negflg = 0;
164 } else {
165 mpihelp_sub_n(prodp, up, up + hsize, hsize);
166 negflg = 1;
167 }
168 if (mpihelp_cmp(vp + hsize, vp, hsize) >= 0) {
169 mpihelp_sub_n(prodp + hsize, vp + hsize, vp, hsize);
170 negflg ^= 1;
171 } else {
172 mpihelp_sub_n(prodp + hsize, vp, vp + hsize, hsize);
173 /* No change of NEGFLG. */
174 }
175 /* Read temporary operands from low part of PROD.
176 * Put result in low part of TSPACE using upper part of TSPACE
177 * as new TSPACE.
178 */
179 MPN_MUL_N_RECURSE(tspace, prodp, prodp + hsize, hsize,
180 tspace + size);
181
182 /* Add/copy product H. */
183 MPN_COPY(prodp + hsize, prodp + size, hsize);
184 cy = mpihelp_add_n(prodp + size, prodp + size,
185 prodp + size + hsize, hsize);
186
187 /* Add product M (if NEGFLG M is a negative number) */
188 if (negflg)
189 cy -=
190 mpihelp_sub_n(prodp + hsize, prodp + hsize, tspace,
191 size);
192 else
193 cy +=
194 mpihelp_add_n(prodp + hsize, prodp + hsize, tspace,
195 size);
196
197 /* Product L. ________________ ________________
198 * |________________||____U0 x V0_____|
199 * Read temporary operands from low part of PROD.
200 * Put result in low part of TSPACE using upper part of TSPACE
201 * as new TSPACE.
202 */
203 MPN_MUL_N_RECURSE(tspace, up, vp, hsize, tspace + size);
204
205 /* Add/copy Product L (twice) */
206
207 cy += mpihelp_add_n(prodp + hsize, prodp + hsize, tspace, size);
208 if (cy)
209 mpihelp_add_1(prodp + hsize + size,
210 prodp + hsize + size, hsize, cy);
211
212 MPN_COPY(prodp, tspace, hsize);
213 cy = mpihelp_add_n(prodp + hsize, prodp + hsize, tspace + hsize,
214 hsize);
215 if (cy)
216 mpihelp_add_1(prodp + size, prodp + size, size, 1);
217 }
218}
219
220void mpih_sqr_n_basecase(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t size)
221{
222 mpi_size_t i;
223 mpi_limb_t cy_limb;
224 mpi_limb_t v_limb;
225
226 /* Multiply by the first limb in V separately, as the result can be
227 * stored (not added) to PROD. We also avoid a loop for zeroing. */
228 v_limb = up[0];
229 if (v_limb <= 1) {
230 if (v_limb == 1)
231 MPN_COPY(prodp, up, size);
232 else
233 MPN_ZERO(prodp, size);
234 cy_limb = 0;
235 } else
236 cy_limb = mpihelp_mul_1(prodp, up, size, v_limb);
237
238 prodp[size] = cy_limb;
239 prodp++;
240
241 /* For each iteration in the outer loop, multiply one limb from
242 * U with one limb from V, and add it to PROD. */
243 for (i = 1; i < size; i++) {
244 v_limb = up[i];
245 if (v_limb <= 1) {
246 cy_limb = 0;
247 if (v_limb == 1)
248 cy_limb = mpihelp_add_n(prodp, prodp, up, size);
249 } else
250 cy_limb = mpihelp_addmul_1(prodp, up, size, v_limb);
251
252 prodp[size] = cy_limb;
253 prodp++;
254 }
255}
256
257void
258mpih_sqr_n(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t size, mpi_ptr_t tspace)
259{
260 if (size & 1) {
261 /* The size is odd, and the code below doesn't handle that.
262 * Multiply the least significant (size - 1) limbs with a recursive
263 * call, and handle the most significant limb of S1 and S2
264 * separately.
265 * A slightly faster way to do this would be to make the Karatsuba
266 * code below behave as if the size were even, and let it check for
267 * odd size in the end. I.e., in essence move this code to the end.
268 * Doing so would save us a recursive call, and potentially make the
269 * stack grow a lot less.
270 */
271 mpi_size_t esize = size - 1; /* even size */
272 mpi_limb_t cy_limb;
273
274 MPN_SQR_N_RECURSE(prodp, up, esize, tspace);
275 cy_limb = mpihelp_addmul_1(prodp + esize, up, esize, up[esize]);
276 prodp[esize + esize] = cy_limb;
277 cy_limb = mpihelp_addmul_1(prodp + esize, up, size, up[esize]);
278
279 prodp[esize + size] = cy_limb;
280 } else {
281 mpi_size_t hsize = size >> 1;
282 mpi_limb_t cy;
283
284 /* Product H. ________________ ________________
285 * |_____U1 x U1____||____U0 x U0_____|
286 * Put result in upper part of PROD and pass low part of TSPACE
287 * as new TSPACE.
288 */
289 MPN_SQR_N_RECURSE(prodp + size, up + hsize, hsize, tspace);
290
291 /* Product M. ________________
292 * |_(U1-U0)(U0-U1)_|
293 */
294 if (mpihelp_cmp(up + hsize, up, hsize) >= 0)
295 mpihelp_sub_n(prodp, up + hsize, up, hsize);
296 else
297 mpihelp_sub_n(prodp, up, up + hsize, hsize);
298
299 /* Read temporary operands from low part of PROD.
300 * Put result in low part of TSPACE using upper part of TSPACE
301 * as new TSPACE. */
302 MPN_SQR_N_RECURSE(tspace, prodp, hsize, tspace + size);
303
304 /* Add/copy product H */
305 MPN_COPY(prodp + hsize, prodp + size, hsize);
306 cy = mpihelp_add_n(prodp + size, prodp + size,
307 prodp + size + hsize, hsize);
308
309 /* Add product M (if NEGFLG M is a negative number). */
310 cy -= mpihelp_sub_n(prodp + hsize, prodp + hsize, tspace, size);
311
312 /* Product L. ________________ ________________
313 * |________________||____U0 x U0_____|
314 * Read temporary operands from low part of PROD.
315 * Put result in low part of TSPACE using upper part of TSPACE
316 * as new TSPACE. */
317 MPN_SQR_N_RECURSE(tspace, up, hsize, tspace + size);
318
319 /* Add/copy Product L (twice). */
320 cy += mpihelp_add_n(prodp + hsize, prodp + hsize, tspace, size);
321 if (cy)
322 mpihelp_add_1(prodp + hsize + size,
323 prodp + hsize + size, hsize, cy);
324
325 MPN_COPY(prodp, tspace, hsize);
326 cy = mpihelp_add_n(prodp + hsize, prodp + hsize, tspace + hsize,
327 hsize);
328 if (cy)
329 mpihelp_add_1(prodp + size, prodp + size, size, 1);
330 }
331}
332
333/* This should be made into an inline function in gmp.h. */
334int mpihelp_mul_n(mpi_ptr_t prodp, mpi_ptr_t up, mpi_ptr_t vp, mpi_size_t size)
335{
336 if (up == vp) {
337 if (size < KARATSUBA_THRESHOLD)
338 mpih_sqr_n_basecase(prodp, up, size);
339 else {
340 mpi_ptr_t tspace;
341 tspace = mpi_alloc_limb_space(2 * size);
342 if (!tspace)
343 return -ENOMEM;
344 mpih_sqr_n(prodp, up, size, tspace);
345 mpi_free_limb_space(tspace);
346 }
347 } else {
348 if (size < KARATSUBA_THRESHOLD)
349 mul_n_basecase(prodp, up, vp, size);
350 else {
351 mpi_ptr_t tspace;
352 tspace = mpi_alloc_limb_space(2 * size);
353 if (!tspace)
354 return -ENOMEM;
355 mul_n(prodp, up, vp, size, tspace);
356 mpi_free_limb_space(tspace);
357 }
358 }
359
360 return 0;
361}
362
363int
364mpihelp_mul_karatsuba_case(mpi_ptr_t prodp,
365 mpi_ptr_t up, mpi_size_t usize,
366 mpi_ptr_t vp, mpi_size_t vsize,
367 struct karatsuba_ctx *ctx)
368{
369 mpi_limb_t cy;
370
371 if (!ctx->tspace || ctx->tspace_size < vsize) {
372 if (ctx->tspace)
373 mpi_free_limb_space(ctx->tspace);
374 ctx->tspace = mpi_alloc_limb_space(2 * vsize);
375 if (!ctx->tspace)
376 return -ENOMEM;
377 ctx->tspace_size = vsize;
378 }
379
380 MPN_MUL_N_RECURSE(prodp, up, vp, vsize, ctx->tspace);
381
382 prodp += vsize;
383 up += vsize;
384 usize -= vsize;
385 if (usize >= vsize) {
386 if (!ctx->tp || ctx->tp_size < vsize) {
387 if (ctx->tp)
388 mpi_free_limb_space(ctx->tp);
389 ctx->tp = mpi_alloc_limb_space(2 * vsize);
390 if (!ctx->tp) {
391 if (ctx->tspace)
392 mpi_free_limb_space(ctx->tspace);
393 ctx->tspace = NULL;
394 return -ENOMEM;
395 }
396 ctx->tp_size = vsize;
397 }
398
399 do {
400 MPN_MUL_N_RECURSE(ctx->tp, up, vp, vsize, ctx->tspace);
401 cy = mpihelp_add_n(prodp, prodp, ctx->tp, vsize);
402 mpihelp_add_1(prodp + vsize, ctx->tp + vsize, vsize,
403 cy);
404 prodp += vsize;
405 up += vsize;
406 usize -= vsize;
407 } while (usize >= vsize);
408 }
409
410 if (usize) {
411 if (usize < KARATSUBA_THRESHOLD) {
412 mpi_limb_t tmp;
413 if (mpihelp_mul(ctx->tspace, vp, vsize, up, usize, &tmp)
414 < 0)
415 return -ENOMEM;
416 } else {
417 if (!ctx->next) {
418 ctx->next = kzalloc(sizeof *ctx, GFP_KERNEL);
419 if (!ctx->next)
420 return -ENOMEM;
421 }
422 if (mpihelp_mul_karatsuba_case(ctx->tspace,
423 vp, vsize,
424 up, usize,
425 ctx->next) < 0)
426 return -ENOMEM;
427 }
428
429 cy = mpihelp_add_n(prodp, prodp, ctx->tspace, vsize);
430 mpihelp_add_1(prodp + vsize, ctx->tspace + vsize, usize, cy);
431 }
432
433 return 0;
434}
435
436void mpihelp_release_karatsuba_ctx(struct karatsuba_ctx *ctx)
437{
438 struct karatsuba_ctx *ctx2;
439
440 if (ctx->tp)
441 mpi_free_limb_space(ctx->tp);
442 if (ctx->tspace)
443 mpi_free_limb_space(ctx->tspace);
444 for (ctx = ctx->next; ctx; ctx = ctx2) {
445 ctx2 = ctx->next;
446 if (ctx->tp)
447 mpi_free_limb_space(ctx->tp);
448 if (ctx->tspace)
449 mpi_free_limb_space(ctx->tspace);
450 kfree(ctx);
451 }
452}
453
454/* Multiply the natural numbers u (pointed to by UP, with USIZE limbs)
455 * and v (pointed to by VP, with VSIZE limbs), and store the result at
456 * PRODP. USIZE + VSIZE limbs are always stored, but if the input
457 * operands are normalized. Return the most significant limb of the
458 * result.
459 *
460 * NOTE: The space pointed to by PRODP is overwritten before finished
461 * with U and V, so overlap is an error.
462 *
463 * Argument constraints:
464 * 1. USIZE >= VSIZE.
465 * 2. PRODP != UP and PRODP != VP, i.e. the destination
466 * must be distinct from the multiplier and the multiplicand.
467 */
468
469int
470mpihelp_mul(mpi_ptr_t prodp, mpi_ptr_t up, mpi_size_t usize,
471 mpi_ptr_t vp, mpi_size_t vsize, mpi_limb_t *_result)
472{
473 mpi_ptr_t prod_endp = prodp + usize + vsize - 1;
474 mpi_limb_t cy;
475 struct karatsuba_ctx ctx;
476
477 if (vsize < KARATSUBA_THRESHOLD) {
478 mpi_size_t i;
479 mpi_limb_t v_limb;
480
481 if (!vsize) {
482 *_result = 0;
483 return 0;
484 }
485
486 /* Multiply by the first limb in V separately, as the result can be
487 * stored (not added) to PROD. We also avoid a loop for zeroing. */
488 v_limb = vp[0];
489 if (v_limb <= 1) {
490 if (v_limb == 1)
491 MPN_COPY(prodp, up, usize);
492 else
493 MPN_ZERO(prodp, usize);
494 cy = 0;
495 } else
496 cy = mpihelp_mul_1(prodp, up, usize, v_limb);
497
498 prodp[usize] = cy;
499 prodp++;
500
501 /* For each iteration in the outer loop, multiply one limb from
502 * U with one limb from V, and add it to PROD. */
503 for (i = 1; i < vsize; i++) {
504 v_limb = vp[i];
505 if (v_limb <= 1) {
506 cy = 0;
507 if (v_limb == 1)
508 cy = mpihelp_add_n(prodp, prodp, up,
509 usize);
510 } else
511 cy = mpihelp_addmul_1(prodp, up, usize, v_limb);
512
513 prodp[usize] = cy;
514 prodp++;
515 }
516
517 *_result = cy;
518 return 0;
519 }
520
521 memset(&ctx, 0, sizeof ctx);
522 if (mpihelp_mul_karatsuba_case(prodp, up, usize, vp, vsize, &ctx) < 0)
523 return -ENOMEM;
524 mpihelp_release_karatsuba_ctx(&ctx);
525 *_result = *prod_endp;
526 return 0;
527}
diff --git a/lib/mpi/mpiutil.c b/lib/mpi/mpiutil.c
new file mode 100644
index 000000000000..26e4ed31e256
--- /dev/null
+++ b/lib/mpi/mpiutil.c
@@ -0,0 +1,211 @@
1/* mpiutil.ac - Utility functions for MPI
2 * Copyright (C) 1998, 1999 Free Software Foundation, Inc.
3 *
4 * This file is part of GnuPG.
5 *
6 * GnuPG is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * GnuPG is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA
19 */
20
21#include "mpi-internal.h"
22
23/****************
24 * Note: It was a bad idea to use the number of limbs to allocate
25 * because on a alpha the limbs are large but we normally need
26 * integers of n bits - So we should chnage this to bits (or bytes).
27 *
28 * But mpi_alloc is used in a lot of places :-)
29 */
30MPI mpi_alloc(unsigned nlimbs)
31{
32 MPI a;
33
34 a = kmalloc(sizeof *a, GFP_KERNEL);
35 if (!a)
36 return a;
37
38 if (nlimbs) {
39 a->d = mpi_alloc_limb_space(nlimbs);
40 if (!a->d) {
41 kfree(a);
42 return NULL;
43 }
44 } else {
45 a->d = NULL;
46 }
47
48 a->alloced = nlimbs;
49 a->nlimbs = 0;
50 a->sign = 0;
51 a->flags = 0;
52 a->nbits = 0;
53 return a;
54}
55EXPORT_SYMBOL_GPL(mpi_alloc);
56
57mpi_ptr_t mpi_alloc_limb_space(unsigned nlimbs)
58{
59 size_t len = nlimbs * sizeof(mpi_limb_t);
60
61 if (!len)
62 return NULL;
63
64 return kmalloc(len, GFP_KERNEL);
65}
66
67void mpi_free_limb_space(mpi_ptr_t a)
68{
69 if (!a)
70 return;
71
72 kfree(a);
73}
74
75void mpi_assign_limb_space(MPI a, mpi_ptr_t ap, unsigned nlimbs)
76{
77 mpi_free_limb_space(a->d);
78 a->d = ap;
79 a->alloced = nlimbs;
80}
81
82/****************
83 * Resize the array of A to NLIMBS. the additional space is cleared
84 * (set to 0) [done by m_realloc()]
85 */
86int mpi_resize(MPI a, unsigned nlimbs)
87{
88 void *p;
89
90 if (nlimbs <= a->alloced)
91 return 0; /* no need to do it */
92
93 if (a->d) {
94 p = kmalloc(nlimbs * sizeof(mpi_limb_t), GFP_KERNEL);
95 if (!p)
96 return -ENOMEM;
97 memcpy(p, a->d, a->alloced * sizeof(mpi_limb_t));
98 kfree(a->d);
99 a->d = p;
100 } else {
101 a->d = kzalloc(nlimbs * sizeof(mpi_limb_t), GFP_KERNEL);
102 if (!a->d)
103 return -ENOMEM;
104 }
105 a->alloced = nlimbs;
106 return 0;
107}
108
109void mpi_clear(MPI a)
110{
111 a->nlimbs = 0;
112 a->nbits = 0;
113 a->flags = 0;
114}
115
116void mpi_free(MPI a)
117{
118 if (!a)
119 return;
120
121 if (a->flags & 4)
122 kfree(a->d);
123 else
124 mpi_free_limb_space(a->d);
125
126 if (a->flags & ~7)
127 pr_info("invalid flag value in mpi\n");
128 kfree(a);
129}
130EXPORT_SYMBOL_GPL(mpi_free);
131
132/****************
133 * Note: This copy function should not interpret the MPI
134 * but copy it transparently.
135 */
136int mpi_copy(MPI *copied, const MPI a)
137{
138 size_t i;
139 MPI b;
140
141 *copied = NULL;
142
143 if (a) {
144 b = mpi_alloc(a->nlimbs);
145 if (!b)
146 return -ENOMEM;
147
148 b->nlimbs = a->nlimbs;
149 b->sign = a->sign;
150 b->flags = a->flags;
151 b->nbits = a->nbits;
152
153 for (i = 0; i < b->nlimbs; i++)
154 b->d[i] = a->d[i];
155
156 *copied = b;
157 }
158
159 return 0;
160}
161
162int mpi_set(MPI w, const MPI u)
163{
164 mpi_ptr_t wp, up;
165 mpi_size_t usize = u->nlimbs;
166 int usign = u->sign;
167
168 if (RESIZE_IF_NEEDED(w, (size_t) usize) < 0)
169 return -ENOMEM;
170
171 wp = w->d;
172 up = u->d;
173 MPN_COPY(wp, up, usize);
174 w->nlimbs = usize;
175 w->nbits = u->nbits;
176 w->flags = u->flags;
177 w->sign = usign;
178 return 0;
179}
180
181int mpi_set_ui(MPI w, unsigned long u)
182{
183 if (RESIZE_IF_NEEDED(w, 1) < 0)
184 return -ENOMEM;
185 w->d[0] = u;
186 w->nlimbs = u ? 1 : 0;
187 w->sign = 0;
188 w->nbits = 0;
189 w->flags = 0;
190 return 0;
191}
192
193MPI mpi_alloc_set_ui(unsigned long u)
194{
195 MPI w = mpi_alloc(1);
196 if (!w)
197 return w;
198 w->d[0] = u;
199 w->nlimbs = u ? 1 : 0;
200 w->sign = 0;
201 return w;
202}
203
204void mpi_swap(MPI a, MPI b)
205{
206 struct gcry_mpi tmp;
207
208 tmp = *a;
209 *a = *b;
210 *b = tmp;
211}
generated by cgit v1.2.2 (git 2.25.0) at 2025-11-11 16:28:23 -0500