LAMMP 4.2.0
Lamina High-Precision Arithmetic Library
载入中...
搜索中...
未找到
mul_toom62.c 文件参考
+ mul_toom62.c 的引用(Include)关系图:

浏览源代码.

结构体

struct  toom62_cache_t
 

宏定义

#define a0   numa
 
#define a0   numa
 
#define a0   numa
 
#define a1   (numa + n)
 
#define a1   (numa + n)
 
#define a1   (numa + n)
 
#define a2   (numa + 2 * n)
 
#define a2   (numa + 2 * n)
 
#define a2   (numa + 2 * n)
 
#define a3   (numa + 3 * n)
 
#define a3   (numa + 3 * n)
 
#define a3   (numa + 3 * n)
 
#define a4   (numa + 4 * n)
 
#define a4   (numa + 4 * n)
 
#define a4   (numa + 4 * n)
 
#define a5   (numa + 5 * n)
 
#define a5   (numa + 5 * n)
 
#define a5   (numa + 5 * n)
 
#define b0   numb
 
#define b0   numb
 
#define b0   numb
 
#define b1   (numb + n)
 
#define b1   (numb + n)
 
#define b1   (numb + n)
 
#define bs1   (cache->bs1)
 
#define bs1   (cache->bs1)
 
#define bs2   (cache->bs2)
 
#define bs2   (cache->bs2)
 
#define bsh   (cache->bsh)
 
#define bsh   (cache->bsh)
 
#define bsm1   (cache->bsm1)
 
#define bsm1   (cache->bsm1)
 
#define bsm2   (cache->bsm2)
 
#define bsm2   (cache->bsm2)
 
#define lmmp_mul_n_(dst, numa, numb, n)
 Copyright (C) 2026 HJimmyK(Jericho Knox)
 
#define n   (cache->n)
 
#define n   (cache->n)
 
#define numb   (cache->numb)
 
#define numb   (cache->numb)
 
#define s   (cache->s)
 
#define s   (cache->s)
 
#define scratch   (cache->scratch)
 
#define scratch   (cache->scratch)
 
#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */
 
#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */
 
#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */
 
#define t   (cache->t)
 
#define t   (cache->t)
 
#define tmp   (cache->tmp)
 
#define tmp   (cache->tmp)
 
#define v0   dst /* 2n */
 
#define v0   dst /* 2n */
 
#define v0   dst /* 2n */
 
#define v1   (dst + 2 * n) /* 2n+1 */
 
#define v1   (dst + 2 * n) /* 2n+1 */
 
#define v1   (dst + 2 * n) /* 2n+1 */
 
#define v2   scratch /* 2n+1 */
 
#define v2   scratch /* 2n+1 */
 
#define v2   scratch /* 2n+1 */
 
#define vh   (scratch + 4 * n + 2) /* 2n+1 */
 
#define vh   (scratch + 4 * n + 2) /* 2n+1 */
 
#define vh   (scratch + 4 * n + 2) /* 2n+1 */
 
#define vinf   (dst + 6 * n) /* s+t */
 
#define vinf   (dst + 6 * n) /* s+t */
 
#define vinf   (dst + 6 * n) /* s+t */
 
#define vm1   (scratch + 6 * n + 3) /* 2n+1 */
 
#define vm1   (scratch + 6 * n + 3) /* 2n+1 */
 
#define vm1   (scratch + 6 * n + 3) /* 2n+1 */
 
#define vm2   (scratch + 2 * n + 1) /* 2n+1 */
 
#define vm2   (scratch + 2 * n + 1) /* 2n+1 */
 
#define vm2   (scratch + 2 * n + 1) /* 2n+1 */
 

函数

void lmmp_mul_toom62_ (mp_ptr restrict dst, mp_srcptr restrict numa, mp_size_t na, mp_srcptr restrict numb, mp_size_t nb)
 
static void lmmp_mul_toom62_cache_ (mp_ptr restrict dst, mp_srcptr restrict numa, const toom62_cache_t *cache, enum toom7_flags bflags)
 
static enum toom7_flags lmmp_mul_toom62_cache_init_ (mp_ptr restrict dst, mp_srcptr restrict numa, toom62_cache_t *cache)
 
void lmmp_mul_toom62_unbalance_ (mp_ptr restrict dst, mp_srcptr restrict numa, mp_size_t na, mp_srcptr restrict numb, mp_size_t nb)
 

结构体说明

◆ toom62_cache_t

struct toom62_cache_t

在文件 mul_toom62.c231 行定义.

+ toom62_cache_t 的协作图:
成员变量
mp_ptr restrict bs1
mp_ptr restrict bs2
mp_ptr restrict bsh
mp_ptr restrict bsm1
mp_ptr restrict bsm2
mp_size_t n
mp_srcptr restrict numb
mp_size_t s
mp_ptr restrict scratch
mp_size_t t
mp_ptr restrict tmp

宏定义说明

◆ a0 [1/3]

#define a0   numa

◆ a0 [2/3]

#define a0   numa

◆ a0 [3/3]

#define a0   numa

◆ a1 [1/3]

#define a1   (numa + n)

◆ a1 [2/3]

#define a1   (numa + n)

◆ a1 [3/3]

#define a1   (numa + n)

◆ a2 [1/3]

#define a2   (numa + 2 * n)

◆ a2 [2/3]

#define a2   (numa + 2 * n)

◆ a2 [3/3]

#define a2   (numa + 2 * n)

◆ a3 [1/3]

#define a3   (numa + 3 * n)

◆ a3 [2/3]

#define a3   (numa + 3 * n)

◆ a3 [3/3]

#define a3   (numa + 3 * n)

◆ a4 [1/3]

#define a4   (numa + 4 * n)

◆ a4 [2/3]

#define a4   (numa + 4 * n)

◆ a4 [3/3]

#define a4   (numa + 4 * n)

◆ a5 [1/3]

#define a5   (numa + 5 * n)

◆ a5 [2/3]

#define a5   (numa + 5 * n)

◆ a5 [3/3]

#define a5   (numa + 5 * n)

◆ b0 [1/3]

#define b0   numb

◆ b0 [2/3]

#define b0   numb

◆ b0 [3/3]

#define b0   numb

◆ b1 [1/3]

#define b1   (numb + n)

◆ b1 [2/3]

#define b1   (numb + n)

◆ b1 [3/3]

#define b1   (numb + n)

◆ bs1 [1/2]

#define bs1   (cache->bs1)

◆ bs1 [2/2]

#define bs1   (cache->bs1)

◆ bs2 [1/2]

#define bs2   (cache->bs2)

◆ bs2 [2/2]

#define bs2   (cache->bs2)

◆ bsh [1/2]

