1 // SPDX-License-Identifier: GPL-2.0-or-later
2 /* mpihelp-mul.c - MPI helper functions
3 * Copyright (C) 1994, 1996, 1998, 1999,
4 * 2000 Free Software Foundation, Inc.
6 * This file is part of GnuPG.
8 * Note: This code is heavily based on the GNU MP Library.
9 * Actually it's the same code with only minor changes in the
10 * way the data is stored; this is to support the abstraction
11 * of an optional secure memory allocation which may be used
12 * to avoid revealing of sensitive data due to paging etc.
13 * The GNU MP Library itself is published under the LGPL;
14 * however I decided to publish this code under the plain GPL.
17 #include <linux/string.h>
18 #include "mpi-internal.h"
21 #define MPN_MUL_N_RECURSE(prodp, up, vp, size, tspace) \
23 if ((size) < KARATSUBA_THRESHOLD) \
24 mul_n_basecase(prodp, up, vp, size); \
26 mul_n(prodp, up, vp, size, tspace); \
29 #define MPN_SQR_N_RECURSE(prodp, up, size, tspace) \
31 if ((size) < KARATSUBA_THRESHOLD) \
32 mpih_sqr_n_basecase(prodp, up, size); \
34 mpih_sqr_n(prodp, up, size, tspace); \
37 /* Multiply the natural numbers u (pointed to by UP) and v (pointed to by VP),
38 * both with SIZE limbs, and store the result at PRODP. 2 * SIZE limbs are
39 * always stored. Return the most significant limb.
41 * Argument constraints:
42 * 1. PRODP != UP and PRODP != VP, i.e. the destination
43 * must be distinct from the multiplier and the multiplicand.
46 * Handle simple cases with traditional multiplication.
48 * This is the most critical code of multiplication. All multiplies rely
49 * on this, both small and huge. Small ones arrive here immediately. Huge
50 * ones arrive here as this is the base case for Karatsuba's recursive
55 mul_n_basecase(mpi_ptr_t prodp
, mpi_ptr_t up
, mpi_ptr_t vp
, mpi_size_t size
)
61 /* Multiply by the first limb in V separately, as the result can be
62 * stored (not added) to PROD. We also avoid a loop for zeroing. */
66 MPN_COPY(prodp
, up
, size
);
68 MPN_ZERO(prodp
, size
);
71 cy
= mpihelp_mul_1(prodp
, up
, size
, v_limb
);
76 /* For each iteration in the outer loop, multiply one limb from
77 * U with one limb from V, and add it to PROD. */
78 for (i
= 1; i
< size
; i
++) {
83 cy
= mpihelp_add_n(prodp
, prodp
, up
, size
);
85 cy
= mpihelp_addmul_1(prodp
, up
, size
, v_limb
);
95 mul_n(mpi_ptr_t prodp
, mpi_ptr_t up
, mpi_ptr_t vp
,
96 mpi_size_t size
, mpi_ptr_t tspace
)
99 /* The size is odd, and the code below doesn't handle that.
100 * Multiply the least significant (size - 1) limbs with a recursive
101 * call, and handle the most significant limb of S1 and S2
103 * A slightly faster way to do this would be to make the Karatsuba
104 * code below behave as if the size were even, and let it check for
105 * odd size in the end. I.e., in essence move this code to the end.
106 * Doing so would save us a recursive call, and potentially make the
107 * stack grow a lot less.
109 mpi_size_t esize
= size
- 1; /* even size */
112 MPN_MUL_N_RECURSE(prodp
, up
, vp
, esize
, tspace
);
113 cy_limb
= mpihelp_addmul_1(prodp
+ esize
, up
, esize
, vp
[esize
]);
114 prodp
[esize
+ esize
] = cy_limb
;
115 cy_limb
= mpihelp_addmul_1(prodp
+ esize
, vp
, size
, up
[esize
]);
116 prodp
[esize
+ size
] = cy_limb
;
118 /* Anatolij Alekseevich Karatsuba's divide-and-conquer algorithm.
120 * Split U in two pieces, U1 and U0, such that
121 * U = U0 + U1*(B**n),
122 * and V in V1 and V0, such that
123 * V = V0 + V1*(B**n).
125 * UV is then computed recursively using the identity
128 * UV = (B + B )U V + B (U -U )(V -V ) + (B + 1)U V
131 * Where B = 2**BITS_PER_MP_LIMB.
133 mpi_size_t hsize
= size
>> 1;
137 /* Product H. ________________ ________________
138 * |_____U1 x V1____||____U0 x V0_____|
139 * Put result in upper part of PROD and pass low part of TSPACE
142 MPN_MUL_N_RECURSE(prodp
+ size
, up
+ hsize
, vp
+ hsize
, hsize
,
145 /* Product M. ________________
148 if (mpihelp_cmp(up
+ hsize
, up
, hsize
) >= 0) {
149 mpihelp_sub_n(prodp
, up
+ hsize
, up
, hsize
);
152 mpihelp_sub_n(prodp
, up
, up
+ hsize
, hsize
);
155 if (mpihelp_cmp(vp
+ hsize
, vp
, hsize
) >= 0) {
156 mpihelp_sub_n(prodp
+ hsize
, vp
+ hsize
, vp
, hsize
);
159 mpihelp_sub_n(prodp
+ hsize
, vp
, vp
+ hsize
, hsize
);
160 /* No change of NEGFLG. */
162 /* Read temporary operands from low part of PROD.
163 * Put result in low part of TSPACE using upper part of TSPACE
166 MPN_MUL_N_RECURSE(tspace
, prodp
, prodp
+ hsize
, hsize
,
169 /* Add/copy product H. */
170 MPN_COPY(prodp
+ hsize
, prodp
+ size
, hsize
);
171 cy
= mpihelp_add_n(prodp
+ size
, prodp
+ size
,
172 prodp
+ size
+ hsize
, hsize
);
174 /* Add product M (if NEGFLG M is a negative number) */
177 mpihelp_sub_n(prodp
+ hsize
, prodp
+ hsize
, tspace
,
181 mpihelp_add_n(prodp
+ hsize
, prodp
+ hsize
, tspace
,
184 /* Product L. ________________ ________________
185 * |________________||____U0 x V0_____|
186 * Read temporary operands from low part of PROD.
187 * Put result in low part of TSPACE using upper part of TSPACE
190 MPN_MUL_N_RECURSE(tspace
, up
, vp
, hsize
, tspace
+ size
);
192 /* Add/copy Product L (twice) */
194 cy
+= mpihelp_add_n(prodp
+ hsize
, prodp
+ hsize
, tspace
, size
);
196 mpihelp_add_1(prodp
+ hsize
+ size
,
197 prodp
+ hsize
+ size
, hsize
, cy
);
199 MPN_COPY(prodp
, tspace
, hsize
);
200 cy
= mpihelp_add_n(prodp
+ hsize
, prodp
+ hsize
, tspace
+ hsize
,
203 mpihelp_add_1(prodp
+ size
, prodp
+ size
, size
, 1);
207 void mpih_sqr_n_basecase(mpi_ptr_t prodp
, mpi_ptr_t up
, mpi_size_t size
)
213 /* Multiply by the first limb in V separately, as the result can be
214 * stored (not added) to PROD. We also avoid a loop for zeroing. */
218 MPN_COPY(prodp
, up
, size
);
220 MPN_ZERO(prodp
, size
);
223 cy_limb
= mpihelp_mul_1(prodp
, up
, size
, v_limb
);
225 prodp
[size
] = cy_limb
;
228 /* For each iteration in the outer loop, multiply one limb from
229 * U with one limb from V, and add it to PROD. */
230 for (i
= 1; i
< size
; i
++) {
235 cy_limb
= mpihelp_add_n(prodp
, prodp
, up
, size
);
237 cy_limb
= mpihelp_addmul_1(prodp
, up
, size
, v_limb
);
239 prodp
[size
] = cy_limb
;
245 mpih_sqr_n(mpi_ptr_t prodp
, mpi_ptr_t up
, mpi_size_t size
, mpi_ptr_t tspace
)
248 /* The size is odd, and the code below doesn't handle that.
249 * Multiply the least significant (size - 1) limbs with a recursive
250 * call, and handle the most significant limb of S1 and S2
252 * A slightly faster way to do this would be to make the Karatsuba
253 * code below behave as if the size were even, and let it check for
254 * odd size in the end. I.e., in essence move this code to the end.
255 * Doing so would save us a recursive call, and potentially make the
256 * stack grow a lot less.
258 mpi_size_t esize
= size
- 1; /* even size */
261 MPN_SQR_N_RECURSE(prodp
, up
, esize
, tspace
);
262 cy_limb
= mpihelp_addmul_1(prodp
+ esize
, up
, esize
, up
[esize
]);
263 prodp
[esize
+ esize
] = cy_limb
;
264 cy_limb
= mpihelp_addmul_1(prodp
+ esize
, up
, size
, up
[esize
]);
266 prodp
[esize
+ size
] = cy_limb
;
268 mpi_size_t hsize
= size
>> 1;
271 /* Product H. ________________ ________________
272 * |_____U1 x U1____||____U0 x U0_____|
273 * Put result in upper part of PROD and pass low part of TSPACE
276 MPN_SQR_N_RECURSE(prodp
+ size
, up
+ hsize
, hsize
, tspace
);
278 /* Product M. ________________
281 if (mpihelp_cmp(up
+ hsize
, up
, hsize
) >= 0)
282 mpihelp_sub_n(prodp
, up
+ hsize
, up
, hsize
);
284 mpihelp_sub_n(prodp
, up
, up
+ hsize
, hsize
);
286 /* Read temporary operands from low part of PROD.
287 * Put result in low part of TSPACE using upper part of TSPACE
289 MPN_SQR_N_RECURSE(tspace
, prodp
, hsize
, tspace
+ size
);
291 /* Add/copy product H */
292 MPN_COPY(prodp
+ hsize
, prodp
+ size
, hsize
);
293 cy
= mpihelp_add_n(prodp
+ size
, prodp
+ size
,
294 prodp
+ size
+ hsize
, hsize
);
296 /* Add product M (if NEGFLG M is a negative number). */
297 cy
-= mpihelp_sub_n(prodp
+ hsize
, prodp
+ hsize
, tspace
, size
);
299 /* Product L. ________________ ________________
300 * |________________||____U0 x U0_____|
301 * Read temporary operands from low part of PROD.
302 * Put result in low part of TSPACE using upper part of TSPACE
304 MPN_SQR_N_RECURSE(tspace
, up
, hsize
, tspace
+ size
);
306 /* Add/copy Product L (twice). */
307 cy
+= mpihelp_add_n(prodp
+ hsize
, prodp
+ hsize
, tspace
, size
);
309 mpihelp_add_1(prodp
+ hsize
+ size
,
310 prodp
+ hsize
+ size
, hsize
, cy
);
312 MPN_COPY(prodp
, tspace
, hsize
);
313 cy
= mpihelp_add_n(prodp
+ hsize
, prodp
+ hsize
, tspace
+ hsize
,
316 mpihelp_add_1(prodp
+ size
, prodp
+ size
, size
, 1);
321 mpihelp_mul_karatsuba_case(mpi_ptr_t prodp
,
322 mpi_ptr_t up
, mpi_size_t usize
,
323 mpi_ptr_t vp
, mpi_size_t vsize
,
324 struct karatsuba_ctx
*ctx
)
328 if (!ctx
->tspace
|| ctx
->tspace_size
< vsize
) {
330 mpi_free_limb_space(ctx
->tspace
);
331 ctx
->tspace
= mpi_alloc_limb_space(2 * vsize
);
334 ctx
->tspace_size
= vsize
;
337 MPN_MUL_N_RECURSE(prodp
, up
, vp
, vsize
, ctx
->tspace
);
342 if (usize
>= vsize
) {
343 if (!ctx
->tp
|| ctx
->tp_size
< vsize
) {
345 mpi_free_limb_space(ctx
->tp
);
346 ctx
->tp
= mpi_alloc_limb_space(2 * vsize
);
349 mpi_free_limb_space(ctx
->tspace
);
353 ctx
->tp_size
= vsize
;
357 MPN_MUL_N_RECURSE(ctx
->tp
, up
, vp
, vsize
, ctx
->tspace
);
358 cy
= mpihelp_add_n(prodp
, prodp
, ctx
->tp
, vsize
);
359 mpihelp_add_1(prodp
+ vsize
, ctx
->tp
+ vsize
, vsize
,
364 } while (usize
>= vsize
);
368 if (usize
< KARATSUBA_THRESHOLD
) {
370 if (mpihelp_mul(ctx
->tspace
, vp
, vsize
, up
, usize
, &tmp
)
375 ctx
->next
= kzalloc(sizeof *ctx
, GFP_KERNEL
);
379 if (mpihelp_mul_karatsuba_case(ctx
->tspace
,
386 cy
= mpihelp_add_n(prodp
, prodp
, ctx
->tspace
, vsize
);
387 mpihelp_add_1(prodp
+ vsize
, ctx
->tspace
+ vsize
, usize
, cy
);
393 void mpihelp_release_karatsuba_ctx(struct karatsuba_ctx
*ctx
)
395 struct karatsuba_ctx
*ctx2
;
398 mpi_free_limb_space(ctx
->tp
);
400 mpi_free_limb_space(ctx
->tspace
);
401 for (ctx
= ctx
->next
; ctx
; ctx
= ctx2
) {
404 mpi_free_limb_space(ctx
->tp
);
406 mpi_free_limb_space(ctx
->tspace
);
411 /* Multiply the natural numbers u (pointed to by UP, with USIZE limbs)
412 * and v (pointed to by VP, with VSIZE limbs), and store the result at
413 * PRODP. USIZE + VSIZE limbs are always stored, but if the input
414 * operands are normalized. Return the most significant limb of the
417 * NOTE: The space pointed to by PRODP is overwritten before finished
418 * with U and V, so overlap is an error.
420 * Argument constraints:
422 * 2. PRODP != UP and PRODP != VP, i.e. the destination
423 * must be distinct from the multiplier and the multiplicand.
427 mpihelp_mul(mpi_ptr_t prodp
, mpi_ptr_t up
, mpi_size_t usize
,
428 mpi_ptr_t vp
, mpi_size_t vsize
, mpi_limb_t
*_result
)
430 mpi_ptr_t prod_endp
= prodp
+ usize
+ vsize
- 1;
432 struct karatsuba_ctx ctx
;
434 if (vsize
< KARATSUBA_THRESHOLD
) {
443 /* Multiply by the first limb in V separately, as the result can be
444 * stored (not added) to PROD. We also avoid a loop for zeroing. */
448 MPN_COPY(prodp
, up
, usize
);
450 MPN_ZERO(prodp
, usize
);
453 cy
= mpihelp_mul_1(prodp
, up
, usize
, v_limb
);
458 /* For each iteration in the outer loop, multiply one limb from
459 * U with one limb from V, and add it to PROD. */
460 for (i
= 1; i
< vsize
; i
++) {
465 cy
= mpihelp_add_n(prodp
, prodp
, up
,
468 cy
= mpihelp_addmul_1(prodp
, up
, usize
, v_limb
);
478 memset(&ctx
, 0, sizeof ctx
);
479 if (mpihelp_mul_karatsuba_case(prodp
, up
, usize
, vp
, vsize
, &ctx
) < 0)
481 mpihelp_release_karatsuba_ctx(&ctx
);
482 *_result
= *prod_endp
;