#define bsh   (cache->bsh)

◆ bsh [2/2]

#define bsh   (cache->bsh)

◆ bsm1 [1/2]

#define bsm1   (cache->bsm1)

◆ bsm1 [2/2]

#define bsm1   (cache->bsm1)

◆ bsm2 [1/2]

#define bsm2   (cache->bsm2)

◆ bsm2 [2/2]

#define bsm2   (cache->bsm2)

◆ lmmp_mul_n_

#define lmmp_mul_n_ (   dst,
  numa,
  numb,
  n 
)
值:
lmmp_mul_basecase_((dst), (numa), (n), (numb), (n)); \
else if ((n) < MUL_TOOM33_THRESHOLD) \
lmmp_mul_toom22_((dst), (numa), (n), (numb), (n)); \
else if ((n) < MUL_TOOM44_THRESHOLD) \
lmmp_mul_toom33_((dst), (numa), (n), (numb), (n)); \
lmmp_mul_toom44_((dst), (numa), (n), (numb), (n))
#define MUL_TOOM22_THRESHOLD
Definition mparam.h:46
#define MUL_TOOM33_THRESHOLD
Definition mparam.h:50
#define MUL_TOOM44_THRESHOLD
Definition mparam.h:52
#define numb
#define n

Copyright (C) 2026 HJimmyK(Jericho Knox)

This file is part of LAMMP.

LAMMP is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License (LGPL) as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.

This program is distributed WITHOUT ANY WARRANTY.

See https://www.gnu.org/licenses/.

在文件 mul_toom62.c21 行定义.

50 {
51 lmmp_param_assert(na >= 3 * nb);
52 lmmp_param_assert(5 * nb >= na);
53
54 mp_size_t n, s, t;
58 mp_ptr gp;
61
62#define a0 numa
63#define a1 (numa + n)
64#define a2 (numa + 2 * n)
65#define a3 (numa + 3 * n)
66#define a4 (numa + 4 * n)
67#define a5 (numa + 5 * n)
68#define b0 numb
69#define b1 (numb + n)
70
71 n = 1 + (na >= 3 * nb ? (na - 1) / (mp_size_t)6 : (nb - 1) >> 1);
72
73 s = na - 5 * n;
74 t = nb - n;
75
76 lmmp_debug_assert(0 < s && s <= n);
77 lmmp_debug_assert(0 < t && t <= n);
78
80
82 as1 = tmp;
83 asm1 = as1 + n + 1;
84 as2 = asm1 + n + 1;
85 asm2 = as2 + n + 1;
86 ash = asm2 + n + 1;
87 bs1 = ash + n + 1;
88 bsm1 = bs1 + n + 1;
89 bs2 = bsm1 + n;
90 bsm2 = bs2 + n + 1;
91 bsh = bsm2 + n + 1;
92
93 gp = dst;
94
95 /* Compute as1 and asm1. */
97
98 /* Compute as2 and asm2. */
100
101 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
102 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
103
104 cy = lmmp_addshl1_n_(ash, a1, a0, n);
105 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
106 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
107 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
108 if (s < n) {
111 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
112 lmmp_inc_1(ash + s, cy2);
113 } else
114 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
115
116 /* Compute bs1 and bsm1. */
117 if (t == n) {
118 if (lmmp_cmp_(b0, b1, n) < 0) {
121 } else {
123 bflags = (enum toom7_flags)0;
124 }
125 bs1[n] = cy >> 1;
126 } else {
127 bs1[n] = lmmp_add_(bs1, b0, n, b1, t);
128 if (lmmp_zero_q_(b0 + t, n - t) && lmmp_cmp_(b0, b1, t) < 0) {
129 lmmp_sub_n_(bsm1, b1, b0, t);
130 lmmp_zero(bsm1 + t, n - t);
132 } else {
133 lmmp_sub_(bsm1, b0, n, b1, t);
134 bflags = (enum toom7_flags)0;
135 }
136 }
137
138 /* Compute bs2 and bsm2. Recycling bs1 and bsm1; bs2=bs1+b1, bsm2 =
139 bsm1 - b1 */
140 lmmp_add_(bs2, bs1, n + 1, b1, t);
141 if (bflags & toom7_w3_neg) {
142 bsm2[n] = lmmp_add_(bsm2, bsm1, n, b1, t);
144 } else {
145 if (t < n) {
146 if (lmmp_zero_q_(bsm1 + t, n - t) && lmmp_cmp_(bsm1, b1, t) < 0) {
148 lmmp_zero(bsm2 + t, n + 1 - t);
150 } else {
151 lmmp_sub_(bsm2, bsm1, n, b1, t);
152 bsm2[n] = 0;
153 }
154 } else {
155 if (lmmp_cmp_(bsm1, b1, n) < 0) {
158 } else {
160 }
161 bsm2[n] = 0;
162 }
163 }
164
165 /* Compute bsh, recycling bs1. bsh=bs1+b0; */
166 bsh[n] = bs1[n] + lmmp_add_n_(bsh, bs1, b0, n);
167
168 lmmp_debug_assert(as1[n] <= 5);
169 lmmp_debug_assert(bs1[n] <= 1);
170 lmmp_debug_assert(asm1[n] <= 2);
171 lmmp_debug_assert(as2[n] <= 62);
172 lmmp_debug_assert(bs2[n] <= 2);
173 lmmp_debug_assert(asm2[n] <= 41);
174 lmmp_debug_assert(bsm2[n] <= 1);
175 lmmp_debug_assert(ash[n] <= 62);
176 lmmp_debug_assert(bsh[n] <= 2);
177
178#define v0 dst /* 2n */
179#define v1 (dst + 2 * n) /* 2n+1 */
180#define vinf (dst + 6 * n) /* s+t */
181#define v2 scratch /* 2n+1 */
182#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
183#define vh (scratch + 4 * n + 2) /* 2n+1 */
184#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
185#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
186 /* Total scratch need: 10*n+5 */
187
188 /* Must be in allocation order, as they overwrite one limb beyond
189 * 2n+1. */
190 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
191 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
192 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
193
194 /* vm1, 2n+1 limbs */
196 cy = 0;
197 if (asm1[n] == 1) {
198 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
199 } else if (asm1[n] == 2) {
200 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
201 }
202 vm1[2 * n] = cy;
203
204 /* v1, 2n+1 limbs */
205 lmmp_mul_n_(v1, as1, bs1, n);
206 if (as1[n] == 1) {
207 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
208 } else if (as1[n] == 2) {
209 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
210 } else if (as1[n] != 0) {
211 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
212 } else
213 cy = 0;
214 if (bs1[n] != 0)
215 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
216 v1[2 * n] = cy;
217
218 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
219
220 /* vinf, s+t limbs */
221 if (s > t)
222 lmmp_mul_(vinf, a5, s, b1, t);
223 else
224 lmmp_mul_(vinf, b1, t, a5, s);
225
227
229}
230
231typedef struct {
233 mp_size_t n;
234 mp_size_t s;
235 mp_size_t t;
244
249) {
250#define numb (cache->numb)
251#define n (cache->n)
252#define s (cache->s)
253#define t (cache->t)
254#define scratch (cache->scratch)
255#define tmp (cache->tmp)
256#define bs1 (cache->bs1)
257#define bsm1 (cache->bsm1)
258#define bs2 (cache->bs2)
259#define bsm2 (cache->bsm2)
260#define bsh (cache->bsh)
261
265
266#define a0 numa
267#define a1 (numa + n)
268#define a2 (numa + 2 * n)
269#define a3 (numa + 3 * n)
270#define a4 (numa + 4 * n)
271#define a5 (numa + 5 * n)
272#define b0 numb
273#define b1 (numb + n)
274
275
276 as1 = tmp;
277 asm1 = as1 + n + 1;
278 as2 = asm1 + n + 1;
279 asm2 = as2 + n + 1;
280 ash = asm2 + n + 1;
281
282
283 /* Compute as1 and asm1. */
285
286 /* Compute as2 and asm2. */
288
289 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
290 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
291
292 cy = lmmp_addshl1_n_(ash, a1, a0, n);
293 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
294 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
295 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
296 if (s < n) {
299 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
300 lmmp_inc_1(ash + s, cy2);
301 } else
302 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
303
304 /* Compute bs1 and bsm1. */
305 if (t == n) {
306 if (lmmp_cmp_(b0, b1, n) < 0) {
309 } else {
311 bflags = (enum toom7_flags)0;
312 }
313 bs1[n] = cy >> 1;
314 } else {
315 bs1[n] = lmmp_add_(bs1, b0, n, b1, t);
316 if (lmmp_zero_q_(b0 + t, n - t) && lmmp_cmp_(b0, b1, t) < 0) {
317 lmmp_sub_n_(bsm1, b1, b0, t);
318 lmmp_zero(bsm1 + t, n - t);
320 } else {
321 lmmp_sub_(bsm1, b0, n, b1, t);
322 bflags = (enum toom7_flags)0;
323 }
324 }
325
326 /* Compute bs2 and bsm2. Recycling bs1 and bsm1; bs2=bs1+b1, bsm2 =
327 bsm1 - b1 */
328 lmmp_add_(bs2, bs1, n + 1, b1, t);
329 if (bflags & toom7_w3_neg) {
330 bsm2[n] = lmmp_add_(bsm2, bsm1, n, b1, t);
332 } else {
333 if (t < n) {
334 if (lmmp_zero_q_(bsm1 + t, n - t) && lmmp_cmp_(bsm1, b1, t) < 0) {
336 lmmp_zero(bsm2 + t, n + 1 - t);
338 } else {
339 lmmp_sub_(bsm2, bsm1, n, b1, t);
340 bsm2[n] = 0;
341 }
342 } else {
343 if (lmmp_cmp_(bsm1, b1, n) < 0) {
346 } else {
348 }
349 bsm2[n] = 0;
350 }
351 }
352
353 /* Compute bsh, recycling bs1. bsh=bs1+b0; */
354 bsh[n] = bs1[n] + lmmp_add_n_(bsh, bs1, b0, n);
355
356 lmmp_debug_assert(as1[n] <= 5);
357 lmmp_debug_assert(bs1[n] <= 1);
358 lmmp_debug_assert(asm1[n] <= 2);
359 lmmp_debug_assert(as2[n] <= 62);
360 lmmp_debug_assert(bs2[n] <= 2);
361 lmmp_debug_assert(asm2[n] <= 41);
362 lmmp_debug_assert(bsm2[n] <= 1);
363 lmmp_debug_assert(ash[n] <= 62);
364 lmmp_debug_assert(bsh[n] <= 2);
365
366#define v0 dst /* 2n */
367#define v1 (dst + 2 * n) /* 2n+1 */
368#define vinf (dst + 6 * n) /* s+t */
369#define v2 scratch /* 2n+1 */
370#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
371#define vh (scratch + 4 * n + 2) /* 2n+1 */
372#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
373#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
374 /* Total scratch need: 10*n+5 */
375
376 /* Must be in allocation order, as they overwrite one limb beyond
377 * 2n+1. */
378 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
379 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
380 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
381
382 /* vm1, 2n+1 limbs */
384 cy = 0;
385 if (asm1[n] == 1) {
386 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
387 } else if (asm1[n] == 2) {
388 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
389 }
390 vm1[2 * n] = cy;
391
392 /* v1, 2n+1 limbs */
393 lmmp_mul_n_(v1, as1, bs1, n);
394 if (as1[n] == 1) {
395 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
396 } else if (as1[n] == 2) {
397 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
398 } else if (as1[n] != 0) {
399 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
400 } else
401 cy = 0;
402 if (bs1[n] != 0)
403 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
404 v1[2 * n] = cy;
405
406 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
407
408 /* vinf, s+t limbs */
409 if (s > t)
410 lmmp_mul_(vinf, a5, s, b1, t);
411 else
412 lmmp_mul_(vinf, b1, t, a5, s);
413
415
416 return bflags;
417#undef numb
418#undef n
419#undef s
420#undef t
421#undef scratch
422#undef tmp
423#undef bs1
424#undef bsm1
425#undef bs2
426#undef bsm2
427#undef bsh
428}
429
430static void lmmp_mul_toom62_cache_(
433 const toom62_cache_t* cache,
435) {
436#define numb (cache->numb)
437#define n (cache->n)
438#define s (cache->s)
439#define t (cache->t)
440#define scratch (cache->scratch)
441#define tmp (cache->tmp)
442#define bs1 (cache->bs1)
443#define bsm1 (cache->bsm1)
444#define bs2 (cache->bs2)
445#define bsm2 (cache->bsm2)
446#define bsh (cache->bsh)
447
450 enum toom7_flags aflags;
451
452#define a0 numa
453#define a1 (numa + n)
454#define a2 (numa + 2 * n)
455#define a3 (numa + 3 * n)
456#define a4 (numa + 4 * n)
457#define a5 (numa + 5 * n)
458#define b0 numb
459#define b1 (numb + n)
460
461 as1 = tmp;
462 asm1 = as1 + n + 1;
463 as2 = asm1 + n + 1;
464 asm2 = as2 + n + 1;
465 ash = asm2 + n + 1;
466
467
468 /* Compute as1 and asm1. */
470
471 /* Compute as2 and asm2. */
473
474 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
475 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
476
477 cy = lmmp_addshl1_n_(ash, a1, a0, n);
478 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
479 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
480 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
481 if (s < n) {
484 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
485 lmmp_inc_1(ash + s, cy2);
486 } else
487 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
488
489 lmmp_debug_assert(as1[n] <= 5);
490 lmmp_debug_assert(asm1[n] <= 2);
491 lmmp_debug_assert(as2[n] <= 62);
492 lmmp_debug_assert(asm2[n] <= 41);
493 lmmp_debug_assert(ash[n] <= 62);
494
495#define v0 dst /* 2n */
496#define v1 (dst + 2 * n) /* 2n+1 */
497#define vinf (dst + 6 * n) /* s+t */
498#define v2 scratch /* 2n+1 */
499#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
500#define vh (scratch + 4 * n + 2) /* 2n+1 */
501#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
502#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
503 /* Total scratch need: 10*n+5 */
504
505 /* Must be in allocation order, as they overwrite one limb beyond
506 * 2n+1. */
507 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
508 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
509 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
510
511 /* vm1, 2n+1 limbs */
513 cy = 0;
514 if (asm1[n] == 1) {
515 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
516 } else if (asm1[n] == 2) {
517 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
518 }
519 vm1[2 * n] = cy;
520
521 /* v1, 2n+1 limbs */
522 lmmp_mul_n_(v1, as1, bs1, n);
523 if (as1[n] == 1) {
524 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
525 } else if (as1[n] == 2) {
526 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
527 } else if (as1[n] != 0) {
528 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
529 } else
530 cy = 0;
531 if (bs1[n] != 0)
532 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
533 v1[2 * n] = cy;
534
535 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
536
537 /* vinf, s+t limbs */
538 if (s > t)
539 lmmp_mul_(vinf, a5, s, b1, t);
540 else
541 lmmp_mul_(vinf, b1, t, a5, s);
542
544
545#undef numb
546#undef n
547#undef s
548#undef t
549#undef scratch
550#undef tmp
551#undef bs1
552#undef bsm1
553#undef bs2
554#undef bsm2
555#undef bsh
556}
557
564) {
565 lmmp_param_assert(na >= 5 * nb);
566 TEMP_DECL;
568
570 cache.numb = numb;
571 cache.n = 1 + (3 * nb - 1) / (mp_size_t)6;
572 cache.s = 3 * nb - 5 * cache.n;
573 cache.t = nb - cache.n;
574 cache.scratch = BALLOC_TYPE(20 * cache.n + 20, mp_limb_t);
575 cache.tmp = cache.scratch + 10 * cache.n + 10;
576 cache.bs1 = cache.tmp + 5 * cache.n + 5;
577 cache.bsm1 = cache.bs1 + cache.n + 1;
578 cache.bs2 = cache.bsm1 + cache.n;
579 cache.bsm2 = cache.bs2 + cache.n + 1;
580 cache.bsh = cache.bsm2 + cache.n + 1;
581
583 dst += 3 * nb;
584 numa += 3 * nb;
585 na -= 3 * nb;
586 lmmp_copy(ws, dst, nb);
587 while (na >= 5 * nb) {
589 if (lmmp_add_n_(dst, dst, ws, nb))
590 lmmp_inc(dst + nb);
591 dst += 3 * nb;
592 numa += 3 * nb;
593 na -= 3 * nb;
594 lmmp_copy(ws, dst, nb);
595 }
596 // 0 <= na < 2 nb
597 if (na >= nb)
598 lmmp_mul_(dst, numa, na, numb, nb);
599 else
600 lmmp_mul_(dst, numb, nb, numa, na);
601 if (lmmp_add_n_(dst, dst, ws, nb))
602 lmmp_inc(dst + nb);
603 TEMP_FREE;
604}
mp_limb_t * mp_ptr
Definition lmmp.h:80
#define lmmp_copy(dst, src, n)
Definition lmmp.h:367
#define lmmp_zero(dst, n)
Definition lmmp.h:369
uint64_t mp_size_t
Definition lmmp.h:77
#define lmmp_debug_assert(x)
Definition lmmp.h:390
const mp_limb_t * mp_srcptr
Definition lmmp.h:81
uint64_t mp_limb_t
Definition lmmp.h:76
#define lmmp_param_assert(x)
Definition lmmp.h:401
static mp_limb_t lmmp_add_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
大数加法静态内联函数 [dst,na]=[numa,na]+[numb,nb]
Definition lmmpn.h:1050
static int lmmp_cmp_(mp_srcptr numa, mp_srcptr numb, mp_size_t n)
大数比较函数(内联)
Definition lmmpn.h:996
#define lmmp_inc(p)
大数加1宏(预期无进位)
Definition lmmpn.h:938
void lmmp_mul_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
不等长大数乘法操作 [dst,na+nb] = [numa,na] * [numb,nb]
mp_limb_t lmmp_addshl1_n_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
加法结合左移1位操作 [dst,n] = [numa,n] + ([numb,n] << 1)
Definition shl.c:66
mp_limb_t lmmp_shl_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shl)
大数左移操作 [dst,na] = [numa,na]<<shl,dst的低shl位填充0
Definition shl.c:19
static mp_limb_t lmmp_sub_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
大数减法静态内联函数 [dst,na]=[numa,na]-[numb,nb]
Definition lmmpn.h:1064
mp_limb_t lmmp_addmul_1_(mp_ptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t b)
大数乘以单limb并累加操作 [numa,n] += [numb,n] * b
mp_limb_t lmmp_add_n_sub_n_(mp_ptr dsta, mp_ptr dstb, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
同时执行n位加法和减法 ([dsta,n],[dstb,n]) = ([numa,n]+[numb,n],[numa,n]-[numb,n])
Definition add_n_sub_n.c:20
mp_limb_t lmmp_sub_n_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
无借位的n位减法 [dst,n] = [numa,n] - [numb,n]
Definition sub_n.c:80
#define lmmp_inc_1(p, inc)
大数加指定值宏(预期无进位)
Definition lmmpn.h:950
mp_limb_t lmmp_add_n_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
无进位的n位加法 [dst,n] = [numa,n] + [numb,n]
Definition add_n.c:81
static int lmmp_zero_q_(mp_srcptr p, mp_size_t n)
大数判零函数(内联)
Definition lmmpn.h:1019
#define asm1
#define as2
#define asm2
#define as1
#define bs1
mp_srcptr restrict numb
Definition mul_toom62.c:232
#define bs2
#define t
#define b0
#define v0
#define a4
#define a3
#define lmmp_mul_n_(dst, numa, numb, n)
Copyright (C) 2026 HJimmyK(Jericho Knox)
Definition mul_toom62.c:21
#define b1
#define bsm1
#define v2
#define vm1
static void lmmp_mul_toom62_cache_(mp_ptr restrict dst, mp_srcptr restrict numa, const toom62_cache_t *cache, enum toom7_flags bflags)
Definition mul_toom62.c:430
#define a5
#define scratch_out
#define scratch
#define bsm2
#define vh
#define a2
#define a0
void lmmp_mul_toom62_unbalance_(mp_ptr restrict dst, mp_srcptr restrict numa, mp_size_t na, mp_srcptr restrict numb, mp_size_t nb)
Definition mul_toom62.c:558
#define a1
#define tmp
#define s
#define vinf
#define bsh
#define v1
static enum toom7_flags lmmp_mul_toom62_cache_init_(mp_ptr restrict dst, mp_srcptr restrict numa, toom62_cache_t *cache)
Definition mul_toom62.c:245
#define vm2
#define TEMP_DECL
Definition tmp_alloc.h:131
#define TEMP_FREE
Definition tmp_alloc.h:150
#define SALLOC_TYPE(n, type)
Definition tmp_alloc.h:144
#define TEMP_S_DECL
Definition tmp_alloc.h:133
#define BALLOC_TYPE(n, type)
Definition tmp_alloc.h:146
#define TEMP_S_FREE
Definition tmp_alloc.h:166
int lmmp_toom_eval_pm2_(mp_ptr xp2, mp_ptr xm2, unsigned k, mp_srcptr xp, mp_size_t n, mp_size_t hn, mp_ptr tp)
通用高阶 Toom 求值:k次多项式在 x = +2 和 x = -2 处求值
toom7_flags
Definition toom_interp.h:24
@ toom7_w1_neg
Definition toom_interp.h:24
@ toom7_w3_neg
Definition toom_interp.h:24
void lmmp_toom_interp7_(mp_ptr dst, mp_size_t n, enum toom7_flags flags, mp_ptr w1, mp_ptr w3, mp_ptr w4, mp_ptr w5, mp_size_t w6n, mp_ptr tp)
Toom插值计算(7点插值):用于Toom-44、Toom-53、Toom-62 乘法算法
int lmmp_toom_eval_pm1_(mp_ptr xp1, mp_ptr xm1, unsigned k, mp_srcptr xp, mp_size_t n, mp_size_t hn, mp_ptr tp)
通用高阶 Toom 求值:k次多项式在 x = +1 和 x = -1 处求值

◆ n [1/2]

#define n   (cache->n)

◆ n [2/2]

#define n   (cache->n)

◆ numb [1/2]

#define numb   (cache->numb)

◆ numb [2/2]

#define numb   (cache->numb)

◆ s [1/2]

#define s   (cache->s)

◆ s [2/2]

#define s   (cache->s)

◆ scratch [1/2]

#define scratch   (cache->scratch)

◆ scratch [2/2]

#define scratch   (cache->scratch)

◆ scratch_out [1/3]

#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */

◆ scratch_out [2/3]

#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */

◆ scratch_out [3/3]

#define scratch_out   (scratch + 8 * n + 4) /* 2n+1 */

◆ t [1/2]

#define t   (cache->t)

◆ t [2/2]

#define t   (cache->t)

◆ tmp [1/2]

#define tmp   (cache->tmp)

◆ tmp [2/2]

#define tmp   (cache->tmp)

◆ v0 [1/3]

#define v0   dst /* 2n */

◆ v0 [2/3]

#define v0   dst /* 2n */

◆ v0 [3/3]

#define v0   dst /* 2n */

◆ v1 [1/3]

#define v1   (dst + 2 * n) /* 2n+1 */

◆ v1 [2/3]

#define v1   (dst + 2 * n) /* 2n+1 */

◆ v1 [3/3]

#define v1   (dst + 2 * n) /* 2n+1 */

◆ v2 [1/3]

#define v2   scratch /* 2n+1 */

◆ v2 [2/3]

#define v2   scratch /* 2n+1 */

◆ v2 [3/3]

#define v2   scratch /* 2n+1 */

◆ vh [1/3]

#define vh   (scratch + 4 * n + 2) /* 2n+1 */

◆ vh [2/3]

#define vh   (scratch + 4 * n + 2) /* 2n+1 */

◆ vh [3/3]

#define vh   (scratch + 4 * n + 2) /* 2n+1 */

◆ vinf [1/3]

#define vinf   (dst + 6 * n) /* s+t */

◆ vinf [2/3]

#define vinf   (dst + 6 * n) /* s+t */

◆ vinf [3/3]

#define vinf   (dst + 6 * n) /* s+t */

◆ vm1 [1/3]

#define vm1   (scratch + 6 * n + 3) /* 2n+1 */

◆ vm1 [2/3]

#define vm1   (scratch + 6 * n + 3) /* 2n+1 */

◆ vm1 [3/3]

#define vm1   (scratch + 6 * n + 3) /* 2n+1 */

◆ vm2 [1/3]

#define vm2   (scratch + 2 * n + 1) /* 2n+1 */

◆ vm2 [2/3]

#define vm2   (scratch + 2 * n + 1) /* 2n+1 */

◆ vm2 [3/3]

#define vm2   (scratch + 2 * n + 1) /* 2n+1 */

函数说明

◆ lmmp_mul_toom62_()

void lmmp_mul_toom62_ ( mp_ptr restrict  dst,
mp_srcptr restrict  numa,
mp_size_t  na,
mp_srcptr restrict  numb,
mp_size_t  nb 
)

在文件 mul_toom62.c50 行定义.

50 {
51 lmmp_param_assert(na >= 3 * nb);
52 lmmp_param_assert(5 * nb >= na);
53
54 mp_size_t n, s, t;
58 mp_ptr gp;
61
62#define a0 numa
63#define a1 (numa + n)
64#define a2 (numa + 2 * n)
65#define a3 (numa + 3 * n)
66#define a4 (numa + 4 * n)
67#define a5 (numa + 5 * n)
68#define b0 numb
69#define b1 (numb + n)
70
71 n = 1 + (na >= 3 * nb ? (na - 1) / (mp_size_t)6 : (nb - 1) >> 1);
72
73 s = na - 5 * n;
74 t = nb - n;
75
76 lmmp_debug_assert(0 < s && s <= n);
77 lmmp_debug_assert(0 < t && t <= n);
78
80
82 as1 = tmp;
83 asm1 = as1 + n + 1;
84 as2 = asm1 + n + 1;
85 asm2 = as2 + n + 1;
86 ash = asm2 + n + 1;
87 bs1 = ash + n + 1;
88 bsm1 = bs1 + n + 1;
89 bs2 = bsm1 + n;
90 bsm2 = bs2 + n + 1;
91 bsh = bsm2 + n + 1;
92
93 gp = dst;
94
95 /* Compute as1 and asm1. */
97
98 /* Compute as2 and asm2. */
100
101 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
102 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
103
104 cy = lmmp_addshl1_n_(ash, a1, a0, n);
105 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
106 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
107 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
108 if (s < n) {
111 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
112 lmmp_inc_1(ash + s, cy2);
113 } else
114 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
115
116 /* Compute bs1 and bsm1. */
117 if (t == n) {
118 if (lmmp_cmp_(b0, b1, n) < 0) {
121 } else {
123 bflags = (enum toom7_flags)0;
124 }
125 bs1[n] = cy >> 1;
126 } else {
127 bs1[n] = lmmp_add_(bs1, b0, n, b1, t);
128 if (lmmp_zero_q_(b0 + t, n - t) && lmmp_cmp_(b0, b1, t) < 0) {
129 lmmp_sub_n_(bsm1, b1, b0, t);
130 lmmp_zero(bsm1 + t, n - t);
132 } else {
133 lmmp_sub_(bsm1, b0, n, b1, t);
134 bflags = (enum toom7_flags)0;
135 }
136 }
137
138 /* Compute bs2 and bsm2. Recycling bs1 and bsm1; bs2=bs1+b1, bsm2 =
139 bsm1 - b1 */
140 lmmp_add_(bs2, bs1, n + 1, b1, t);
141 if (bflags & toom7_w3_neg) {
142 bsm2[n] = lmmp_add_(bsm2, bsm1, n, b1, t);
144 } else {
145 if (t < n) {
146 if (lmmp_zero_q_(bsm1 + t, n - t) && lmmp_cmp_(bsm1, b1, t) < 0) {
148 lmmp_zero(bsm2 + t, n + 1 - t);
150 } else {
151 lmmp_sub_(bsm2, bsm1, n, b1, t);
152 bsm2[n] = 0;
153 }
154 } else {
155 if (lmmp_cmp_(bsm1, b1, n) < 0) {
158 } else {
160 }
161 bsm2[n] = 0;
162 }
163 }
164
165 /* Compute bsh, recycling bs1. bsh=bs1+b0; */
166 bsh[n] = bs1[n] + lmmp_add_n_(bsh, bs1, b0, n);
167
168 lmmp_debug_assert(as1[n] <= 5);
169 lmmp_debug_assert(bs1[n] <= 1);
170 lmmp_debug_assert(asm1[n] <= 2);
171 lmmp_debug_assert(as2[n] <= 62);
172 lmmp_debug_assert(bs2[n] <= 2);
173 lmmp_debug_assert(asm2[n] <= 41);
174 lmmp_debug_assert(bsm2[n] <= 1);
175 lmmp_debug_assert(ash[n] <= 62);
176 lmmp_debug_assert(bsh[n] <= 2);
177
178#define v0 dst /* 2n */
179#define v1 (dst + 2 * n) /* 2n+1 */
180#define vinf (dst + 6 * n) /* s+t */
181#define v2 scratch /* 2n+1 */
182#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
183#define vh (scratch + 4 * n + 2) /* 2n+1 */
184#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
185#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
186 /* Total scratch need: 10*n+5 */
187
188 /* Must be in allocation order, as they overwrite one limb beyond
189 * 2n+1. */
190 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
191 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
192 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
193
194 /* vm1, 2n+1 limbs */
196 cy = 0;
197 if (asm1[n] == 1) {
198 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
199 } else if (asm1[n] == 2) {
200 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
201 }
202 vm1[2 * n] = cy;
203
204 /* v1, 2n+1 limbs */
205 lmmp_mul_n_(v1, as1, bs1, n);
206 if (as1[n] == 1) {
207 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
208 } else if (as1[n] == 2) {
209 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
210 } else if (as1[n] != 0) {
211 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
212 } else
213 cy = 0;
214 if (bs1[n] != 0)
215 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
216 v1[2 * n] = cy;
217
218 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
219
220 /* vinf, s+t limbs */
221 if (s > t)
222 lmmp_mul_(vinf, a5, s, b1, t);
223 else
224 lmmp_mul_(vinf, b1, t, a5, s);
225
227
229}

引用了 a0, a1, a2, a3, a4, a5, as1, as2, asm1, asm2, b0, b1, bs1, bs2, bsh, bsm1, bsm2, lmmp_add_(), lmmp_add_n_(), lmmp_add_n_sub_n_(), lmmp_addmul_1_(), lmmp_addshl1_n_(), lmmp_cmp_(), lmmp_debug_assert, lmmp_inc_1, lmmp_mul_(), lmmp_mul_n_, lmmp_param_assert, lmmp_shl_(), lmmp_sub_(), lmmp_sub_n_(), lmmp_toom_eval_pm1_(), lmmp_toom_eval_pm2_(), lmmp_toom_interp7_(), lmmp_zero, lmmp_zero_q_(), n, s, SALLOC_TYPE, scratch, scratch_out, t, TEMP_S_DECL, TEMP_S_FREE, tmp, toom7_w1_neg, toom7_w3_neg, v0, v1, v2, vh, vinf, vm1 , 以及 vm2.

+ 函数调用图:

◆ lmmp_mul_toom62_cache_()

static void lmmp_mul_toom62_cache_ ( mp_ptr restrict  dst,
mp_srcptr restrict  numa,
const toom62_cache_t cache,
enum toom7_flags  bflags 
)
static

在文件 mul_toom62.c430 行定义.

435 {
436#define numb (cache->numb)
437#define n (cache->n)
438#define s (cache->s)
439#define t (cache->t)
440#define scratch (cache->scratch)
441#define tmp (cache->tmp)
442#define bs1 (cache->bs1)
443#define bsm1 (cache->bsm1)
444#define bs2 (cache->bs2)
445#define bsm2 (cache->bsm2)
446#define bsh (cache->bsh)
447
450 enum toom7_flags aflags;
451
452#define a0 numa
453#define a1 (numa + n)
454#define a2 (numa + 2 * n)
455#define a3 (numa + 3 * n)
456#define a4 (numa + 4 * n)
457#define a5 (numa + 5 * n)
458#define b0 numb
459#define b1 (numb + n)
460
461 as1 = tmp;
462 asm1 = as1 + n + 1;
463 as2 = asm1 + n + 1;
464 asm2 = as2 + n + 1;
465 ash = asm2 + n + 1;
466
467
468 /* Compute as1 and asm1. */
470
471 /* Compute as2 and asm2. */
473
474 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
475 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
476
477 cy = lmmp_addshl1_n_(ash, a1, a0, n);
478 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
479 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
480 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
481 if (s < n) {
484 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
485 lmmp_inc_1(ash + s, cy2);
486 } else
487 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
488
489 lmmp_debug_assert(as1[n] <= 5);
490 lmmp_debug_assert(asm1[n] <= 2);
491 lmmp_debug_assert(as2[n] <= 62);
492 lmmp_debug_assert(asm2[n] <= 41);
493 lmmp_debug_assert(ash[n] <= 62);
494
495#define v0 dst /* 2n */
496#define v1 (dst + 2 * n) /* 2n+1 */
497#define vinf (dst + 6 * n) /* s+t */
498#define v2 scratch /* 2n+1 */
499#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
500#define vh (scratch + 4 * n + 2) /* 2n+1 */
501#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
502#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
503 /* Total scratch need: 10*n+5 */
504
505 /* Must be in allocation order, as they overwrite one limb beyond
506 * 2n+1. */
507 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
508 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
509 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
510
511 /* vm1, 2n+1 limbs */
513 cy = 0;
514 if (asm1[n] == 1) {
515 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
516 } else if (asm1[n] == 2) {
517 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
518 }
519 vm1[2 * n] = cy;
520
521 /* v1, 2n+1 limbs */
522 lmmp_mul_n_(v1, as1, bs1, n);
523 if (as1[n] == 1) {
524 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
525 } else if (as1[n] == 2) {
526 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
527 } else if (as1[n] != 0) {
528 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
529 } else
530 cy = 0;
531 if (bs1[n] != 0)
532 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
533 v1[2 * n] = cy;
534
535 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
536
537 /* vinf, s+t limbs */
538 if (s > t)
539 lmmp_mul_(vinf, a5, s, b1, t);
540 else
541 lmmp_mul_(vinf, b1, t, a5, s);
542
544
545#undef numb
546#undef n
547#undef s
548#undef t
549#undef scratch
550#undef tmp
551#undef bs1
552#undef bsm1
553#undef bs2
554#undef bsm2
555#undef bsh
556}

引用了 a0, a1, a2, a3, a4, a5, as1, as2, asm1, asm2, b0, b1, bs1, bs2, bsh, bsm1, bsm2, lmmp_add_n_(), lmmp_addmul_1_(), lmmp_addshl1_n_(), lmmp_debug_assert, lmmp_inc_1, lmmp_mul_(), lmmp_mul_n_, lmmp_shl_(), lmmp_toom_eval_pm1_(), lmmp_toom_eval_pm2_(), lmmp_toom_interp7_(), n, s, scratch_out, t, tmp, toom7_w1_neg, toom7_w3_neg, v0, v1, v2, vh, vinf, vm1 , 以及 vm2.

被这些函数引用 lmmp_mul_toom62_unbalance_().

+ 函数调用图:
+ 这是这个函数的调用关系图:

◆ lmmp_mul_toom62_cache_init_()

static enum toom7_flags lmmp_mul_toom62_cache_init_ ( mp_ptr restrict  dst,
mp_srcptr restrict  numa,
toom62_cache_t cache 
)
static

在文件 mul_toom62.c245 行定义.

249 {
250#define numb (cache->numb)
251#define n (cache->n)
252#define s (cache->s)
253#define t (cache->t)
254#define scratch (cache->scratch)
255#define tmp (cache->tmp)
256#define bs1 (cache->bs1)
257#define bsm1 (cache->bsm1)
258#define bs2 (cache->bs2)
259#define bsm2 (cache->bsm2)
260#define bsh (cache->bsh)
261
265
266#define a0 numa
267#define a1 (numa + n)
268#define a2 (numa + 2 * n)
269#define a3 (numa + 3 * n)
270#define a4 (numa + 4 * n)
271#define a5 (numa + 5 * n)
272#define b0 numb
273#define b1 (numb + n)
274
275
276 as1 = tmp;
277 asm1 = as1 + n + 1;
278 as2 = asm1 + n + 1;
279 asm2 = as2 + n + 1;
280 ash = asm2 + n + 1;
281
282
283 /* Compute as1 and asm1. */
285
286 /* Compute as2 and asm2. */
288
289 /* Compute ash = 32 a0 + 16 a1 + 8 a2 + 4 a3 + 2 a4 + a5
290 = 2*(2*(2*(2*(2*a0 + a1) + a2) + a3) + a4) + a5 */
291
292 cy = lmmp_addshl1_n_(ash, a1, a0, n);
293 cy = 2 * cy + lmmp_addshl1_n_(ash, a2, ash, n);
294 cy = 2 * cy + lmmp_addshl1_n_(ash, a3, ash, n);
295 cy = 2 * cy + lmmp_addshl1_n_(ash, a4, ash, n);
296 if (s < n) {
299 ash[n] = 2 * cy + lmmp_shl_(ash + s, ash + s, n - s, 1);
300 lmmp_inc_1(ash + s, cy2);
301 } else
302 ash[n] = 2 * cy + lmmp_addshl1_n_(ash, a5, ash, n);
303
304 /* Compute bs1 and bsm1. */
305 if (t == n) {
306 if (lmmp_cmp_(b0, b1, n) < 0) {
309 } else {
311 bflags = (enum toom7_flags)0;
312 }
313 bs1[n] = cy >> 1;
314 } else {
315 bs1[n] = lmmp_add_(bs1, b0, n, b1, t);
316 if (lmmp_zero_q_(b0 + t, n - t) && lmmp_cmp_(b0, b1, t) < 0) {
317 lmmp_sub_n_(bsm1, b1, b0, t);
318 lmmp_zero(bsm1 + t, n - t);
320 } else {
321 lmmp_sub_(bsm1, b0, n, b1, t);
322 bflags = (enum toom7_flags)0;
323 }
324 }
325
326 /* Compute bs2 and bsm2. Recycling bs1 and bsm1; bs2=bs1+b1, bsm2 =
327 bsm1 - b1 */
328 lmmp_add_(bs2, bs1, n + 1, b1, t);
329 if (bflags & toom7_w3_neg) {
330 bsm2[n] = lmmp_add_(bsm2, bsm1, n, b1, t);
332 } else {
333 if (t < n) {
334 if (lmmp_zero_q_(bsm1 + t, n - t) && lmmp_cmp_(bsm1, b1, t) < 0) {
336 lmmp_zero(bsm2 + t, n + 1 - t);
338 } else {
339 lmmp_sub_(bsm2, bsm1, n, b1, t);
340 bsm2[n] = 0;
341 }
342 } else {
343 if (lmmp_cmp_(bsm1, b1, n) < 0) {
346 } else {
348 }
349 bsm2[n] = 0;
350 }
351 }
352
353 /* Compute bsh, recycling bs1. bsh=bs1+b0; */
354 bsh[n] = bs1[n] + lmmp_add_n_(bsh, bs1, b0, n);
355
356 lmmp_debug_assert(as1[n] <= 5);
357 lmmp_debug_assert(bs1[n] <= 1);
358 lmmp_debug_assert(asm1[n] <= 2);
359 lmmp_debug_assert(as2[n] <= 62);
360 lmmp_debug_assert(bs2[n] <= 2);
361 lmmp_debug_assert(asm2[n] <= 41);
362 lmmp_debug_assert(bsm2[n] <= 1);
363 lmmp_debug_assert(ash[n] <= 62);
364 lmmp_debug_assert(bsh[n] <= 2);
365
366#define v0 dst /* 2n */
367#define v1 (dst + 2 * n) /* 2n+1 */
368#define vinf (dst + 6 * n) /* s+t */
369#define v2 scratch /* 2n+1 */
370#define vm2 (scratch + 2 * n + 1) /* 2n+1 */
371#define vh (scratch + 4 * n + 2) /* 2n+1 */
372#define vm1 (scratch + 6 * n + 3) /* 2n+1 */
373#define scratch_out (scratch + 8 * n + 4) /* 2n+1 */
374 /* Total scratch need: 10*n+5 */
375
376 /* Must be in allocation order, as they overwrite one limb beyond
377 * 2n+1. */
378 lmmp_mul_n_(v2, as2, bs2, n + 1); /* v2, 2n+1 limbs */
379 lmmp_mul_n_(vm2, asm2, bsm2, n + 1); /* vm2, 2n+1 limbs */
380 lmmp_mul_n_(vh, ash, bsh, n + 1); /* vh, 2n+1 limbs */
381
382 /* vm1, 2n+1 limbs */
384 cy = 0;
385 if (asm1[n] == 1) {
386 cy = lmmp_add_n_(vm1 + n, vm1 + n, bsm1, n);
387 } else if (asm1[n] == 2) {
388 cy = lmmp_addshl1_n_(vm1 + n, vm1 + n, bsm1, n);
389 }
390 vm1[2 * n] = cy;
391
392 /* v1, 2n+1 limbs */
393 lmmp_mul_n_(v1, as1, bs1, n);
394 if (as1[n] == 1) {
395 cy = bs1[n] + lmmp_add_n_(v1 + n, v1 + n, bs1, n);
396 } else if (as1[n] == 2) {
397 cy = 2 * bs1[n] + lmmp_addshl1_n_(v1 + n, v1 + n, bs1, n);
398 } else if (as1[n] != 0) {
399 cy = as1[n] * bs1[n] + lmmp_addmul_1_(v1 + n, bs1, n, as1[n]);
400 } else
401 cy = 0;
402 if (bs1[n] != 0)
403 cy += lmmp_add_n_(v1 + n, v1 + n, as1, n);
404 v1[2 * n] = cy;
405
406 lmmp_mul_n_(v0, a0, b0, n); /* v0, 2n limbs */
407
408 /* vinf, s+t limbs */
409 if (s > t)
410 lmmp_mul_(vinf, a5, s, b1, t);
411 else
412 lmmp_mul_(vinf, b1, t, a5, s);
413
415
416 return bflags;
417#undef numb
418#undef n
419#undef s
420#undef t
421#undef scratch
422#undef tmp
423#undef bs1
424#undef bsm1
425#undef bs2
426#undef bsm2
427#undef bsh
428}

引用了 a0, a1, a2, a3, a4, a5, as1, as2, asm1, asm2, b0, b1, bs1, bs2, bsh, bsm1, bsm2, lmmp_add_(), lmmp_add_n_(), lmmp_add_n_sub_n_(), lmmp_addmul_1_(), lmmp_addshl1_n_(), lmmp_cmp_(), lmmp_debug_assert, lmmp_inc_1, lmmp_mul_(), lmmp_mul_n_, lmmp_shl_(), lmmp_sub_(), lmmp_sub_n_(), lmmp_toom_eval_pm1_(), lmmp_toom_eval_pm2_(), lmmp_toom_interp7_(), lmmp_zero, lmmp_zero_q_(), n, s, scratch_out, t, tmp, toom7_w1_neg, toom7_w3_neg, v0, v1, v2, vh, vinf, vm1 , 以及 vm2.

被这些函数引用 lmmp_mul_toom62_unbalance_().

+ 函数调用图:
+ 这是这个函数的调用关系图:

◆ lmmp_mul_toom62_unbalance_()

void lmmp_mul_toom62_unbalance_ ( mp_ptr restrict  dst,
mp_srcptr restrict  numa,
mp_size_t  na,
mp_srcptr restrict  numb,
mp_size_t  nb 
)

在文件 mul_toom62.c558 行定义.

564 {
565 lmmp_param_assert(na >= 5 * nb);
566 TEMP_DECL;
568
570 cache.numb = numb;
571 cache.n = 1 + (3 * nb - 1) / (mp_size_t)6;
572 cache.s = 3 * nb - 5 * cache.n;
573 cache.t = nb - cache.n;
574 cache.scratch = BALLOC_TYPE(20 * cache.n + 20, mp_limb_t);
575 cache.tmp = cache.scratch + 10 * cache.n + 10;
576 cache.bs1 = cache.tmp + 5 * cache.n + 5;
577 cache.bsm1 = cache.bs1 + cache.n + 1;
578 cache.bs2 = cache.bsm1 + cache.n;
579 cache.bsm2 = cache.bs2 + cache.n + 1;
580 cache.bsh = cache.bsm2 + cache.n + 1;
581
583 dst += 3 * nb;
584 numa += 3 * nb;
585 na -= 3 * nb;
586 lmmp_copy(ws, dst, nb);
587 while (na >= 5 * nb) {
589 if (lmmp_add_n_(dst, dst, ws, nb))
590 lmmp_inc(dst + nb);
591 dst += 3 * nb;
592 numa += 3 * nb;
593 na -= 3 * nb;
594 lmmp_copy(ws, dst, nb);
595 }
596 // 0 <= na < 2 nb
597 if (na >= nb)
598 lmmp_mul_(dst, numa, na, numb, nb);
599 else
600 lmmp_mul_(dst, numb, nb, numa, na);
601 if (lmmp_add_n_(dst, dst, ws, nb))
602 lmmp_inc(dst + nb);
603 TEMP_FREE;
604}

引用了 BALLOC_TYPE, lmmp_add_n_(), lmmp_copy, lmmp_inc, lmmp_mul_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_param_assert, n, toom62_cache_t::numb, numb, SALLOC_TYPE, TEMP_DECL , 以及 TEMP_FREE.

+ 函数调用图: