LAMMP 4.1.0
Lamina High-Precision Arithmetic Library
载入中...
搜索中...
未找到
宏定义 | 函数
lmmpn.h 文件参考
#include <stdbool.h>
#include <stdio.h>
#include "lmmp.h"
+ lmmpn.h 的引用(Include)关系图:
+ 此图展示该文件直接或间接的被哪些文件引用了:

浏览源代码.

宏定义

#define INLINE_   static inline
 
#define LMMP_ADDCB_(r, x, y)   ((r) < (y))
 
#define LMMP_AORS_(FUNCTION, TEST)
 
#define LMMP_AORS_1_(OP, CB)
 
#define lmmp_dec(p)
 大数减1宏（预期无借位）
 
#define lmmp_dec_1(p, dec)
 大数减指定值宏（预期无借位）
 
#define lmmp_inc(p)
 大数加1宏（预期无进位）
 
#define lmmp_inc_1(p, inc)
 大数加指定值宏（预期无进位）
 
#define LMMP_SUBCB_(r, x, y)   ((x) < (y))
 

函数

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]
 
static mp_limb_t lmmp_add_1_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
 大数加单精度数静态内联函数 [dst,na]=[numa,na]+x
 
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]
 
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])
 
mp_limb_t lmmp_add_nc_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
 带进位的n位加法 [dst,n] = [numa,n] + [numb,n] + c
 
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_addshl1_n_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 加法结合左移1位操作 [dst,n] = [numa,n] + ([numb,n] << 1)
 
void lmmp_bninv_ (mp_ptr dstq, mp_srcptr numa, mp_size_t na, mp_size_t ni)
 精确逆元计算 [dstq,na+ni+2] = B^(2*(na+ni)) / ([numa,na] * B^ni)
 
static int lmmp_cmp_ (mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 大数比较函数（内联）
 
void lmmp_div_ (mp_ptr dstq, mp_ptr dstr, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 大数除法和取模操作
 
mp_limb_t lmmp_div_1_ (mp_ptr dstq, mp_srcptr numa, mp_size_t na, mp_limb_t x)
 单精度数除法
 
mp_limb_t lmmp_div_1_s_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_limb_t x)
 单精度数除法（除数为1个limb）
 
void lmmp_div_2_ (mp_ptr dstq, mp_srcptr numa, mp_size_t na, mp_ptr numb)
 双精度数除法 (除数为2个limb)
 
mp_limb_t lmmp_div_2_s_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb)
 双精度数除法（除数为2个limb）
 
mp_limb_t lmmp_div_3_2_ (mp_ptr numa, mp_srcptr numb, mp_limb_t inv21)
 3/2位除法运算 [numa,2]=[numa,3] mod [numb,2]
 
mp_limb_t lmmp_div_basecase_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_limb_t inv21)
 基础除法运算
 
mp_limb_t lmmp_div_divide_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_limb_t inv21)
 分治除法运算
 
static mp_size_t lmmp_div_inv_size_ (mp_size_t nq, mp_size_t nb)
 计算预计算逆元的尺寸
 
mp_limb_t lmmp_div_mulinv_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_srcptr invappr, mp_size_t ni)
 乘法逆元除法
 
mp_limb_t lmmp_div_s_ (mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 除法运算
 
static bool lmmp_endian (void)
 运行时判断端序
 
mp_bitcnt_t lmmp_extract_bits_ (mp_srcptr num, mp_size_t n, mp_limb_t *ext, int bits)
 提取高位指定位数，并返回低位bits位数
 
mp_size_t lmmp_fft_next_size_ (mp_size_t n)
 计算满足 >=n 的最小费马/梅森乘法可行尺寸
 
mp_size_t lmmp_from_str_ (mp_ptr dst, const mp_byte_t *src, mp_size_t len, int base)
 字符串转大数操作 [src,len,base] to [dst,return value,B]
 
mp_size_t lmmp_from_str_len_ (const mp_byte_t *src, mp_size_t len, int base)
 计算字符串转大数所需的 limb 缓冲区长度
 
void lmmp_inv_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t nf)
 大数求逆操作 [dst,na+nf+1] = (B^(2*(na+nf)) - 1) / ([numa,na]*B^nf) + [0|-1]
 
mp_limb_t lmmp_inv_1_ (mp_limb_t x)
 1阶逆元计算 (inv1)
 
mp_limb_t lmmp_inv_2_1_ (mp_limb_t xh, mp_limb_t xl)
 2-1阶逆元计算 (inv21)
 
void lmmp_inv_basecase_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 近似逆元计算
 
void lmmp_inv_prediv_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t ni)
 除法前的逆元预计算，[dst,ni] = invappr( (ni+1 MSLs of numa) + 1 ) / B
 
void lmmp_invappr_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 近似逆元计算 (invappr)
 
void lmmp_invappr_newton_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 近似逆元计算（牛顿迭代法）
 
int lmmp_leading_zeros_ (mp_limb_t x)
 计算一个单精度数(limb)中前导零的个数
 
int lmmp_limb_bits_ (mp_limb_t x)
 计算满足 2^k > x 的最小自然数k
 
int lmmp_limb_popcnt_ (mp_limb_t x)
 计算一个64位无符号整数中1的个数
 
mp_limb_t lmmp_mod_1_ (mp_srcptr numa, mp_size_t na, mp_limb_t x)
 单精度数取余
 
void lmmp_mod_2_ (mp_srcptr numa, mp_size_t na, mp_ptr numb)
 双精度数取余 (除数为2个limb)
 
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_mul_1_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
 大数乘以单limb操作 [dst,na] = [numa,na] * x
 
void lmmp_mul_basecase_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 基础乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_fermat_ (mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 费马数模乘法 [dst,rn+1]=[numa,na]*[numb,nb] mod B^rn+1
 
void lmmp_mul_fft_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 FFT乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_fft_unbalance_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 FFT不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_mersenne_ (mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 梅森数模乘法 [dst,rn] = [numa,na]*[numb,nb] mod B^rn-1
 
void lmmp_mul_n_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 等长大数乘法操作 [dst,2*n] = [numa,n] * [numb,n]
 
void lmmp_mul_toom22_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-22乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom32_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-32乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom33_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-33乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom42_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-42乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom42_unbalance_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-42不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom43_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-43乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom44_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-44乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom52_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-52乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom53_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-53乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom62_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-62乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
void lmmp_mul_toom62_unbalance_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
 Toom-62不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
 
mp_limb_t lmmp_mulh_ (mp_limb_t a, mp_limb_t b)
 计算两个64位无符号整数相乘的高位结果 (a*b)/2^64
 
void lmmp_mullh_ (mp_limb_t a, mp_limb_t b, mp_ptr dst)
 计算两个64位无符号整数相乘的128位结果 (a*b)
 
void lmmp_mullo_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 低位乘法 [dst,n] = [numa,n] * [numb,n] mod B^n
 
void lmmp_mullo_dc_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_ptr tp, mp_size_t n)
 低位乘法 [dst,n] = [numa,n] * [numb,n] mod B^n
 
void lmmp_mullo_fft_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_ptr scratch)
 低位FFT乘法 [dst,n] = [numa,n] * [numb,n] mod B^n
 
void lmmp_not_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 大数按位取反操作 [dst,na] = ~[numa,na] (对每个limb执行按位非操作)
 
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
 
mp_limb_t lmmp_shl_c_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shl, mp_limb_t c)
 带进位的大数左移操作 [dst,na] = [numa,na]<<shl，dst的低shl位填充c的低shl位
 
mp_limb_t lmmp_shlnot_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shl)
 左移后按位取反操作 [dst,na] = ~([numa,na] << shl)，dst的低shl位填充1
 
mp_limb_t lmmp_shr1add_n_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 加法后右移1位 [dst,n] = ([numa,n] + [numb,n]) >> 1
 
mp_limb_t lmmp_shr1add_nc_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
 带进位加法后右移1位 [dst,n] = ([numa,n] + [numb,n] + c) >> 1
 
mp_limb_t lmmp_shr1sub_n_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 减法后右移1位 [dst,n] = ([numa,n] - [numb,n]) >> 1
 
mp_limb_t lmmp_shr1sub_nc_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
 带借位减法后右移1位 [dst,n] = ([numa,n] - [numb,n] - c) >> 1
 
mp_limb_t lmmp_shr_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shr)
 大数右移操作 [dst,na] = [numa,na] >> shr，dst的高shr位填充0
 
mp_limb_t lmmp_shr_c_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shr, mp_limb_t c)
 带进位的大数右移操作 [dst,na] = [numa,na]>>shr，dst的高shr位填充c的高shr位
 
void lmmp_sqr_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 大数平方操作 [dst,2*na] = [numa,na]^2
 
void lmmp_sqr_basecase_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 基础平方运算 [dst,2*na] = [numa,na]^2
 
void lmmp_sqr_toom2_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 Toom-2平方运算 [dst,2*na] = [numa,na]^2
 
void lmmp_sqr_toom3_ (mp_ptr dst, mp_srcptr numa, mp_size_t na)
 Toom-3平方运算 [dst,2*na] = [numa,na]^2
 
void lmmp_sqr_toom4_ (mp_ptr pp, mp_srcptr ap, mp_size_t an)
 Toom-4平方运算 [dst,2*na] = [numa,na]^2
 
void lmmp_sqrlo_dc_ (mp_ptr dst, mp_srcptr numa, mp_ptr tp, mp_size_t n)
 低位平方 [dst,n] = [numa,n]^2 mod B^n
 
void lmmp_sqrt_ (mp_ptr dsts, mp_ptr dstr, mp_srcptr numa, mp_size_t na, mp_size_t nf)
 大数平方根和取余操作
 
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]
 
static mp_limb_t lmmp_sub_1_ (mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
 大数减单精度数静态内联函数 [dst,na]=[numa,na]-x
 
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]
 
mp_limb_t lmmp_sub_nc_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
 带借位的n位减法 [dst,n] = [numa,n] - [numb,n] - c
 
mp_limb_t lmmp_submul_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_subshl1_n_ (mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
 减法结合左移1位操作 [dst,n] = [numa,n] - ([numb,n] << 1)
 
int lmmp_tailing_zeros_ (mp_limb_t x)
 计算一个单精度数(limb)中末尾零的个数
 
mp_size_t lmmp_to_str_ (mp_byte_t *dst, mp_srcptr numa, mp_size_t na, int base)
 大数转字符串操作 [numa,na,B] to [dst,return value,base]
 
mp_size_t lmmp_to_str_len_ (mp_srcptr numa, mp_size_t na, int base)
 计算大数转换为字符串，字符串需要的缓冲区长度
 
static int lmmp_zero_q_ (mp_srcptr p, mp_size_t n)
 大数判零函数（内联）
 

宏定义说明

◆ INLINE_

#define INLINE_   static inline

在文件 lmmpn.h 第 59 行定义.

◆ LMMP_ADDCB_

#define LMMP_ADDCB_ (   r,
  x,
  y 
)    ((r) < (y))

在文件 lmmpn.h 第 1098 行定义.

◆ LMMP_AORS_

#define LMMP_AORS_ (   FUNCTION,
  TEST 
)
值:
mp_limb_t _x_; \
if (FUNCTION(dst, numa, numb, nb)) { \
do { \
if (nb >= na) \
return 1; \
_x_ = numa[nb]; \
} while (TEST); \
} \
if (dst != numa && na != nb) \
lmmp_copy(dst + nb, numa + nb, na - nb); \
return 0
mp_limb_t
uint64_t mp_limb_t
Definition lmmp.h:211

在文件 lmmpn.h 第 1035 行定义.

1037 { \
1038 do { \
1039 if (nb >= na) \
1040 return 1; \
1041 _x_ = numa[nb]; \
1042 } while (TEST); \
1043 } \
1044 if (dst != numa && na != nb) \
1045 lmmp_copy(dst + nb, numa + nb, na - nb); \
1046 return 0

◆ LMMP_AORS_1_

#define LMMP_AORS_1_ (   OP,
  CB 
)
值:
mp_size_t _i_ = 1; \
mp_limb_t _x_ = numa[0], _r_ = _x_ OP x; \
dst[0] = _r_; \
if (CB(_r_, _x_, x)) { \
do { \
if (_i_ >= na) \
return 1; \
_x_ = numa[_i_]; \
_r_ = _x_ OP 1; \
dst[_i_] = _r_; \
++_i_; \
} while (CB(_r_, _x_, 1)); \
} \
if (numa != dst && na != _i_) \
lmmp_copy(dst + _i_, numa + _i_, na - _i_); \
return 0
mp_size_t
uint64_t mp_size_t
Definition lmmp.h:212

在文件 lmmpn.h 第 1079 行定义.

1083 { \
1084 do { \
1085 if (_i_ >= na) \
1086 return 1; \
1087 _x_ = numa[_i_]; \
1088 _r_ = _x_ OP 1; \
1089 dst[_i_] = _r_; \
1090 ++_i_; \
1091 } while (CB(_r_, _x_, 1)); \
1092 } \
1093 if (numa != dst && na != _i_) \
1094 lmmp_copy(dst + _i_, numa + _i_, na - _i_); \
1095 return 0

◆ lmmp_dec

#define lmmp_dec (   p)
值:
do { \
mp_ptr _p_ = (p); \
while ((*(_p_++))-- == 0); \
} while (0)
mp_ptr
mp_limb_t * mp_ptr
Definition lmmp.h:215

大数减1宏（预期无借位）

参数
p指向大数起始位置的指针
注解
从最低位开始减1，直到遇到非零值（预期无借位溢出）

在文件 lmmpn.h 第 973 行定义.

974 { \
975 mp_ptr _p_ = (p); \
976 while ((*(_p_++))-- == 0); \
977 } while (0)

◆ lmmp_dec_1

#define lmmp_dec_1 (   p,
  dec 
)
值:
do { \
mp_ptr _p_ = (p); \
mp_limb_t _dec_ = (dec), _x_; \
_x_ = *_p_; \
*_p_ = _x_ - _dec_; \
if (_x_ < _dec_) \
while ((*(++_p_))-- == 0); \
} while (0)

大数减指定值宏（预期无借位）

参数
p指向大数起始位置的指针
dec要减的单精度数值
注解
先减最低位，若产生借位则逐位减1，直到无借位（预期无溢出）

在文件 lmmpn.h 第 985 行定义.

986 { \
987 mp_ptr _p_ = (p); \
988 mp_limb_t _dec_ = (dec), _x_; \
989 _x_ = *_p_; \
990 *_p_ = _x_ - _dec_; \
991 if (_x_ < _dec_) \
992 while ((*(++_p_))-- == 0); \
993 } while (0)

◆ lmmp_inc

#define lmmp_inc (   p)
值:
do { \
mp_ptr _p_ = (p); \
while (++(*(_p_++)) == 0); \
} while (0)

大数加1宏（预期无进位）

参数
p指向大数起始位置的指针
注解
从最低位开始加1，直到遇到非零值（预期无进位溢出）

在文件 lmmpn.h 第 946 行定义.

947 { \
948 mp_ptr _p_ = (p); \
949 while (++(*(_p_++)) == 0); \
950 } while (0)

◆ lmmp_inc_1

#define lmmp_inc_1 (   p,
  inc 
)
值:
do { \
mp_ptr _p_ = (p); \
mp_limb_t _inc_ = (inc), _x_; \
_x_ = *_p_ + _inc_; \
*_p_ = _x_; \
if (_x_ < _inc_) \
while (++(*(++_p_)) == 0); \
} while (0)

大数加指定值宏（预期无进位）

参数
p指向大数起始位置的指针
inc要加的单精度数值
注解
先加最低位，若产生进位则逐位加1，直到无进位（预期无溢出）

在文件 lmmpn.h 第 958 行定义.

959 { \
960 mp_ptr _p_ = (p); \
961 mp_limb_t _inc_ = (inc), _x_; \
962 _x_ = *_p_ + _inc_; \
963 *_p_ = _x_; \
964 if (_x_ < _inc_) \
965 while (++(*(++_p_)) == 0); \
966 } while (0)

◆ LMMP_SUBCB_

#define LMMP_SUBCB_ (   r,
  x,
  y 
)    ((x) < (y))

在文件 lmmpn.h 第 1100 行定义.

函数说明

◆ lmmp_add_()

static mp_limb_t lmmp_add_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)
inlinestatic

大数加法静态内联函数 [dst,na]=[numa,na]+[numb,nb]

参数
dst输出结果缓冲区，存储numa + numb
numa第一个加数，长度为na
na第一个加数的 limb 长度
numb第二个加数，长度为nb
nb第二个加数的 limb 长度
返回
进位标志（1表示有进位，0表示无进位）
警告
0<nb<=na, eqsep(dst,[numa|numb])

在文件 lmmpn.h 第 1058 行定义.

1058 {
1059 LMMP_AORS_(lmmp_add_n_, ((dst[nb++] = _x_ + 1) == 0));
1060}
LMMP_AORS_
#define LMMP_AORS_(FUNCTION, TEST)
Definition lmmpn.h:1035
lmmp_add_n_
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:71

引用了 lmmp_add_n_() , 以及 LMMP_AORS_.

被这些函数引用 lmmp_add_signed_(), lmmp_invsqrt_newton_(), lmmp_mul_fermat_recombine_(), lmmp_mul_fft_(), lmmp_mul_fft_cache_(), lmmp_mul_mersenne_(), lmmp_mul_mersenne_single_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom43_(), lmmp_mul_toom52_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_(), lmmp_mullo_fft_(), lmmp_sqr_toom2_(), lmmp_sqr_toom3_(), lmmp_toom_eval_dgr3_pm1_(), lmmp_toom_eval_dgr3_pm2_() , 以及 lmmp_toom_eval_pm1_().

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

◆ lmmp_add_1_()

static mp_limb_t lmmp_add_1_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_limb_t  x 
)
inlinestatic

大数加单精度数静态内联函数 [dst,na]=[numa,na]+x

参数
dst输出结果缓冲区，存储numa + x
numa被加数，长度为na
na被加数的 limb 长度
x加数（单个 limb ）
返回
进位标志（1表示有进位，0表示无进位）
警告
na>0, eqsep(dst,numa)

在文件 lmmpn.h 第 1111 行定义.

1111{ LMMP_AORS_1_(+, LMMP_ADDCB_); }
LMMP_AORS_1_
#define LMMP_AORS_1_(OP, CB)
Definition lmmpn.h:1079
LMMP_ADDCB_
#define LMMP_ADDCB_(r, x, y)
Definition lmmpn.h:1098

引用了 LMMP_ADDCB_ , 以及 LMMP_AORS_1_.

被这些函数引用 lmmp_from_str_basecase_(), lmmp_inv_prediv_(), lmmp_mul_fermat_recombine_(), lmmp_mul_mersenne_(), lmmp_mul_mersenne_single_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom43_(), lmmp_sqr_toom3_(), lmmp_sqrt_divide_() , 以及 lmmp_toom_eval_pm2_().

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

◆ lmmp_add_n_()

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]

参数
dst结果输出指针
numa第一个加数指针
numb第二个加数指针
nlimb长度
警告
n>0, eqsep(dst,[numa|numb])
返回
运算后的最终进位值 [0|1]

在文件 add_n.c 第 71 行定义.

71 {
72 mp_size_t i = 0;
73 mp_limb_t cy = 0;
74
75 for (; i + 4 <= n; i += 4) {
76 mp_limb_t a0, b0, r0;
77 mp_limb_t a1, b1, r1;
78 mp_limb_t a2, b2, r2;
79 mp_limb_t a3, b3, r3;
80
81 a0 = numa[i + 0];
82 b0 = numb[i + 0];
83
84 a1 = numa[i + 1];
85 b1 = numb[i + 1];
86
87 a2 = numa[i + 2];
88 b2 = numb[i + 2];
89
90 a3 = numa[i + 3];
91 b3 = numb[i + 3];
92
93 r0 = a0 + cy;
94 cy = (r0 < cy);
95 r0 += b0;
96 cy += (r0 < b0);
97
98 r1 = a1 + cy;
99 cy = (r1 < cy);
100 r1 += b1;
101 cy += (r1 < b1);
102
103 r2 = a2 + cy;
104 cy = (r2 < cy);
105 r2 += b2;
106 cy += (r2 < b2);
107
108 r3 = a3 + cy;
109 cy = (r3 < cy);
110 r3 += b3;
111 cy += (r3 < b3);
112
113 dst[i + 0] = r0;
114 dst[i + 1] = r1;
115 dst[i + 2] = r2;
116 dst[i + 3] = r3;
117 }
118
119 for (; i < n; i++) {
120 mp_limb_t a, b, r;
121 a = numa[i];
122 b = numb[i];
123 r = a + cy;
124 cy = (r < cy);
125 r += b;
126 cy += (r < b);
127 dst[i] = r;
128 }
129
130 return cy;
131}
b0
#define b0
b1
#define b1
a0
#define a0
a1
#define a1
r2
#define r2
r1
#define r1
a2
#define a2
r3
#define r3
r0
#define r0
b2
#define b2
a3
#define a3
b3
#define b3

引用了 a0, a1, a2, a3, b0, b1, b2, b3, r0, r1, r2 , 以及 r3.

被这些函数引用 lmmp_add_(), lmmp_binvert_n_dc_(), lmmp_div_(), lmmp_div_basecase_(), lmmp_div_divide_n_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_from_str_divide_(), lmmp_invappr_newton_(), lmmp_mul_(), lmmp_mul_fft_unbalance_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom42_unbalance_(), lmmp_mul_toom43_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_mul_toom62_unbalance_(), lmmp_mullo_dc_(), lmmp_sqr_toom2_(), lmmp_sqrlo_dc_(), lmmp_sqrt_divide_(), lmmp_toom_eval_dgr3_pm1_(), lmmp_toom_eval_dgr3_pm2_(), lmmp_toom_eval_pm1_(), lmmp_toom_interp5_(), lmmp_toom_interp6_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_add_n_sub_n_()

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])

参数
dsta加法结果输出指针
dstb减法结果输出指针
numa第一个操作数指针（被加数/被减数）
numb第二个操作数指针（加数/减数）
nlimb长度
警告
n>0, eqsep(dsta,[numa|numb]), eqsep(dstb,[numa|numb])
返回
组合返回值 cb = 2*c + b (c为加法进位, b为减法借位) 返回值范围: 0(无进位无借位),1(无进位有借位),2(有进位无借位),3(有进位有借位)

在文件 add_n_sub_n.c 第 10 行定义.

10 {
11 /*
12 这段代码看起来有点奇怪的原因是，对于使用x64汇编时，我们会使用带进位的加法和减法，而x64中
13 只能使用同一个进位寄存器，所以我们需要将两条指令分开执行。
14 而不使用汇编时，编译器通常不会使用进位寄存器。因此我们可以同时读取两路内存，以减少读写次数。
15 */
16#ifdef USE_ASM
17 mp_limb_t acyo = 0, scyo = 0;
18 mp_size_t off, this_n;
19
20 if (dsta != numa && dsta != numb) {
21 for (off = 0; off < n; off += PART_SIZE) {
22 this_n = LMMP_MIN(n - off, PART_SIZE);
23 acyo = lmmp_add_nc_(dsta + off, numa + off, numb + off, this_n, acyo);
24 scyo = lmmp_sub_nc_(dstb + off, numa + off, numb + off, this_n, scyo);
25 }
26 } else if (dstb != numa && dstb != numb) {
27 for (off = 0; off < n; off += PART_SIZE) {
28 this_n = LMMP_MIN(n - off, PART_SIZE);
29 scyo = lmmp_sub_nc_(dstb + off, numa + off, numb + off, this_n, scyo);
30 acyo = lmmp_add_nc_(dsta + off, numa + off, numb + off, this_n, acyo);
31 }
32 } else {
33 mp_limb_t tp[PART_SIZE];
34 for (off = 0; off < n; off += PART_SIZE) {
35 this_n = LMMP_MIN(n - off, PART_SIZE);
36 acyo = lmmp_add_nc_(tp, numa + off, numb + off, this_n, acyo);
37 scyo = lmmp_sub_nc_(dstb + off, numa + off, numb + off, this_n, scyo);
38 lmmp_copy(dsta + off, tp, this_n);
39 }
40 }
41 return 2 * acyo + scyo;
42#else
43 mp_size_t i;
44 mp_limb_t acyo, scyo;
45
46 for (i = 0, acyo = 0, scyo = 0; i < n; i++) {
47 mp_limb_t a, b, r;
48 a = numa[i];
49 b = numb[i];
50 r = a + acyo;
51 acyo = (r < acyo);
52 r += b;
53 acyo += (r < b);
54 dsta[i] = r;
55
56 b += scyo;
57 scyo = (b < scyo);
58 scyo += (a < b);
59 dstb[i] = a - b;
60 }
61 return 2 * acyo + scyo;
62#endif
63}
lmmp_copy
#define lmmp_copy(dst, src, n)
Definition lmmp.h:364
LMMP_MIN
#define LMMP_MIN(l, o)
Definition lmmp.h:348
lmmp_add_nc_
mp_limb_t lmmp_add_nc_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
带进位的n位加法 [dst,n] = [numa,n] + [numb,n] + c
Definition add_n.c:9
lmmp_sub_nc_
mp_limb_t lmmp_sub_nc_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
带借位的n位减法 [dst,n] = [numa,n] - [numb,n] - c
Definition sub_n.c:9
PART_SIZE
#define PART_SIZE
Definition mparam.h:89
tp
#define tp

引用了 lmmp_add_nc_(), lmmp_copy, LMMP_MIN, lmmp_sub_nc_(), PART_SIZE , 以及 tp.

被这些函数引用 lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom43_(), lmmp_mul_toom52_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_(), lmmp_sqr_toom3_(), lmmp_toom_eval_dgr3_pm1_(), lmmp_toom_eval_dgr3_pm2_(), lmmp_toom_eval_pm1_() , 以及 lmmp_toom_eval_pm2_().

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

◆ lmmp_add_nc_()

mp_limb_t lmmp_add_nc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_limb_t  c 
)

带进位的n位加法 [dst,n] = [numa,n] + [numb,n] + c

参数
dst结果输出指针
numa第一个加数指针
numb第二个加数指针
nlimb长度
c初始进位值 [0|1]
警告
c=[0|1], n>0, eqsep(dst,[numa|numb])
返回
运算后的最终进位值 [0|1]

在文件 add_n.c 第 9 行定义.

9 {
10 mp_size_t i = 0;
11 mp_limb_t cy = c;
12
13 for (; i + 4 <= n; i += 4) {
14 mp_limb_t a0, b0, r0;
15 mp_limb_t a1, b1, r1;
16 mp_limb_t a2, b2, r2;
17 mp_limb_t a3, b3, r3;
18
19 a0 = numa[i + 0];
20 b0 = numb[i + 0];
21
22 a1 = numa[i + 1];
23 b1 = numb[i + 1];
24
25 a2 = numa[i + 2];
26 b2 = numb[i + 2];
27
28 a3 = numa[i + 3];
29 b3 = numb[i + 3];
30
31 r0 = a0 + cy;
32 cy = (r0 < cy);
33 r0 += b0;
34 cy += (r0 < b0);
35
36 r1 = a1 + cy;
37 cy = (r1 < cy);
38 r1 += b1;
39 cy += (r1 < b1);
40
41 r2 = a2 + cy;
42 cy = (r2 < cy);
43 r2 += b2;
44 cy += (r2 < b2);
45
46 r3 = a3 + cy;
47 cy = (r3 < cy);
48 r3 += b3;
49 cy += (r3 < b3);
50
51 dst[i + 0] = r0;
52 dst[i + 1] = r1;
53 dst[i + 2] = r2;
54 dst[i + 3] = r3;
55 }
56
57 for (; i < n; i++) {
58 mp_limb_t a, b, r;
59 a = numa[i];
60 b = numb[i];
61 r = a + cy;
62 cy = (r < cy);
63 r += b;
64 cy += (r < b);
65 dst[i] = r;
66 }
67
68 return cy;
69}

引用了 a0, a1, a2, a3, b0, b1, b2, b3, r0, r1, r2 , 以及 r3.

被这些函数引用 lmmp_add_n_sub_n_(), lmmp_fft_bfy_(), lmmp_ifft_bfy_() , 以及 lmmp_invappr_newton_().

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

◆ lmmp_addmul_1_()

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

参数
numa被加数指针（结果也存储在此）
numb乘数指针
nlimb长度
b乘数
警告
n>0, eqsep(numa,numb))
返回
运算后的进位limb值

被这些函数引用 lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_sqrt_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_addshl1_n_()

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)

参数
dst结果输出指针
numa被加数指针
numb加数指针（先左移1位）
nlimb长度
警告
n>0, eqsep(dst,[numa|numb])
返回
运算后的进位值 [0|1|2]

在文件 shl.c 第 56 行定义.

56 {
57 mp_size_t i, c = 0, mb = 0;
58
59 for (i = 0; i < n; i++) {
60 mp_limb_t a, b, r;
61 a = numa[i];
62 b = (numb[i] << 1) + mb;
63 mb = numb[i] >> (LIMB_BITS - 1);
64 r = a + c;
65 c = (r < c);
66 r += b;
67 c += (r < b);
68 dst[i] = r;
69 }
70 return c + mb;
71}
LIMB_BITS
#define LIMB_BITS
Definition lmmp.h:221

引用了 LIMB_BITS.

被这些函数引用 lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom44_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_sqr_toom3_(), lmmp_sqr_toom4_(), lmmp_sqrlo_dc_() , 以及 lmmp_sqrt_divide_().

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

◆ lmmp_bninv_()

void lmmp_bninv_ ( mp_ptr  dstq,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  ni 
)

精确逆元计算 [dstq,na+ni+2] = B^(2*(na+ni)) / ([numa,na] * B^ni)

参数
dstq输出商的缓冲区，长度至少为na+ni+2
numa输入被除数（长度na）
na被除数的 limb 长度
ni精度因子
警告
na>0, sep(dstq,numa), dstq!=NULL, numa[na-1]!=0
注解
也就是计算 B^(2*na+ni) div ([numa,na]

◆ lmmp_cmp_()

static int lmmp_cmp_ ( mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)
inlinestatic

大数比较函数（内联）

参数
numa第一个大数，长度为n
numb第二个大数，长度为n
n大数的单精度数(limb)长度
返回
1(numa>numb) / 0(numa==numb) / -1(numa<numb)
警告
n>0, numa!=NULL, numb!=NULL
注解
从最高位开始逐位比较，直到找到不同位

在文件 lmmpn.h 第 1004 行定义.

1004 {
1005 lmmp_param_assert(n > 0);
1006 lmmp_param_assert(numa != NULL);
1007 lmmp_param_assert(numb != NULL);
1008 mp_ssize_t i = n;
1009 mp_limb_t x, y;
1010 while (--i >= 0) {
1011 x = numa[i];
1012 y = numb[i];
1013 if (x != y)
1014 return (x > y ? 1 : -1);
1015 }
1016 return 0;
1017}
mp_ssize_t
int64_t mp_ssize_t
Definition lmmp.h:214
lmmp_param_assert
#define lmmp_param_assert(x)
Definition lmmp.h:398

引用了 lmmp_param_assert.

被这些函数引用 lmmp_add_signed_(), lmmp_div_(), lmmp_div_basecase_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_gcd_basecase_(), lmmp_gcd_lehmer_(), lmmp_invappr_newton_(), lmmp_lehmer_mul_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom43_(), lmmp_mul_toom52_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_(), lmmp_sqr_toom2_(), lmmp_sqr_toom3_(), lmmp_toom_eval_dgr3_pm1_(), lmmp_toom_eval_dgr3_pm2_(), lmmp_toom_eval_pm1_(), lmmp_toom_eval_pm2_() , 以及 try_div_().

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

◆ lmmp_div_()

void lmmp_div_ ( mp_ptr  dstq,
mp_ptr  dstr,
mp_srcptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)

大数除法和取模操作

注解
如果dstq不为NULL: [dstq,na-nb+1] = [numa,na] / [numb,nb] (商) 如果dstr不为NULL: [dstr,nb] = [numa,na] mod [numb,nb] (余数)
警告
0<nb<=na, numb[nb-1]!=0, sep(dstq,[numa|numb]), eqsep(dstr,[numa|numb])) 特殊情况: nb==1时, dstq>=numa-1 是允许的 nb==2时, dstq>=numa 是允许的
参数
dstq商结果输出指针（NULL表示不计算商）
dstr余数结果输出指针（NULL表示不计算余数）
numa被除数指针
na被除数的 limb 长度
numb除数指针
nb除数的 limb 长度

在文件 div.c 第 57 行定义.

57 {
58 if (nb == 1) {
59 mp_limb_t rem = lmmp_div_1_(dstq, numa, na, *numb);
60 if (dstr)
61 *dstr = rem;
62 } else if (nb == 2) {
63 mp_limb_t brem[2];
64 brem[0] = numb[0];
65 brem[1] = numb[1];
66 lmmp_div_2_(dstq, numa, na, brem);
67 if (dstr) {
68 dstr[0] = brem[0];
69 dstr[1] = brem[1];
70 }
71 } else {
72 int adjust = numa[na - 1] >= numb[nb - 1];
73 int cnt = lmmp_leading_zeros_(numb[nb - 1]);
74 mp_size_t nq = na + adjust - nb;
75 if (nq == 0) {
76 if (dstr && dstr != numa)
77 lmmp_copy(dstr, numa, nb);
78 if (dstq)
79 dstq[0] = 0;
80 return;
81 }
82 TEMP_DECL;
83
84 if (!dstq)
85 dstq = TALLOC_TYPE(na - nb + 1, mp_limb_t);
86 dstq[na - nb] = 0;
87
88 if (nq >= nb) {
89 mp_ptr numa2 = TALLOC_TYPE(na + 1, mp_limb_t);
90 mp_ptr numb2;
91 if (cnt) {
92 numa2[na] = lmmp_shl_(numa2, numa, na, cnt);
93 numb2 = TALLOC_TYPE(nb, mp_limb_t);
94 lmmp_shl_(numb2, numb, nb, cnt);
95 } else {
96 numa2[na] = 0;
97 lmmp_copy(numa2, numa, na);
98 numb2 = (mp_ptr)numb;
99 }
100
101 mp_limb_t inv21 = lmmp_inv_2_1_(numb2[nb - 1], numb2[nb - 2]);
102 na += adjust;
103
104 if (nb < DIV_DIVIDE_THRESHOLD)
105 lmmp_div_basecase_(dstq, numa2, na, numb2, nb, inv21);
106 else if (nb < DIV_MULINV_L_THRESHOLD || na < 2 * DIV_MULINV_N_THRESHOLD)
107 lmmp_div_divide_(dstq, numa2, na, numb2, nb, inv21);
108 else {
109 mp_limb_t ni = lmmp_div_inv_size_(nq, nb);
110 mp_ptr invappr = TALLOC_TYPE(ni, mp_limb_t);
111 lmmp_inv_prediv_(invappr, numb2, nb, ni);
112 lmmp_div_mulinv_(dstq, numa2, na, numb2, nb, invappr, ni);
113 }
114
115 if (dstr) {
116 if (cnt)
117 lmmp_shr_(dstr, numa2, nb, cnt);
118 else
119 lmmp_copy(dstr, numa2, nb);
120 }
121 } else {
122 // nq=na-nb+adj<nb
123 //-> na+adj>=2nq+1
124 mp_size_t ni = nb - nq;
125 mp_ptr numa2, numb2;
126 mp_ptr tp = TALLOC_TYPE(nb, mp_limb_t);
127 mp_limb_t cy;
128
129 numa2 = TALLOC_TYPE(nq * 2 + 1, mp_limb_t);
130 if (cnt) {
131 numb2 = TALLOC_TYPE(nq, mp_limb_t);
132 lmmp_shl_(numb2, numb + ni, nq, cnt);
133 numb2[0] |= numb[ni - 1] >> (LIMB_BITS - cnt);
134 cy = lmmp_shl_(numa2, numa + na - 2 * nq, 2 * nq, cnt);
135 if (adjust) {
136 numa2[2 * nq] = cy;
137 ++numa2; // numa2[0] is as significant as numa[ni=na-2nq+adjust]
138 } else
139 numa2[0] |= numa[na - 2 * nq - 1] >> (LIMB_BITS - cnt);
140 } else {
141 numb2 = (mp_ptr)numb + ni;
142 lmmp_copy(numa2, numa + na - 2 * nq, 2 * nq);
143 if (adjust) {
144 numa2[2 * nq] = 0;
145 ++numa2;
146 }
147 }
148
149 // now: 0<=numa2<B^2nq, B^nq/2<=numb2<B^nq, and 0<=numa2/numb2<B^nq
150 // ignored bits could be seen as fraction part of numa and numb
151 // we can prove: Q<=Qh<=Q+2
152 // where Q=floor(numa/numb) is the real quotient
153 // Qh=floor(floor(numa)/floor(numb)) as below
154
155 if (nq == 1) {
156 lmmp_div_1_s_(dstq, numa2, 2, *numb2);
157 } else if (nq == 2) {
158 lmmp_div_2_s_(dstq, numa2, 4, numb2);
159 } else {
160 mp_limb_t inv21 = lmmp_inv_2_1_(numb2[nq - 1], numb2[nq - 2]);
161
162 if (nq < DIV_DIVIDE_THRESHOLD)
163 lmmp_div_basecase_(dstq, numa2, 2 * nq, numb2, nq, inv21);
164 else if (nq < DIV_MULINV_N_THRESHOLD)
165 lmmp_div_divide_(dstq, numa2, 2 * nq, numb2, nq, inv21);
166 else {
167 mp_limb_t ni = lmmp_div_inv_size_(nq, nq);
168 mp_ptr invappr = tp;
169 lmmp_inv_prediv_(invappr, numb2, nq, ni);
170 lmmp_div_mulinv_(dstq, numa2, 2 * nq, numb2, nq, invappr, ni);
171 }
172 }
173 /*
174 true remainder = partial remainder - quotient * ignored divisor limbs
175
176 Multiply the first ignored divisor limb by the most significant
177 quotient limb. If that product is > the partial remainder's
178 most significant limb, we know the quotient is too large. This
179 test quickly catches most cases where the quotient is too large;
180 it catches all cases where the quotient is 2 too large.*/
181
182 mp_limb_t x;
183 if (cnt) {
184 mp_limb_t dl;
185 if (ni < 2)
186 dl = 0;
187 else
188 dl = numb[ni - 2];
189 x = (numb[ni - 1] << cnt) | (dl >> (LIMB_BITS - cnt));
190 } else
191 x = numb[ni - 1];
192 mp_limb_t h = (x >> LIMB_BITS / 2) * (dstq[nq - 1] >> LIMB_BITS / 2);
193 mp_limb_t rnb = 0; // remainder[nb]
194 mp_size_t nr = nq; // remainder=rnb:[numa2,nr]:[...,ni]
195
196 if (h > numa2[nq - 1]) {
197 lmmp_dec(dstq);
198 rnb = lmmp_add_n_(numa2, numa2, numb2, nq);
199 }
200
201 // if cnt, recover the shift of partial remainder
202 // and remove the effect of the partial-ignored numa[ni-1] and numb[ni-1]
203 if (cnt) {
204 numa2[nq] = rnb;
205 ++nr;
206 --ni;
207 lmmp_shl_(numa2, numa2, nr, LIMB_BITS - cnt);
208 numa2[0] |= numa[ni] & (LIMB_MAX >> cnt);
209 cy = lmmp_submul_1_(numa2, dstq, nq, numb[ni] & (LIMB_MAX >> cnt));
210 rnb = -(numa2[nq] < cy);
211 numa2[nq] -= cy;
212 }
213
214 if (ni == 0) {
215 if (dstr) {
216 if (rnb)
217 lmmp_add_n_(dstr, numa2, numb, nr);
218 else
219 lmmp_copy(dstr, numa2, nr);
220 }
221 } else {
222 tp[nb - 1] = 0;
223 if (ni < nq)
224 lmmp_mul_(tp, dstq, nq, numb, ni);
225 else
226 lmmp_mul_(tp, numb, ni, dstq, nq);
227
228 if (dstr) {
229 mp_ptr remptr = dstr == numb ? tp : dstr;
230 cy = lmmp_sub_n_(remptr, numa, tp, ni);
231 rnb -= lmmp_sub_nc_(remptr + ni, numa2, tp + ni, nr, cy);
232 if (rnb)
233 lmmp_add_n_(dstr, remptr, numb, nb);
234 else if (dstr != remptr)
235 lmmp_copy(dstr, remptr, nb);
236 } else {
237 int hcmp = lmmp_cmp_(numa2, tp + ni, nr);
238 if (hcmp < 0)
239 --rnb;
240 else if (hcmp == 0)
241 rnb -= (lmmp_cmp_(numa, tp, ni) < 0);
242 }
243 }
244
245 if (rnb)
246 lmmp_dec(dstq);
247 }
248
249 TEMP_FREE;
250 }
251}
LIMB_MAX
#define LIMB_MAX
Definition lmmp.h:224
lmmp_div_inv_size_
static mp_size_t lmmp_div_inv_size_(mp_size_t nq, mp_size_t nb)
计算预计算逆元的尺寸
Definition lmmpn.h:812
lmmp_div_1_s_
mp_limb_t lmmp_div_1_s_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_limb_t x)
单精度数除法（除数为1个limb）
lmmp_div_1_
mp_limb_t lmmp_div_1_(mp_ptr dstq, mp_srcptr numa, mp_size_t na, mp_limb_t x)
单精度数除法
Definition div.c:66
lmmp_leading_zeros_
int lmmp_leading_zeros_(mp_limb_t x)
计算一个单精度数(limb)中前导零的个数
Definition tiny.c:35
lmmp_cmp_
static int lmmp_cmp_(mp_srcptr numa, mp_srcptr numb, mp_size_t n)
大数比较函数（内联）
Definition lmmpn.h:1004
lmmp_dec
#define lmmp_dec(p)
大数减1宏（预期无借位）
Definition lmmpn.h:973
lmmp_inv_prediv_
void lmmp_inv_prediv_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t ni)
除法前的逆元预计算，[dst,ni] = invappr( (ni+1 MSLs of numa) + 1 ) / B
Definition div_mulinv.c:11
lmmp_div_2_
void lmmp_div_2_(mp_ptr dstq, mp_srcptr numa, mp_size_t na, mp_ptr numb)
双精度数除法 (除数为2个limb)
Definition div.c:223
lmmp_shr_
mp_limb_t lmmp_shr_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shr)
大数右移操作 [dst,na] = [numa,na] >> shr，dst的高shr位填充0
Definition shr.c:9
lmmp_mul_
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]
lmmp_shl_
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:9
lmmp_div_mulinv_
mp_limb_t lmmp_div_mulinv_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_srcptr invappr, mp_size_t ni)
乘法逆元除法
Definition div_mulinv.c:36
lmmp_div_2_s_
mp_limb_t lmmp_div_2_s_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb)
双精度数除法（除数为2个limb）
lmmp_submul_1_
mp_limb_t lmmp_submul_1_(mp_ptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t b)
大数乘以单limb并累减操作 [numa,n] -= [numb,n] * b
lmmp_sub_n_
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:70
lmmp_inv_2_1_
mp_limb_t lmmp_inv_2_1_(mp_limb_t xh, mp_limb_t xl)
2-1阶逆元计算 (inv21)
Definition inv.c:10
lmmp_div_basecase_
mp_limb_t lmmp_div_basecase_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_limb_t inv21)
基础除法运算
Definition div_basecase.c:10
lmmp_div_divide_
mp_limb_t lmmp_div_divide_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb, mp_limb_t inv21)
分治除法运算
Definition div_divide.c:53
DIV_DIVIDE_THRESHOLD
#define DIV_DIVIDE_THRESHOLD
Definition mparam.h:26
DIV_MULINV_N_THRESHOLD
#define DIV_MULINV_N_THRESHOLD
Definition mparam.h:30
DIV_MULINV_L_THRESHOLD
#define DIV_MULINV_L_THRESHOLD
Definition mparam.h:28
TEMP_DECL
#define TEMP_DECL
Definition tmp_alloc.h:72
TEMP_FREE
#define TEMP_FREE
Definition tmp_alloc.h:93
TALLOC_TYPE
#define TALLOC_TYPE(n, type)
Definition tmp_alloc.h:91

引用了 DIV_DIVIDE_THRESHOLD, DIV_MULINV_L_THRESHOLD, DIV_MULINV_N_THRESHOLD, LIMB_BITS, LIMB_MAX, lmmp_add_n_(), lmmp_cmp_(), lmmp_copy, lmmp_dec, lmmp_div_1_(), lmmp_div_1_s_(), lmmp_div_2_(), lmmp_div_2_s_(), lmmp_div_basecase_(), lmmp_div_divide_(), lmmp_div_inv_size_(), lmmp_div_mulinv_(), lmmp_inv_2_1_(), lmmp_inv_prediv_(), lmmp_leading_zeros_(), lmmp_mul_(), lmmp_shl_(), lmmp_shr_(), lmmp_sub_n_(), lmmp_sub_nc_(), lmmp_submul_1_(), TALLOC_TYPE, TEMP_DECL, TEMP_FREE , 以及 tp.

被这些函数引用 lmmp_bninv_(), lmmp_gcd_basecase_(), lmmp_gcd_lehmer_(), lmmp_trialdiv_() , 以及 try_div_().

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

◆ lmmp_div_1_()

mp_limb_t lmmp_div_1_ ( mp_ptr  dstq,
mp_srcptr  numa,
mp_size_t  na,
mp_limb_t  x 
)

单精度数除法

参数
dstq输出商的缓冲区（可为NULL，此时仅计算余数）
numa输入被除数，长度为na
na被除数的 limb 长度
x除数（单个 limb ）
返回
除法余数（单个 limb ）
警告
na>0, x!=0, eqsep(dstq,numa), dstq>=numa-1 是可以接受的
注解
if (dstq!=NULL) [dstq,na] = [numa,na] div x

在文件 div.c 第 66 行定义.

66 {
67 mp_limb_t ah, al;
68 if (na == 1) {
69 ah = numa[0];
70 if (dstq)
71 dstq[0] = ah / x;
72 return ah % x;
73 }
74 if (dstq) {
75 mp_limb_t t = numa[na - 2], q = 0, r = 0;
76 const int shift = lmmp_leading_zeros_(x);
77 if (shift > 0) {
78 /*
79 ah al
80 X|XXXtttX|XXXmmmX|XXXnnnX|XXX----|
81 |000XXXX|tttXXXX|mmmXXXX|nnnXXXX|
82 t numa[na]
83
84 ah al
85 X|XXXtttX|XXXmmmX|XXXnnnX|XXX----|
86 |000XXXX|tttXXXX|mmmXXXX|nnnXXXX|
87 t
88 numa[na]
89 */
90 const int rshift = LIMB_BITS - shift;
91 ah = numa[na - 1] >> rshift;
92 t = numa[na - 2];
93 al = (numa[na - 1] << shift) | (t >> rshift);
94 x <<= shift;
95 const mp_limb_t inv = lmmp_inv_1_(x);
96 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
97 dstq[na - 1] = q;
98 na -= 2;
99 while (na-- > 0) {
100 ah = r;
101 al = t << shift;
102 t = numa[na];
103 al |= t >> rshift;
104 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
105 dstq[na + 1] = q;
106 }
107 ah = r;
108 al = t << shift;
109 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
110 dstq[0] = q;
111 return r >> shift;
112 } else {
113 /*
114 ah al
115 |000XXXX|tttXXXX|mmmXXXX|nnnXXXX|
116 t numa[na]
117 */
118 ah = 0;
119 t = numa[na - 2];
120 al = numa[na - 1];
121 const mp_limb_t inv = lmmp_inv_1_(x);
122 q = al / x;
123 r = al % x;
124 dstq[na - 1] = q;
125 na -= 2;
126 while (na-- > 0) {
127 ah = r;
128 al = t;
129 t = numa[na];
130 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
131 dstq[na + 1] = q;
132 }
133 ah = r;
134 al = t;
135 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
136 dstq[0] = q;
137 return r;
138 }
139 } else {
140 return lmmp_mod_1_(numa, na, x);
141 }
142}
lmmp_mod_1_
mp_limb_t lmmp_mod_1_(mp_srcptr numa, mp_size_t na, mp_limb_t x)
单精度数取余
Definition div.c:20
lmmp_inv_1_
mp_limb_t lmmp_inv_1_(mp_limb_t x)
1阶逆元计算 (inv1)
Definition inv.c:107
_udiv_qrnnd_preinv
#define _udiv_qrnnd_preinv(q, r, nh, nl, d, di)
Definition longlong.h:424

引用了 _udiv_qrnnd_preinv, LIMB_BITS, lmmp_inv_1_(), lmmp_leading_zeros_() , 以及 lmmp_mod_1_().

被这些函数引用 lmmp_bninv_(), lmmp_div_(), lmmp_odd_nCr_uint_(), lmmp_odd_nCr_ushort_() , 以及 lmmp_to_str_basecase_().

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

◆ lmmp_div_1_s_()

mp_limb_t lmmp_div_1_s_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_limb_t  x 
)

单精度数除法（除数为1个limb）

参数
dstq输出商的缓冲区，长度至少为na-1
numa输入被除数（长度na），运算后存储余数（长度1）
na被除数的 limb 长度
x除数（单个 limb ）
返回
商的最高位（qh）
警告
na>1, MSB(x)=1, sep(dstq,numa)
注解
qh:[dstq,na-1]=[numa,na] div x, [numa,1]=[numa,na] mod x, return qh

被这些函数引用 lmmp_div_(), lmmp_div_s_(), mont63_R2() , 以及 mont64_R2().

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

◆ lmmp_div_2_()

void lmmp_div_2_ ( mp_ptr  dstq,
mp_srcptr  numa,
mp_size_t  na,
mp_ptr  numb 
)

双精度数除法 (除数为2个limb)

参数
dstq输出商的缓冲区，长度至少为na-1
numa输入被除数（长度na）
na被除数的 limb 长度
numb输入除数（长度2）[numb,2]=[numa,na] mod [numb,2]
警告
na>=2, numb[1]!=0, eqsep(dstq,numa), dstq>=numa 是可以接受的
注解
if (dstq!=NULL) [dstq,na-1]=[numa,na] div [numb,2]

在文件 div.c 第 223 行定义.

223 {
224 mp_limb_t q, r1, r0, a2, a1, a0, b1, b0;
225 b1 = numb[1];
226 b0 = numb[0];
227 if (na == 2) {
228 int shift = lmmp_leading_zeros_(b1);
229 if (shift > 0) {
230 const int rshift = LIMB_BITS - shift;
231 b1 = (b1 << shift) | (b0 >> rshift);
232 b0 <<= shift;
233 a2 = numa[1] >> rshift;
234 a1 = (numa[1] << shift) | (numa[0] >> rshift);
235 a0 = (numa[0] << shift);
236 mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
237 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
238 if (dstq)
239 dstq[0] = q;
240 numb[0] = (r0 >> shift) | (r1 << rshift);
241 numb[1] = r1 >> shift;
242 return;
243 } else {
244 if (_u128cmp(numa, numb)) {
245 numb[0] = numa[0];
246 numb[1] = numa[1];
247 if (dstq)
248 dstq[0] = 0;
249 return;
250 } else {
251 _u128sub(numb, numa, numb);
252 if (dstq)
253 dstq[0] = 1;
254 return;
255 }
256 }
257 }
258 if (dstq) {
259 int shift = lmmp_leading_zeros_(b1);
260 if (shift > 0) {
261 /*
262 a2 a1 a0
263 X|XXXtttX|XXXmmmX|XXXnnnX|XXX----|
264 |000XXXX|tttXXXX|mmmXXXX|nnnXXXX|
265 numa[na]
266 */
267 const int rshift = LIMB_BITS - shift;
268 b1 = (b1 << shift) | (b0 >> rshift);
269 b0 <<= shift;
270 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
271 a2 = numa[na - 1] >> rshift;
272 a1 = (numa[na - 1] << shift) | (numa[na - 2] >> rshift);
273 a0 = (numa[na - 2] << shift) | (numa[na - 3] >> rshift);
274 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
275 dstq[na - 2] = q;
276 na -= 2;
277 while (na-- > 1) {
278 a2 = r1;
279 a1 = r0;
280 a0 = (numa[na] << shift) | (numa[na - 1] >> rshift);
281 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
282 dstq[na] = q;
283 }
284
285 a2 = r1;
286 a1 = r0;
287 a0 = (numa[na] << shift);
288 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
289 dstq[0] = q;
290 numb[0] = (r0 >> shift) | (r1 << rshift);
291 numb[1] = r1 >> shift;
292 return;
293 } else {
294 /*
295 a2 a1 a0
296 |000XXXX|tttXXXX|mmmXXXX|nnnXXXX|
297 numa[na]
298 */
299 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
300 a2 = 0;
301 a1 = numa[na - 1];
302 a0 = numa[na - 2];
303 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
304 dstq[na - 2] = q;
305 na -= 2;
306 while (na-- > 1) {
307 a2 = r1;
308 a1 = r0;
309 a0 = numa[na];
310 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
311 dstq[na] = q;
312 }
313 a2 = r1;
314 a1 = r0;
315 a0 = numa[na];
316 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
317 dstq[0] = q;
318 numb[0] = r0;
319 numb[1] = r1;
320 return;
321 }
322 } else {
323 int shift = lmmp_leading_zeros_(b1);
324 if (shift > 0) {
325 const int rshift = LIMB_BITS - shift;
326 b1 = (b1 << shift) | (b0 >> rshift);
327 b0 <<= shift;
328 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
329 a2 = numa[na - 1] >> rshift;
330 a1 = (numa[na - 1] << shift) | (numa[na - 2] >> rshift);
331 a0 = (numa[na - 2] << shift) | (numa[na - 3] >> rshift);
332 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
333 na -= 2;
334 while (na-- > 1) {
335 a2 = r1;
336 a1 = r0;
337 a0 = (numa[na] << shift) | (numa[na - 1] >> rshift);
338 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
339 }
340
341 a2 = r1;
342 a1 = r0;
343 a0 = (numa[na] << shift);
344 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
345 numb[0] = (r0 >> shift) | (r1 << rshift);
346 numb[1] = r1 >> shift;
347 return;
348 } else {
349 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
350 a2 = 0;
351 a1 = numa[na - 1];
352 a0 = numa[na - 2];
353 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
354 na -= 2;
355 while (na-- > 1) {
356 a2 = r1;
357 a1 = r0;
358 a0 = numa[na];
359 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
360 }
361 a2 = r1;
362 a1 = r0;
363 a0 = numa[na];
364 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
365 numb[0] = r0;
366 numb[1] = r1;
367 return;
368 }
369 }
370}
_u128sub
#define _u128sub(r, x, y)
Definition longlong.h:282
_u128cmp
#define _u128cmp(x, y)
Definition longlong.h:280
_udiv_qr_3by2
#define _udiv_qr_3by2(q, r1, r0, n2, n1, n0, d1, d0, dinv)
Definition longlong.h:441

引用了 _u128cmp, _u128sub, _udiv_qr_3by2, a0, a1, a2, b0, b1, LIMB_BITS, lmmp_inv_2_1_(), lmmp_leading_zeros_(), r0 , 以及 r1.

被这些函数引用 lmmp_bninv_() , 以及 lmmp_div_().

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

◆ lmmp_div_2_s_()

mp_limb_t lmmp_div_2_s_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_srcptr  numb 
)

双精度数除法（除数为2个limb）

参数
dstq输出商的缓冲区，长度至少为na-2
numa输入被除数（长度na），运算后存储余数（长度2）
na被除数的 limb 长度
numb输入除数，长度为2
返回
商的最高位（qh）
警告
na>2, MSB(numb)=1, sep(dstq,numa,numb)
注解
qh:[dstq,na-2]=[numa,na] div [numb,2], [numa,2]=[numa,na] mod [numb,2], return qh

被这些函数引用 lmmp_div_(), lmmp_div_s_() , 以及 lmmp_inv_basecase_().

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

◆ lmmp_div_3_2_()

mp_limb_t lmmp_div_3_2_ ( mp_ptr  numa,
mp_srcptr  numb,
mp_limb_t  inv21 
)

3/2位除法运算 [numa,2]=[numa,3] mod [numb,2]

参数
numa输入被除数（长度3），运算后存储余数（长度2）
numb输入除数（长度2）
inv21除数的2-1阶逆元（提前计算好的inv21([numb,2])）
返回
商值（单精度数）
警告
[numa,3]<[numb,2]*B, MSB(numb)=1, inv21=inv21([numb,2]), eqsep(numa,numb)

被这些函数引用 lmmp_div_basecase_().

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

◆ lmmp_div_basecase_()

mp_limb_t lmmp_div_basecase_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb,
mp_limb_t  inv21 
)

基础除法运算

参数
dstq输出商的缓冲区，长度至少为na-nb
numa输入被除数（长度na），运算后存储余数（长度nb）
na被除数的单精度数(limb)长度
numb输入除数，长度为nb
nb除数的单精度数(limb)长度
inv21除数的2-1阶逆元（inv21([numb+nb-2,2])）
返回
商的最高位（qh）
警告
na>=nb>=3, MSB(numb)=1, inv21=inv21([numb+nb-2,2]), sep(dstq,numa,numb)
注解
qh:[dstq,na-nb]=[numa,na] div [numb,nb], [numa,na-nb]=[numa,na] mod [numb,nb], return qh

在文件 div_basecase.c 第 10 行定义.

10 {
11 lmmp_param_assert(na >= nb);
12 lmmp_param_assert(nb >= 3);
13 lmmp_param_assert(numb[nb - 1] >= LIMB_B_2);
14 mp_size_t nq = na - nb;
15
16 numa += na;
17
18 mp_limb_t qh = lmmp_cmp_(numa - nb, numb, nb) >= 0;
19 if (qh) {
20 lmmp_sub_n_(numa - nb, numa - nb, numb, nb);
21 }
22
23 nb -= 2;
24 numa -= 2;
25
26 mp_limb_t d1 = numb[nb + 1], d0 = numb[nb + 0];
27
28 while (nq) {
29 mp_limb_t q;
30 --numa;
31 if (numa[2] == d1 && numa[1] == d0) {
32 q = LIMB_MAX;
33 lmmp_submul_1_(numa - nb, numb, nb + 2, q);
34 } else {
35 mp_limb_t cy, cy1;
36 q = lmmp_div_3_2_(numa, numb + nb, inv21);
37 cy = lmmp_submul_1_(numa - nb, numb, nb, q);
38 cy1 = numa[0] < cy;
39 numa[0] -= cy;
40 cy = numa[1] < cy1;
41 numa[1] -= cy1;
42 if (cy) {
43 lmmp_add_n_(numa - nb, numa - nb, numb, nb + 2);
44 --q;
45 }
46 }
47 dstq[--nq] = q;
48 }
49 return qh;
50}
lmmp_div_3_2_
mp_limb_t lmmp_div_3_2_(mp_ptr numa, mp_srcptr numb, mp_limb_t inv21)
3/2位除法运算 [numa,2]=[numa,3] mod [numb,2]
LIMB_B_2
#define LIMB_B_2
Definition mparam.h:160

引用了 LIMB_B_2, LIMB_MAX, lmmp_add_n_(), lmmp_cmp_(), lmmp_div_3_2_(), lmmp_param_assert, lmmp_sub_n_() , 以及 lmmp_submul_1_().

被这些函数引用 lmmp_bninv_appr_newton_(), lmmp_div_(), lmmp_div_divide_n_(), lmmp_div_s_() , 以及 lmmp_inv_basecase_().

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

◆ lmmp_div_divide_()

mp_limb_t lmmp_div_divide_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb,
mp_limb_t  inv21 
)

分治除法运算

参数
dstq输出商的缓冲区，长度至少为na-nb
numa输入被除数（长度na），运算后存储余数（长度nb）
na被除数的单精度数(limb)长度
numb输入除数，长度为nb
nb除数的单精度数(limb)长度
inv21除数的2-1阶逆元（inv21([numb+nb-2,2])）
返回
商的最高位（qh）
警告
na>=2*nb, nb>=6, MSB(numb)=1, inv21=inv21([numb+nb-2,2]), sep(dstq,numa,numb)
注解
qh:[dstq,na-nb]=[numa,na] div [numb,nb], [numa,na-nb]=[numa,na] mod [numb,nb], return qh

在文件 div_divide.c 第 53 行定义.

53 {
54 lmmp_param_assert(na >= 2 * nb);
55 lmmp_param_assert(nb >= 6);
56 lmmp_param_assert(numb[nb - 1] >= LIMB_B_2);
57 mp_size_t nq = na - nb;
58
59 dstq += nq;
60 numa += nq;
61
62 do {
63 nq -= nb;
64 } while (nq >= nb);
65
66 dstq -= nq;
67 numa -= nq;
68
69 /* Perform the typically smaller block first. */
70 mp_limb_t qh = lmmp_div_s_(dstq, numa, nb + nq, numb, nb);
71
72 TEMP_DECL;
73 mp_ptr tp = TALLOC_TYPE(nb, mp_limb_t);
74 nq = na - nb - nq;
75
76 do {
77 dstq -= nb;
78 numa -= nb;
79 lmmp_div_divide_n_(dstq, numa, numb, nb, inv21, tp);
80 nq -= nb;
81 } while (nq > 0);
82
83 TEMP_FREE;
84 return qh;
85}
lmmp_div_divide_n_
static mp_limb_t lmmp_div_divide_n_(mp_ptr dstq, mp_ptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t inv21, mp_ptr tp)
Definition div_divide.c:13
lmmp_div_s_
mp_limb_t lmmp_div_s_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
除法运算
Definition div.c:11

引用了 LIMB_B_2, lmmp_div_divide_n_(), lmmp_div_s_(), lmmp_param_assert, TALLOC_TYPE, TEMP_DECL, TEMP_FREE , 以及 tp.

被这些函数引用 lmmp_div_(), lmmp_div_s_() , 以及 lmmp_inv_basecase_().

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

◆ lmmp_div_inv_size_()

static mp_size_t lmmp_div_inv_size_ ( mp_size_t  nq,
mp_size_t  nb 
)
inlinestatic

计算预计算逆元的尺寸

参数
nq商的 limb 长度
nb除数的 limb 长度
返回
计算需要预计算逆元尺寸ni（ni<=nb）
注解
用于已归一化除法([nq+nb]/[nb]=[nq])的逆元 ni 尺寸

在文件 lmmpn.h 第 812 行定义.

812 {
813 mp_size_t ni, b;
814 if (nq > nb) {
815 b = (nq - 1) / nb + 1; // ceil(nq/nb), number of blocks
816 ni = (nq - 1) / b + 1; // ceil(nq/b)
817 } else if (3 * nq > nb) {
818 ni = (nq - 1) / 2 + 1; // b=2
819 } else {
820 ni = (nq - 1) / 1 + 1; // b=1
821 }
822 return ni;
823}

被这些函数引用 lmmp_div_(), lmmp_div_s_() , 以及 lmmp_to_str_().

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

◆ lmmp_div_mulinv_()

mp_limb_t lmmp_div_mulinv_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb,
mp_srcptr  invappr,
mp_size_t  ni 
)

乘法逆元除法

参数
dstq输出商的缓冲区，长度至少为na-nb
numa输入被除数（长度na），运算后存储余数（长度nb）
na被除数的 limb 长度
numb输入除数，长度为nb
nb除数的 limb 长度
invappr预计算的近似逆元，长度为ni
ni预计算逆元的 limb 长度
返回
商的最高位（qh）
警告
na>=nb>=ni>0, MSB(numb)=1, [invappr,ni]=inv_prediv([numb,nb]), sep(dstq,numa,numb,invappr))
注解
qh:[dstq,na-1]=[numa,na] div x, [numa,1]=[numa,na] mod x, return qh

在文件 div_mulinv.c 第 36 行定义.

44 {
45 lmmp_param_assert(na >= nb && nb >= ni);
46 lmmp_param_assert(ni > 0);
47 lmmp_param_assert(numb[nb - 1] >= LIMB_B_2);
48 mp_size_t nq = na - nb, ntp = LMMP_MIN(ni, nq) + nb;
49 mp_limb_t qh;
50 TEMP_DECL;
51 mp_ptr tp = TALLOC_TYPE(ntp, mp_limb_t);
52
53 numa += nq;
54 dstq += nq;
55
56 qh = lmmp_cmp_(numa, numb, nb) >= 0;
57 if (qh) {
58 lmmp_sub_n_(numa, numa, numb, nb);
59 }
60 while (nq) {
61 if (nq < ni) {
62 invappr += ni - nq;
63 ni = nq;
64 }
65 numa -= ni;
66 dstq -= ni;
67 nq -= ni;
68
69 lmmp_mul_n_(tp, numa + nb, invappr, ni);
70 lmmp_assert(lmmp_add_n_(dstq, tp + ni, numa + nb, ni) == 0);
71
72 mp_size_t mn, wn;
73 mp_limb_t cy;
74
75 if (nb < DIV_MULINV_MODM_THRESHOLD || (mn = lmmp_fft_next_size_(nb + 1)) >= nb + ni) {
76 lmmp_mul_(tp, numb, nb, dstq, ni); // nb+ni limbs, high 'ni' cancels
77 } else {
78 // 0<wn<ni<=nb<mn<nb+ni
79 wn = nb + ni - mn;
80
81 // x=b*q
82 // tp=x mod 2^mn-1
83 lmmp_mul_mersenne_(tp, mn, numb, nb, dstq, ni);
84
85 // tp-=ah:0 mod B^mn-1, if result=0, represent it as B^mn-1
86 cy = lmmp_sub_nc_(tp, tp, numa + mn, wn, 1);
87 if (cy)
88 cy = lmmp_sub_1_(tp + wn, tp + wn, mn - wn, 1);
89 if (!cy)
90 lmmp_inc(tp);
91
92 // if al<<tp,
93 if (lmmp_cmp_(numa + nb, tp + nb, mn - nb) < 0) {
94 // maybe ah=xh+1 and al<<xl,
95 // so we subtracted 1 too much when tp-=ah,
96 // now tp=xl-1 mod B^mn-1, and 0<=al<<xl-1<B^mn-1, so tp=xl-1
97 // or ah=xh and al>=xl,
98 // tp=xl mod B^mn-1, the only possibility is we represented xl=0 as tp=B^mn-1
99 // whatever, just inc and then tp=xl
100 tp[mn] = 0; // set a limit
101 lmmp_inc(tp);
102 }
103 }
104
105 mp_limb_t r = numa[nb] - tp[nb];
106 cy = lmmp_sub_n_(numa, numa, tp, nb);
107
108 while ((r -= cy) || lmmp_cmp_(numa, numb, nb) >= 0) {
109 lmmp_inc(dstq);
110 cy = lmmp_sub_n_(numa, numa, numb, nb);
111 }
112 }
113 TEMP_FREE;
114 return qh;
115}
lmmp_assert
#define lmmp_assert(x)
Definition lmmp.h:370
lmmp_mul_mersenne_
void lmmp_mul_mersenne_(mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
梅森数模乘法 [dst,rn] = [numa,na]*[numb,nb] mod B^rn-1
Definition mul_fft.c:752
lmmp_inc
#define lmmp_inc(p)
大数加1宏（预期无进位）
Definition lmmpn.h:946
lmmp_mul_n_
void lmmp_mul_n_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n)
等长大数乘法操作 [dst,2*n] = [numa,n] * [numb,n]
Definition mul.c:99
lmmp_fft_next_size_
mp_size_t lmmp_fft_next_size_(mp_size_t n)
计算满足 >=n 的最小费马/梅森乘法可行尺寸
Definition mul_fft.c:84
lmmp_sub_1_
static mp_limb_t lmmp_sub_1_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
大数减单精度数静态内联函数 [dst,na]=[numa,na]-x
Definition lmmpn.h:1122
DIV_MULINV_MODM_THRESHOLD
#define DIV_MULINV_MODM_THRESHOLD
Definition mparam.h:38

引用了 DIV_MULINV_MODM_THRESHOLD, LIMB_B_2, lmmp_add_n_(), lmmp_assert, lmmp_cmp_(), lmmp_fft_next_size_(), lmmp_inc, LMMP_MIN, lmmp_mul_(), lmmp_mul_mersenne_(), lmmp_mul_n_(), lmmp_param_assert, lmmp_sub_1_(), lmmp_sub_n_(), lmmp_sub_nc_(), TALLOC_TYPE, TEMP_DECL, TEMP_FREE , 以及 tp.

被这些函数引用 lmmp_div_(), lmmp_div_s_() , 以及 lmmp_to_str_divide_().

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

◆ lmmp_div_s_()

mp_limb_t lmmp_div_s_ ( mp_ptr  dstq,
mp_ptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)

除法运算

参数
dstq输出商的缓冲区，长度至少为na-nb
numa输入被除数（长度na），运算后存储余数（长度nb）
na被除数的 limb 长度
numb输入除数，长度为nb
nb除数的 limb 长度
返回
商的最高位（qh）
警告
na>=nb>0, MSB(numb)=1, sep(dstq,numa,numb)
注解
qh:[dstq,na-nb]=[numa,na] div [numb,nb], [numa,nb]=[numa,na] mod [numb,nb], return qh

在文件 div.c 第 11 行定义.

11 {
12 TEMP_DECL;
13 mp_limb_t nq = na - nb;
14 mp_limb_t qh;
15 if (nq == 0) {
16 qh = lmmp_cmp_(numa, numb, nb) >= 0;
17 if (qh)
18 lmmp_sub_n_(numa, numa, numb, nb);
19 } else if (nb == 1) {
20 qh = lmmp_div_1_s_(dstq, numa, na, *numb);
21 } else if (nb == 2) {
22 qh = lmmp_div_2_s_(dstq, numa, na, numb);
23 } else if (nq < nb) {
24 qh = lmmp_div_s_(dstq, numa + na - 2 * nq, 2 * nq, numb + nb - nq, nq);
25
26 mp_ptr tp = TALLOC_TYPE(nb, mp_limb_t);
27 if (nq > nb - nq)
28 lmmp_mul_(tp, dstq, nq, numb, nb - nq);
29 else
30 lmmp_mul_(tp, numb, nb - nq, dstq, nq);
31
32 mp_limb_t cy = lmmp_sub_n_(numa, numa, tp, nb);
33 if (qh)
34 cy += lmmp_sub_n_(numa + nq, numa + nq, numb, nb - nq);
35
36 while (cy) {
37 qh -= lmmp_sub_1_(dstq, dstq, nq, 1);
38 cy -= lmmp_add_n_(numa, numa, numb, nb);
39 }
40 } else {
41 mp_limb_t inv21 = lmmp_inv_2_1_(numb[nb - 1], numb[nb - 2]);
42 if (nb < DIV_DIVIDE_THRESHOLD)
43 qh = lmmp_div_basecase_(dstq, numa, na, numb, nb, inv21);
44 else if (nb < DIV_MULINV_L_THRESHOLD || na < 2 * DIV_MULINV_N_THRESHOLD)
45 qh = lmmp_div_divide_(dstq, numa, na, numb, nb, inv21);
46 else {
47 mp_limb_t ni = lmmp_div_inv_size_(nq, nb);
48 mp_ptr invappr = TALLOC_TYPE(ni, mp_limb_t);
49 lmmp_inv_prediv_(invappr, numb, nb, ni);
50 qh = lmmp_div_mulinv_(dstq, numa, na, numb, nb, invappr, ni);
51 }
52 }
53 TEMP_FREE;
54 return qh;
55}
lmmp_div_s_
mp_limb_t lmmp_div_s_(mp_ptr dstq, mp_ptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
除法运算
Definition div.c:11

引用了 DIV_DIVIDE_THRESHOLD, DIV_MULINV_L_THRESHOLD, DIV_MULINV_N_THRESHOLD, lmmp_add_n_(), lmmp_cmp_(), lmmp_div_1_s_(), lmmp_div_2_s_(), lmmp_div_basecase_(), lmmp_div_divide_(), lmmp_div_inv_size_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_inv_2_1_(), lmmp_inv_prediv_(), lmmp_mul_(), lmmp_sub_1_(), lmmp_sub_n_(), TALLOC_TYPE, TEMP_DECL, TEMP_FREE , 以及 tp.

被这些函数引用 lmmp_div_divide_(), lmmp_div_s_(), lmmp_invsqrt_newton_(), lmmp_sqrt_divide_() , 以及 lmmp_to_str_divide_().

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

◆ lmmp_endian()

static bool lmmp_endian ( void  )
inlinestatic

运行时判断端序

返回
true 表示小端序，false 表示大端序

在文件 lmmpn.h 第 69 行定义.

69 {
70 int num = 1;
71 return (*(char*)&num) == 0;
72}

◆ lmmp_extract_bits_()

mp_bitcnt_t lmmp_extract_bits_ ( mp_srcptr  num,
mp_size_t  n,
mp_limb_t *  ext,
int  bits 
)

提取高位指定位数，并返回低位bits位数

参数
num待提取的大数指针
nnum的 limb 长度
bits待提取的位数(1-64)
ext提取结果输出指针
警告
n>0, 1<=bits<=64, ext!=NULL
注解
如果bits大于num的实际位数，则不会保证ext有效位数为bits位； 如果bits小于等于num的实际位数，则ext将会有bits位有效位数。
返回
剩余的低位bits数量

◆ lmmp_fft_next_size_()

mp_size_t lmmp_fft_next_size_ ( mp_size_t  n)

计算满足 >=n 的最小费马/梅森乘法可行尺寸

参数
n输入的目标尺寸
返回
满足条件的SSA乘法最小尺寸

计算满足 >=n 的最小费马/梅森乘法可行尺寸

参数
n- 原始长度
返回
规整后的长度（为2^k的倍数）

在文件 mul_fft.c 第 84 行定义.

84 {
85 mp_size_t k = lmmp_fft_best_k_(n);
86 lmmp_debug_assert(k >= LOG2_LIMB_BITS);
87 k -= LOG2_LIMB_BITS;
88 n = (((n - 1) >> k) + 1) << k;
89 return n;
90}
k
#define k
lmmp_debug_assert
#define lmmp_debug_assert(x)
Definition lmmp.h:387
LOG2_LIMB_BITS
#define LOG2_LIMB_BITS
Definition lmmp.h:223
lmmp_fft_best_k_
static mp_size_t lmmp_fft_best_k_(mp_size_t n)
查找对于 m>=n 的模 B^m+1 FFT运算的最优k值
Definition mul_fft.c:73

引用了 k, lmmp_debug_assert, lmmp_fft_best_k_() , 以及 LOG2_LIMB_BITS.

被这些函数引用 binvert_mulhi_(), lmmp_div_mulinv_(), lmmp_invappr_newton_(), lmmp_invsqrt_newton_(), lmmp_mul_fft_(), lmmp_mul_fft_unbalance_() , 以及 lmmp_mullo_fft_().

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

◆ lmmp_from_str_()

mp_size_t lmmp_from_str_ ( mp_ptr  dst,
const mp_byte_t *  src,
mp_size_t  len,
int  base 
)

字符串转大数操作 [src,len,base] to [dst,return value,B]

警告
len>=0, 2<=base<=256
参数
dst大数结果输出指针
src字符串源指针
len字符串长度
base字符串的进制基数
返回
转换后的大数 limb 长度

在文件 from_str.c 第 128 行定义.

128 {
129 do {
130 if (len == 0)
131 return 0;
132 } while (src[--len] == 0);
133 ++len;
134
135 mp_size_t limbs;
136 if (LMMP_POW2_Q(base)) {
137 mp_limb_t curlimb = 0;
138 const mp_byte_t* srcend = src + len;
139 int bitspd = lmmp_bases_table[base - 2].large_base;
140 int bitpos = 0;
141 limbs = 0;
142
143 do {
144 mp_limb_t curdigit = *src;
145 curlimb |= curdigit << bitpos;
146 bitpos += bitspd;
147 if (bitpos >= LIMB_BITS) {
148 dst[limbs] = curlimb;
149 ++limbs;
150 bitpos -= LIMB_BITS;
151 curlimb = curdigit >> (bitspd - bitpos);
152 }
153 } while (++src != srcend);
154 if (curlimb) {
155 dst[limbs] = curlimb;
156 ++limbs;
157 }
158 } else if (lmmp_from_str_len_(0, len, base) < FROM_STR_BASEPOW_THRESHOLD) {
159 limbs = lmmp_from_str_basecase_(dst, src, len, base);
160 } else {
161 TEMP_DECL;
162 mp_basepow_t powers[LIMB_BITS];
163 mp_limb_t lbase = lmmp_bases_table[base - 2].large_base;
164 mp_size_t digitspl = lmmp_bases_table[base - 2].digits_in_limb;
165 mp_size_t bexp, lexp = (len - 1) / digitspl + 1;
166 mp_size_t tzbit = lmmp_tailing_zeros_(lbase);
167 // need 1 extra limb to store result
168 mp_size_t alloc_size = lmmp_from_str_len_(0, len, base) + 1;
169 mp_limb_t cy;
170 mp_ptr tp;
171
172 int cpow = lmmp_limb_bits_(lexp - 1);
173
174 for (int i = cpow; i > 0; --i) {
175 // we will calculate lbase^bexp
176 bexp = ((lexp - 1) >> i) + 1;
177 // we will calculate lbase^(bexp-1) first, and trim it s. t.
178 // it contains at most 2 tailing 0 limb, then multiply it by lbase,
179 // so we need npow limbs to store lbase^bexp
180 mp_size_t npow = lmmp_from_str_len_(0, (bexp - 1) * digitspl + 1, base) + 1;
181
182 if (tzbit) {
183 mp_size_t tzlimb = tzbit * (bexp - 1) / LIMB_BITS;
184 if (tzlimb >= 2)
185 npow -= tzlimb - 2;
186 }
187
188 // space needed for a trimmed npow-limb lbase^bexp
189 alloc_size += npow;
190 }
191
192 tp = BALLOC_TYPE(alloc_size, mp_limb_t);
193
194 for (int i = 0; i < 2; ++i) {
195 tp[0] = lbase;
196 powers[i].p = tp;
197 powers[i].np = 1;
198 tp += i + 1;
199 powers[i].zeros = 0;
200 powers[i].digits = digitspl * (i + 1);
201 powers[i].base = base;
202 }
203
204 mp_ptr p = powers[1].p;
205 mp_size_t zeros = 0, np = 1;
206 for (int i = 2; i < cpow; ++i) {
207 lmmp_sqr_(tp, p, np);
208 np *= 2;
209 np -= tp[np - 1] == 0;
210 bexp = (lexp - 1) >> (cpow - i);
211 if (bexp & 1) {
212 cy = lmmp_mul_1_(tp, tp, np, lbase);
213 tp[np] = cy;
214 np += cy != 0;
215 }
216 zeros *= 2;
217 while (tp[0] == 0) {
218 // at most 2 tailing 0 limb here
219 ++zeros;
220 ++tp;
221 --np;
222 }
223 p = tp;
224 powers[i].p = p;
225 powers[i].np = np;
226 powers[i].zeros = zeros;
227 powers[i].digits = digitspl * (bexp + 1);
228 powers[i].base = base;
229 tp += np + 1;
230 }
231
232 for (int i = 1; i < cpow; ++i) {
233 p = powers[i].p;
234 np = powers[i].np;
235 cy = lmmp_mul_1_(p, p, np, lbase);
236 p[np] = cy;
237 np += cy != 0;
238 if (p[0] == 0) {
239 ++powers[i].zeros;
240 ++p;
241 --np;
242 }
243
244 powers[i].p = p;
245 powers[i].np = np;
246 }
247
248 limbs = lmmp_from_str_divide_(tp, src, len, powers + cpow - 1, dst);
249 lmmp_copy(dst, tp, limbs);
250
251 TEMP_FREE;
252 }
253 return limbs;
254}
mp_basepow_t::base
int base
Definition base_table.h:58
mp_base_t::digits_in_limb
int digits_in_limb
Definition base_table.h:37
mp_basepow_t::digits
mp_size_t digits
Definition base_table.h:54
mp_basepow_t::np
mp_size_t np
Definition base_table.h:46
lmmp_bases_table
const mp_base_t lmmp_bases_table[255]
Definition base_table.c:9
mp_base_t::large_base
mp_limb_t large_base
Definition base_table.h:28
mp_basepow_t::zeros
mp_size_t zeros
Definition base_table.h:52
mp_basepow_t::p
mp_ptr p
Definition base_table.h:44
mp_basepow_t
Definition base_table.h:42
lmmp_from_str_basecase_
static mp_size_t lmmp_from_str_basecase_(mp_ptr dst, const mp_byte_t *src, mp_size_t len, int base)
Definition from_str.c:26
lmmp_from_str_divide_
static mp_size_t lmmp_from_str_divide_(mp_ptr dst, const mp_byte_t *src, mp_size_t len, mp_basepow_t *pow, mp_ptr tp)
Definition from_str.c:69
lmmp_from_str_len_
mp_size_t lmmp_from_str_len_(const mp_byte_t *src, mp_size_t len, int base)
计算字符串转大数所需的 limb 缓冲区长度
Definition from_str.c:13
mp_byte_t
uint8_t mp_byte_t
Definition lmmp.h:210
LMMP_POW2_Q
#define LMMP_POW2_Q(n)
Definition lmmp.h:359
lmmp_tailing_zeros_
int lmmp_tailing_zeros_(mp_limb_t x)
计算一个单精度数(limb)中末尾零的个数
Definition tiny.c:54
lmmp_sqr_
void lmmp_sqr_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
大数平方操作 [dst,2*na] = [numa,na]^2
Definition sqr.c:10
lmmp_limb_bits_
int lmmp_limb_bits_(mp_limb_t x)
计算满足 2^k > x 的最小自然数k
Definition tiny.c:11
lmmp_mul_1_
mp_limb_t lmmp_mul_1_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
大数乘以单limb操作 [dst,na] = [numa,na] * x
FROM_STR_BASEPOW_THRESHOLD
#define FROM_STR_BASEPOW_THRESHOLD
Definition mparam.h:74
BALLOC_TYPE
#define BALLOC_TYPE(n, type)
Definition tmp_alloc.h:89

引用了 BALLOC_TYPE, mp_basepow_t::base, mp_basepow_t::digits, mp_base_t::digits_in_limb, FROM_STR_BASEPOW_THRESHOLD, mp_base_t::large_base, LIMB_BITS, lmmp_bases_table, lmmp_copy, lmmp_from_str_basecase_(), lmmp_from_str_divide_(), lmmp_from_str_len_(), lmmp_limb_bits_(), lmmp_mul_1_(), LMMP_POW2_Q, lmmp_sqr_(), lmmp_tailing_zeros_(), mp_basepow_t::np, mp_basepow_t::p, TEMP_DECL, TEMP_FREE, tp , 以及 mp_basepow_t::zeros.

+ 函数调用图:

◆ lmmp_from_str_len_()

mp_size_t lmmp_from_str_len_ ( const mp_byte_t *  src,
mp_size_t  len,
int  base 
)

计算字符串转大数所需的 limb 缓冲区长度

参数
src输入字符串指针
len字符串长度
base字符串的基数（2~256）
返回
存储该字符串数值所需的 limb 缓冲区长度
警告
len>=0, 2<=base<=256
注解
将会忽略前导零，
  1. if (src!=NULL) 返回的长度可能会多分配一个 limb 空间
  2. if (src==NULL) 返回len位base进制数的最大可能 limb 长度（最坏情况）

在文件 from_str.c 第 13 行定义.

13 {
14 lmmp_param_assert(base >= 2 && base <= 256);
15 if (src) {
16 do {
17 if (len == 0)
18 return 1;
19 } while (src[--len] == 0);
20 ++len;
21 }
22 return lmmp_mulh_(len, lmmp_bases_table[base - 2].lg_base) + 1;
23}
lmmp_mulh_
mp_limb_t lmmp_mulh_(mp_limb_t a, mp_limb_t b)
计算两个64位无符号整数相乘的高位结果 (a*b)/2^64
Definition tiny.c:73

引用了 lmmp_bases_table, lmmp_mulh_() , 以及 lmmp_param_assert.

被这些函数引用 lmmp_from_str_(), lmmp_from_str_divide_() , 以及 lmmp_to_str_().

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

◆ lmmp_inv_()

void lmmp_inv_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  nf 
)

大数求逆操作 [dst,na+nf+1] = (B^(2*(na+nf)) - 1) / ([numa,na]*B^nf) + [0|-1]

参数
dst逆元结果输出指针
numa源操作数指针
na操作数的 limb 长度
nf精度因子
警告
na>0, numa[na-1]!=0, eqsep(dst,numa)

在文件 inv.c 第 152 行定义.

152 {
153 lmmp_param_assert(na > 0);
154 lmmp_param_assert(numa[na - 1] != 0);
155 mp_limb_t high = numa[na - 1];
156 int nsh = lmmp_leading_zeros_(high);
157 TEMP_DECL;
158 if (dst == numa || nsh || nf) {
159 nf += nsh != 0;
160 mp_ptr numa2 = TALLOC_TYPE(na + nf, mp_limb_t);
161 lmmp_zero(numa2, nf);
162 if (nsh)
163 lmmp_shl_(numa2 + nf, numa, na, nsh);
164 else
165 lmmp_copy(numa2 + nf, numa, na);
166 numa = numa2;
167 }
168 lmmp_invappr_(dst, numa, na + nf);
169 if (nsh)
170 lmmp_shr_c_(dst, dst, na + nf, LIMB_BITS - nsh, (mp_limb_t)1 << nsh);
171 else
172 dst[na + nf] = 1;
173 TEMP_FREE;
174}
lmmp_invappr_
void lmmp_invappr_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
近似逆元计算 (invappr)
Definition inv.c:176
lmmp_zero
#define lmmp_zero(dst, n)
Definition lmmp.h:366
lmmp_shr_c_
mp_limb_t lmmp_shr_c_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_size_t shr, mp_limb_t c)
带进位的大数右移操作 [dst,na] = [numa,na]>>shr，dst的高shr位填充c的高shr位
Definition shr.c:30

引用了 LIMB_BITS, lmmp_copy, lmmp_invappr_(), lmmp_leading_zeros_(), lmmp_param_assert, lmmp_shl_(), lmmp_shr_c_(), lmmp_zero, TALLOC_TYPE, TEMP_DECL , 以及 TEMP_FREE.

+ 函数调用图:

◆ lmmp_inv_1_()

mp_limb_t lmmp_inv_1_ ( mp_limb_t  x)

1阶逆元计算 (inv1)

参数
x输入的64位无符号整数，最高位为1(MSB(x)=1)
返回
计算结果：(B^2-1)/x - B
警告
MSB(x)=1, 即x>=2^63

在文件 inv.c 第 107 行定义.

107 {
108 mp_limb_t r, m;
109 {
110 mp_limb_t p, ql;
111 unsigned ul, uh, qh;
112
113 ul = x & LLIMB_MASK;
114 uh = x >> (LIMB_BITS / 2);
115 qh = (x ^ LIMB_MAX) / uh;
116
117 r = ((~x - (mp_limb_t)qh * uh) << (LIMB_BITS / 2)) | LLIMB_MASK;
118 p = (mp_limb_t)qh * ul;
119 if (r < p) {
120 qh--;
121 r += x;
122 if (r >= x)
123 if (r < p) {
124 qh--;
125 r += x;
126 }
127 }
128 r -= p;
129 p = (r >> (LIMB_BITS / 2)) * qh + r;
130 ql = (p >> (LIMB_BITS / 2)) + 1;
131 r = (r << (LIMB_BITS / 2)) + LLIMB_MASK - ql * x;
132 if (r >= (LIMB_MAX & (p << (LIMB_BITS / 2)))) {
133 ql--;
134 r += x;
135 }
136 m = ((mp_limb_t)qh << (LIMB_BITS / 2)) + ql;
137 if (r >= x) {
138 m++;
139 r -= x;
140 }
141 }
142 return m;
143}
LLIMB_MASK
#define LLIMB_MASK
Definition lmmp.h:227

引用了 LIMB_BITS, LIMB_MAX , 以及 LLIMB_MASK.

被这些函数引用 lmmp_div_1_(), lmmp_div_1_s_(), lmmp_inv_basecase_(), lmmp_mod_1_() , 以及 lmmp_mulmod_ulong_().

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

◆ lmmp_inv_2_1_()

mp_limb_t lmmp_inv_2_1_ ( mp_limb_t  xh,
mp_limb_t  xl 
)

2-1阶逆元计算 (inv21)

参数
xh输入数的高64位部分
xl输入数的低64位部分
返回
计算结果：(B^3-1)/(xh*B+xl) - B
警告
MSB(xh)=1, 即xh>=2^63

在文件 inv.c 第 10 行定义.

10 {
11 mp_limb_t r, m;
12 {
13 mp_limb_t p, ql;
14 unsigned ul, uh, qh;
15
16 /* For notation, let b denote the half-limb base, so that B = b^2.
17 Split u1 = b uh + ul. */
18 ul = xh & LLIMB_MASK;
19 uh = xh >> (LIMB_BITS / 2);
20
21 /* Approximation of the high half of quotient. Differs from the 2/1
22 inverse of the half limb uh, since we have already subtracted
23 u0. */
24 qh = (xh ^ LIMB_MAX) / uh;
25
26 /* Adjust to get a half-limb 3/2 inverse, i.e., we want
27
28 qh' = floor( (b^3 - 1) / u) - b = floor ((b^3 - b u - 1) / u
29 = floor( (b (~u) + b-1) / u),
30
31 and the remainder
32
33 r = b (~u) + b-1 - qh (b uh + ul)
34 = b (~u - qh uh) + b-1 - qh ul
35
36 Subtraction of qh ul may underflow, which implies adjustments.
37 But by normalization, 2 u >= B > qh ul, so we need to adjust by
38 at most 2.
39 */
40
41 r = ((~xh - (mp_limb_t)qh * uh) << (LIMB_BITS / 2)) | LLIMB_MASK;
42
43 p = (mp_limb_t)qh * ul;
44 /* Adjustment steps taken from udiv_qrnnd_c */
45 if (r < p) {
46 qh--;
47 r += xh;
48 if (r >= xh) /* i.e. we didn't get carry when adding to r */
49 if (r < p) {
50 qh--;
51 r += xh;
52 }
53 }
54 r -= p;
55
56 /* Low half of the quotient is
57
58 ql = floor ( (b r + b-1) / u1).
59
60 This is a 3/2 division (on half-limbs), for which qh is a
61 suitable inverse. */
62
63 p = (r >> (LIMB_BITS / 2)) * qh + r;
64 /* Unlike full-limb 3/2, we can add 1 without overflow. For this to
65 work, it is essential that ql is a full mp_limb_t. */
66 ql = (p >> (LIMB_BITS / 2)) + 1;
67
68 /* By the 3/2 trick, we don't need the high half limb. */
69 r = (r << (LIMB_BITS / 2)) + LLIMB_MASK - ql * xh;
70
71 if (r >= (LIMB_MAX & (p << (LIMB_BITS / 2)))) {
72 ql--;
73 r += xh;
74 }
75 m = ((mp_limb_t)qh << (LIMB_BITS / 2)) + ql;
76 if (r >= xh) {
77 m++;
78 r -= xh;
79 }
80 }
81
82 /* Now m is the 2/1 inverse of u1. If u0 > 0, adjust it to become a
83 3/2 inverse. */
84 if (xl > 0) {
85 mp_limb_t th, tl;
86 r = ~r;
87 r += xl;
88 if (r < xl) {
89 m--;
90 if (r >= xh) {
91 m--;
92 r -= xh;
93 }
94 r -= xh;
95 }
96 _umul64to128_(xl, m, &tl, &th);
97 r += th;
98 if (r < th) {
99 m--;
100 m -= ((r > xh) | ((r == xh) & (tl > xl)));
101 }
102 }
103
104 return m;
105}
_umul64to128_
static void _umul64to128_(uint64_t a, uint64_t b, uint64_t *low, uint64_t *high)
Definition longlong.h:31

引用了 _umul64to128_(), LIMB_BITS, LIMB_MAX , 以及 LLIMB_MASK.

被这些函数引用 lmmp_bninv_appr_newton_(), lmmp_div_(), lmmp_div_2_(), lmmp_div_2_s_(), lmmp_div_s_(), lmmp_inv_basecase_() , 以及 lmmp_mod_2_().

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

◆ lmmp_inv_basecase_()

void lmmp_inv_basecase_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

近似逆元计算

参数
dst输出结果缓冲区，长度为na
numa输入操作数，长度为na
na输入操作数的 limb 长度
警告
na>0, MSB(numa)=1, sep(dst,numa)
返回
无返回值，结果存储在dst中，[dst,na]=(B^(2*na)-1)/[numa,na] - B^na

在文件 inv.c 第 11 行定义.

11 {
12 lmmp_param_assert(na > 0);
13 lmmp_param_assert(numa[na - 1] >= LIMB_B_2);
14 if (na == 1)
15 *dst = lmmp_inv_1_(*numa);
16 else {
17 TEMP_DECL;
18 mp_ptr xp = TALLOC_TYPE(2 * na, mp_limb_t);
19 mp_size_t i = na;
20 do {
21 xp[--i] = LIMB_MAX;
22 } while (i);
23 lmmp_not_(xp + na, numa, na);
24 //[xp,2*na] = B^(2*na)-1 - [numa,na]*B^na
25
26 if (na == 2) {
27 lmmp_div_2_s_(dst, xp, 4, numa);
28 } else {
29 mp_limb_t inv21 = lmmp_inv_2_1_(numa[na - 1], numa[na - 2]);
30 if (na < DIV_DIVIDE_THRESHOLD) {
31 lmmp_div_basecase_(dst, xp, 2 * na, numa, na, inv21);
32 } else {
33 lmmp_div_divide_(dst, xp, 2 * na, numa, na, inv21);
34 }
35 }
36 TEMP_FREE;
37 }
38}
lmmp_inv_1_
mp_limb_t lmmp_inv_1_(mp_limb_t x)
1阶逆元计算 (inv1)
Definition inv.c:107
lmmp_inv_2_1_
mp_limb_t lmmp_inv_2_1_(mp_limb_t xh, mp_limb_t xl)
2-1阶逆元计算 (inv21)
Definition inv.c:10
lmmp_not_
void lmmp_not_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
大数按位取反操作 [dst,na] = ~[numa,na] (对每个limb执行按位非操作)

引用了 DIV_DIVIDE_THRESHOLD, LIMB_B_2, LIMB_MAX, lmmp_div_2_s_(), lmmp_div_basecase_(), lmmp_div_divide_(), lmmp_inv_1_(), lmmp_inv_2_1_(), lmmp_not_(), lmmp_param_assert, TALLOC_TYPE, TEMP_DECL , 以及 TEMP_FREE.

被这些函数引用 lmmp_invappr_() , 以及 lmmp_invappr_newton_().

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

◆ lmmp_inv_prediv_()

void lmmp_inv_prediv_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  ni 
)

除法前的逆元预计算，[dst,ni] = invappr( (ni+1 MSLs of numa) + 1 ) / B

参数
dst输出预计算逆元的缓冲区，长度为ni
numa输入操作数，长度为na
na输入操作数的 limb 长度
ni预计算逆元的目标尺寸
警告
na>=ni>0, MSB(numa)=1, eqsep(dst,numa)
注解
if (ni=na) [dst,na] = (B^(2*na)-1) / [numa,na] - B^na

在文件 div_mulinv.c 第 11 行定义.

11 {
12 lmmp_param_assert(na >= ni);
13 lmmp_param_assert(ni > 0);
14 lmmp_param_assert(numa[na - 1] >= LIMB_B_2);
15 TEMP_DECL;
16 mp_limb_t cy;
17 mp_ptr tp = TALLOC_TYPE(ni + 1, mp_limb_t);
18
19 if (na == ni) {
20 lmmp_copy(tp + 1, numa, ni);
21 tp[0] = 1;
22 cy = 0;
23 } else {
24 cy = lmmp_add_1_(tp, numa + na - (ni + 1), ni + 1, 1);
25 }
26 if (cy)
27 lmmp_zero(dst, ni);
28 else {
29 mp_ptr invappr = TALLOC_TYPE(ni + 1, mp_limb_t);
30 lmmp_invappr_(invappr, tp, ni + 1);
31 lmmp_copy(dst, invappr + 1, ni);
32 }
33 TEMP_FREE;
34}
lmmp_invappr_
void lmmp_invappr_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
近似逆元计算 (invappr)
Definition inv.c:176
lmmp_add_1_
static mp_limb_t lmmp_add_1_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_limb_t x)
大数加单精度数静态内联函数 [dst,na]=[numa,na]+x
Definition lmmpn.h:1111

引用了 LIMB_B_2, lmmp_add_1_(), lmmp_copy, lmmp_invappr_(), lmmp_param_assert, lmmp_zero, TALLOC_TYPE, TEMP_DECL, TEMP_FREE , 以及 tp.

被这些函数引用 lmmp_div_(), lmmp_div_s_() , 以及 lmmp_to_str_().

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

◆ lmmp_invappr_()

void lmmp_invappr_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

近似逆元计算 (invappr)

参数
dst输出结果缓冲区，长度为na
numa输入操作数，长度为na
na输入操作数的 limb 长度
警告
na>0, MSB(numa)=1, sep(dst,numa)
返回
无返回值，结果存储在dst中，[dst,na] = (B^(2*na)-1)/[numa,na] - B^na + [0|-1]

在文件 inv.c 第 176 行定义.

176 {
177 if (na < INV_NEWTON_THRESHOLD)
178 lmmp_inv_basecase_(dst, numa, na);
179 else
180 lmmp_invappr_newton_(dst, numa, na);
181}
lmmp_invappr_newton_
void lmmp_invappr_newton_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
近似逆元计算（牛顿迭代法）
Definition inv.c:40
lmmp_inv_basecase_
void lmmp_inv_basecase_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
近似逆元计算
Definition inv.c:11
INV_NEWTON_THRESHOLD
#define INV_NEWTON_THRESHOLD
Definition mparam.h:33

引用了 INV_NEWTON_THRESHOLD, lmmp_inv_basecase_() , 以及 lmmp_invappr_newton_().

被这些函数引用 lmmp_inv_() , 以及 lmmp_inv_prediv_().

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

◆ lmmp_invappr_newton_()

void lmmp_invappr_newton_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

近似逆元计算（牛顿迭代法）

参数
dst输出结果缓冲区，长度为na
numa输入操作数，长度为na
na输入操作数的 limb 长度
警告
na>4, MSB(numa)=1, sep(dst,numa)
返回
无返回值，结果存储在dst中，[dst,na]=(B^(2*na)-1)/[numa,na]-B^na+[0|-1]

在文件 inv.c 第 40 行定义.

40 {
41 lmmp_param_assert(na > 4);
42 lmmp_param_assert(numa[na - 1] >= LIMB_B_2);
43
44 mp_limb_t cy;
45 mp_size_t nr = na, mn;
46 mp_size_t sizes[LIMB_BITS], *sizp = sizes;
47
48 do {
49 *sizp = nr;
50 nr = (nr >> 1) + 1;
51 ++sizp;
52 } while (nr >= INV_NEWTON_THRESHOLD);
53
54 numa += na;
55 dst += na;
56
57 lmmp_inv_basecase_(dst - nr, numa - nr, nr);
58
59 TEMP_DECL;
60 mp_ptr xp = TALLOC_TYPE(3 * (na >> 1) + 3, mp_limb_t);
61 do {
62 na = *--sizp;
63
64 // ar = 0:[numa-nr,nr]
65 // an = 0:[numa-na,na]
66 // ir = 1:[dst-nr,nr] = (B^(2*nr)-1)/ar - [0|1]
67 // rem = ir*an-B^(na+nr)
68 //-2*B^na < rem < 2*B^na
69
70 //[xp] = rem
71 if (na < INV_MODM_THRESHOLD || (mn = lmmp_fft_next_size_(na + 1)) >= na + nr) {
72 lmmp_mul_(xp, numa - na, na, dst - nr, nr);
73 lmmp_add_n_(xp + nr, xp + nr, numa - na, na + 1 - nr);
74 cy = 1; // for mod B^(na+1)
75 } else { // nr < na < mn < na+nr
76
77 //[xp,mn] = [dst,nr] * [numa,na] mod (B^mn-1)
78 lmmp_mul_mersenne_(xp, mn, numa - na, na, dst - nr, nr);
79
80 //[xp,mn] += [numa,na]*B^nr mod (B^mn-1)
81 cy = lmmp_add_n_(xp + nr, xp + nr, numa - na, mn - nr);
82 cy = lmmp_add_nc_(xp, xp, numa - (na - (mn - nr)), na - (mn - nr), cy);
83
84 //[xp,mn] -= B^(na+nr) mod (B^mn-1)
85 xp[mn] = 1;
86 lmmp_dec_1(xp + na + nr - mn, 1 - cy);
87 lmmp_dec_1(xp, 1 - xp[mn]);
88
89 cy = 0; // for mod (B^mn-1)
90 }
91
92 // adjust ir,rem s.t.
93 // -B^na < rem = ir*an - B^(na+nr) < 0
94 // use this we can prove B^nr <= ir < 2*B^nr
95 // so inc/dec ir won't overflow
96 if (xp[na] < 2) { // rem>=0
97
98 // rem-=cy*an s.t. rem[na]=0
99 if (cy = xp[na]) {
100 if (!lmmp_sub_n_(xp, xp, numa - na, na)) {
101 ++cy;
102 lmmp_sub_n_(xp, xp, numa - na, na);
103 }
104 }
105
106 // rem-=cy*an s.t. 0<=rem<an
107 if (lmmp_cmp_(xp, numa - na, na) >= 0) {
108 lmmp_sub_n_(xp, xp, numa - na, na);
109 ++cy;
110 }
111
112 // 0 < an-rem <= an < B^na , trunc to nr limbs
113 lmmp_sub_nc_(xp + 2 * nr, numa - nr, xp + na - nr, nr, lmmp_cmp_(xp, numa - na, na - nr) > 0);
114 ++cy;
115
116 lmmp_dec_1(dst - nr, cy);
117 } else { // rem<0
118 if (cy)
119 lmmp_dec(xp); // for neg to not
120 // else (neg to not) compensate (mod transfer)
121
122 if (xp[na] != LIMB_MAX) {
123 lmmp_assert(xp[na] + lmmp_add_n_(xp, xp, numa - na, na) == LIMB_MAX);
124 lmmp_inc(dst - nr);
125 }
126
127 //-rem
128 lmmp_not_(xp + 2 * nr, xp + na - nr, nr);
129 }
130
131 // in = 1:[dst-na,na]
132 // in = ir*B^(na-nr) + ir*(-rem/B^(na-nr))/B^(3*nr-na)
133 // use inequality an*ir!=B^(na+nr),
134 //(otherwise obviously contradictory),
135 // we can prove
136 // an*in <= an*ir * ( 2*B^(na+nr) - an*ir ) * B^(-2*nr) < B^(2*na)
137 // so in < B^(2*na)/an <= 2*B^(na),
138 // inc below won't overflow
139
140 // and via inequality -B^na < an*ir - B^(na+nr) < 0
141 // we can prove in = (B^(2*na)-1)/an - [0|1]
142 lmmp_mul_n_(xp, xp + 2 * nr, dst - nr, nr);
143 cy = lmmp_add_n_(xp + nr, xp + nr, xp + 2 * nr, 2 * nr - na);
144 if (lmmp_add_nc_(dst - na, xp + 3 * nr - na, xp + 4 * nr - na, na - nr, cy))
145 lmmp_inc(dst - nr);
146
147 nr = na;
148 } while (sizp != sizes);
149 TEMP_FREE;
150}
lmmp_dec_1
#define lmmp_dec_1(p, dec)
大数减指定值宏（预期无借位）
Definition lmmpn.h:985
INV_MODM_THRESHOLD
#define INV_MODM_THRESHOLD
Definition mparam.h:35

引用了 INV_MODM_THRESHOLD, INV_NEWTON_THRESHOLD, LIMB_B_2, LIMB_BITS, LIMB_MAX, lmmp_add_n_(), lmmp_add_nc_(), lmmp_assert, lmmp_cmp_(), lmmp_dec, lmmp_dec_1, lmmp_fft_next_size_(), lmmp_inc, lmmp_inv_basecase_(), lmmp_mul_(), lmmp_mul_mersenne_(), lmmp_mul_n_(), lmmp_not_(), lmmp_param_assert, lmmp_sub_n_(), lmmp_sub_nc_(), TALLOC_TYPE, TEMP_DECL , 以及 TEMP_FREE.

被这些函数引用 lmmp_invappr_().

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

◆ lmmp_leading_zeros_()

int lmmp_leading_zeros_ ( mp_limb_t  x)

计算一个单精度数(limb)中前导零的个数

参数
x输入的64位无符号整数
返回
前导零的位数（范围：0~64）

在文件 tiny.c 第 35 行定义.

35 {
36 if (x == 0) return 64;
37#ifdef __GNUC__
38 return __builtin_clzll(x);
39#elif _MSC_VER
40 unsigned long index;
41 _BitScanReverse64(&index, x);
42 return 63 - (int)index;
43#else
44 int n = 0;
45 if (x <= 0x00000000FFFFFFFF) { n += 32; x <<= 32; }
46 if (x <= 0x0000FFFFFFFFFFFF) { n += 16; x <<= 16; }
47 if (x <= 0x00FFFFFFFFFFFFFF) { n += 8; x <<= 8; }
48 if (x <= 0x0FFFFFFFFFFFFFFF) { n += 4; x <<= 4; }
49 if (x <= 0x3FFFFFFFFFFFFFFF) { n += 2; x <<= 2; }
50 if (x <= 0x7FFFFFFFFFFFFFFF) { n += 1; x <<= 1; }
51#endif
52}

被这些函数引用 lmmp_bninv_(), lmmp_div_(), lmmp_div_1_(), lmmp_div_2_(), lmmp_inv_(), lmmp_mod_1_(), lmmp_mod_2_(), lmmp_mulmod_ulong_(), lmmp_pow_basecase_(), lmmp_sqrt_(), lmmp_sqrt_newton_(), lmmp_to_str_() , 以及 mont63_R2().

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

◆ lmmp_limb_bits_()

int lmmp_limb_bits_ ( mp_limb_t  x)

计算满足 2^k > x 的最小自然数k

参数
x输入的64位无符号整数
返回
满足条件的最小自然数k

在文件 tiny.c 第 11 行定义.

11 {
12 int k = 0;
13 while (x) {
14 x >>= 1;
15 k++;
16 }
17 return k;
18}

引用了 k.

被这些函数引用 lmmp_extract_bits_(), lmmp_from_str_(), lmmp_lehmer_extract_(), lmmp_pow_(), lmmp_to_str_() , 以及 lmmp_to_str_len_().

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

◆ lmmp_limb_popcnt_()

int lmmp_limb_popcnt_ ( mp_limb_t  x)

计算一个64位无符号整数中1的个数

参数
x输入的64位无符号整数
返回
1的个数

在文件 tiny.c 第 20 行定义.

20 {
21#ifdef __GNUC__
22 return __builtin_popcountll(x);
23#elif _MSC_VER
24 return (int)__popcnt64(x);
25#else
26 int k = 0;
27 while (x) {
28 k += x & 1;
29 x >>= 1;
30 }
31 return k;
32#endif
33}

引用了 k.

被这些函数引用 lmmp_factorial_size_(), lmmp_multinomial_(), lmmp_nCr_size_(), lmmp_nPr_size_() , 以及 lmmp_pow_().

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

◆ lmmp_mod_1_()

mp_limb_t lmmp_mod_1_ ( mp_srcptr  numa,
mp_size_t  na,
mp_limb_t  x 
)

单精度数取余

参数
numa输入被除数，长度为na
na被除数的 limb 长度
x除数（单个 limb ）
返回
除法余数（单个 limb ）
警告
na>0, x!=0, eqsep(dstq,numa), dstq>=numa-1 是可以接受的

在文件 div.c 第 20 行定义.

20 {
21 mp_limb_t ah, al;
22 // q: assigned for macro reuse, unused in this logic (known warning)
23 mp_limb_t t = numa[na - 2], q = 0, r = 0;
24 const int shift = lmmp_leading_zeros_(x);
25 if (shift > 0) {
26 const int rshift = LIMB_BITS - shift;
27 ah = numa[na - 1] >> rshift;
28 t = numa[na - 2];
29 al = (numa[na - 1] << shift) | (t >> rshift);
30 x <<= shift;
31 const mp_limb_t inv = lmmp_inv_1_(x);
32 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
33 na -= 2;
34 while (na-- > 0) {
35 ah = r;
36 al = t << shift;
37 t = numa[na];
38 al |= t >> rshift;
39 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
40 }
41 ah = r;
42 al = t << shift;
43 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
44 return r >> shift;
45 } else {
46 ah = 0;
47 t = numa[na - 2];
48 al = numa[na - 1];
49 const mp_limb_t inv = lmmp_inv_1_(x);
50 q = al / x;
51 r = al % x;
52 na -= 2;
53 while (na-- > 0) {
54 ah = r;
55 al = t;
56 t = numa[na];
57 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
58 }
59 ah = r;
60 al = t;
61 _udiv_qrnnd_preinv(q, r, ah, al, x, inv);
62 return r;
63 }
64}

引用了 _udiv_qrnnd_preinv, LIMB_BITS, lmmp_inv_1_() , 以及 lmmp_leading_zeros_().

被这些函数引用 lmmp_div_1_(), lmmp_gcd_1_() , 以及 lmmp_trialdiv_short_().

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

◆ lmmp_mod_2_()

void lmmp_mod_2_ ( mp_srcptr  numa,
mp_size_t  na,
mp_ptr  numb 
)

双精度数取余 (除数为2个limb)

参数
numa输入被除数（长度na）
na被除数的 limb 长度
numb输入除数（长度2）[numb,2]=[numa,na] mod [numb,2]
警告
na>=2, numb[1]!=0, eqsep(dstq,numa), dstq>=numa 是可以接受的

在文件 div.c 第 144 行定义.

144 {
145 mp_limb_t q, r1, r0, a2, a1, a0, b1, b0;
146 b1 = numb[1];
147 b0 = numb[0];
148 if (na == 2) {
149 int shift = lmmp_leading_zeros_(b1);
150 if (shift > 0) {
151 const int rshift = LIMB_BITS - shift;
152 b1 = (b1 << shift) | (b0 >> rshift);
153 b0 <<= shift;
154 a2 = numa[1] >> rshift;
155 a1 = (numa[1] << shift) | (numa[0] >> rshift);
156 a0 = (numa[0] << shift);
157 mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
158 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
159 numb[0] = (r0 >> shift) | (r1 << rshift);
160 numb[1] = r1 >> shift;
161 return;
162 } else {
163 if (_u128cmp(numa, numb)) {
164 numb[0] = numa[0];
165 numb[1] = numa[1];
166 return;
167 } else {
168 _u128sub(numb, numa, numb);
169 return;
170 }
171 }
172 } else {
173 int shift = lmmp_leading_zeros_(b1);
174 if (shift > 0) {
175 const int rshift = LIMB_BITS - shift;
176 b1 = (b1 << shift) | (b0 >> rshift);
177 b0 <<= shift;
178 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
179 a2 = numa[na - 1] >> rshift;
180 a1 = (numa[na - 1] << shift) | (numa[na - 2] >> rshift);
181 a0 = (numa[na - 2] << shift) | (numa[na - 3] >> rshift);
182 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
183 na -= 2;
184 while (na-- > 1) {
185 a2 = r1;
186 a1 = r0;
187 a0 = (numa[na] << shift) | (numa[na - 1] >> rshift);
188 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
189 }
190
191 a2 = r1;
192 a1 = r0;
193 a0 = (numa[na] << shift);
194 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
195 numb[0] = (r0 >> shift) | (r1 << rshift);
196 numb[1] = r1 >> shift;
197 return;
198 } else {
199 const mp_limb_t inv = lmmp_inv_2_1_(b1, b0);
200 a2 = 0;
201 a1 = numa[na - 1];
202 a0 = numa[na - 2];
203 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
204 na -= 2;
205 while (na-- > 1) {
206 a2 = r1;
207 a1 = r0;
208 a0 = numa[na];
209 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
210 }
211 a2 = r1;
212 a1 = r0;
213 a0 = numa[na];
214 _udiv_qr_3by2(q, r1, r0, a2, a1, a0, b1, b0, inv);
215 numb[0] = r0;
216 numb[1] = r1;
217 return;
218 }
219 }
220}

引用了 _u128cmp, _u128sub, _udiv_qr_3by2, a0, a1, a2, b0, b1, LIMB_BITS, lmmp_inv_2_1_(), lmmp_leading_zeros_(), r0 , 以及 r1.

被这些函数引用 lmmp_gcd_2_().

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

◆ lmmp_mul_()

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]

警告
0<nb<=na, sep(dst,[numa|numb]) 特殊情况: nb==1时dst<=numa+1是允许的 nb==2时dst<=numa是允许的
参数
dst乘积结果输出指针（需要 na+nb 的 limb 长度）
numa第一个乘数指针（较长的操作数）
na第一个操作数的 limb 长度
numb第二个乘数指针（较短的操作数）
nb第二个操作数的 limb 长度

被这些函数引用 lmmp_bninv_appr_newton_(), lmmp_div_(), lmmp_div_divide_n_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_elem_mul_ulong_(), lmmp_from_str_divide_(), lmmp_invappr_newton_(), lmmp_invsqrt_newton_(), lmmp_mul_fft_unbalance_(), lmmp_mul_signed_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom42_unbalance_(), lmmp_mul_toom43_(), lmmp_mul_toom44_(), lmmp_mul_toom52_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_mul_toom62_unbalance_(), lmmp_num_heap_mul_(), lmmp_pow_win2_(), lmmp_sqrt_newton_() , 以及 pow_nPr_().

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

◆ lmmp_mul_1_()

mp_limb_t lmmp_mul_1_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_limb_t  x 
)

大数乘以单limb操作 [dst,na] = [numa,na] * x

参数
dst结果输出指针
numa被乘数指针
na操作数的位数（limb数量）
x单个limb乘数
警告
na>0, eqsep(dst,numa) 支持 dst<=numa+1 的内存布局
返回
运算后的进位limb值

被这些函数引用 lmmp_3pow_1_(), lmmp_elem_mul_ulong_(), lmmp_from_str_(), lmmp_from_str_basecase_(), lmmp_lehmer_mul_(), lmmp_mullo_dc_(), lmmp_odd_nCr_uint_(), lmmp_odd_nCr_ushort_(), lmmp_pow_(), lmmp_sqrlo_dc_(), lmmp_to_str_(), lmmp_u32_pow_1_(), lmmp_u64_pow_1_() , 以及 lmmp_u8_pow_1_().

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

◆ lmmp_mul_basecase_()

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

基础乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
0<nb<=na, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_(), lmmp_mul_n_(), lmmp_u16_pow_1_() , 以及 lmmp_vec_elem_mul_().

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

◆ lmmp_mul_fermat_()

void lmmp_mul_fermat_ ( mp_ptr  dst,
mp_size_t  rn,
mp_srcptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)

费马数模乘法 [dst,rn+1]=[numa,na]*[numb,nb] mod B^rn+1

参数
dst输出结果缓冲区，长度至少为 rn+1
rn模运算的阶数参数，rn = lmmp_fft_next_size_((na + nb + 1) >> 1)
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
0<=[numa,na]<2*B^rn, 0<=[numb,nb]<2*B^rn, rn = lmmp_fft_next_size_((na + nb + 1) >> 1)
返回
无返回值，结果存储在dst中

在文件 mul_fft.c 第 677 行定义.

677 {
678 int nsqr = numa != numb || na != nb; // 判断是否为平方运算
679 mp_size_t N = rn * LIMB_BITS; // 结果总比特数
680 mp_size_t k = lmmp_fft_best_k_(rn); // 最优FFT层数
681 mp_size_t K = ((mp_size_t)1) << k; // FFT块数（2^k）
682 lmmp_debug_assert(!(N & (K - 1)));
683 mp_size_t M = N >> k; // 每个块的比特数
684 mp_size_t n = 2 * M + k + 2; // 扩展系数长度
685
686 n = (n + LIMB_BITS - 1) & (-LIMB_BITS);
687 n = (((n - 1) >> k) + 1) << k;
688
689 // 初始化内存栈
690 fft_memstack msr;
691 msr.maxdepth = -1;
692 msr.tempdepth = -1;
693 msr.lenw = n / LIMB_BITS; // 系数长度（机器字）
694 mp_size_t nlen = msr.lenw + 1; // 系数总长度
695
696 msr.temp_coef = (mp_ptr)lmmp_fft_memstack_(&msr, (((nlen + 1) << (k + nsqr)) + nlen) * LIMB_BYTES);
697
698 mp_ptr *pfca = (mp_ptr*)(msr.temp_coef + nlen), *pfcb = pfca;
699 mp_size_t narest = na * LIMB_BITS, nbrest = nb * LIMB_BITS;
700
701 for (mp_size_t i = 0; i < K; ++i) {
702 mp_size_t coeflen;
703 pfca[i] = (mp_ptr)(pfca + K) + i * nlen;
704 if (narest > 0) {
705 coeflen = M + (i == K - 1);
706 coeflen = LMMP_MIN(narest, coeflen);
707 narest -= coeflen;
708 lmmp_fft_extract_coef_(pfca[i], numa, M * i, coeflen, msr.lenw);
709 // 非第一个块：左移补偿
710 if (i > 0)
711 lmmp_fft_shl_coef_(&msr, pfca + i, i * n >> k);
712 } else {
713 lmmp_zero(pfca[i], nlen);
714 }
715 }
716 lmmp_fft_(&msr, pfca, k, n >> (k - 1));
717
718 if (nsqr) {
719 pfcb += (nlen + 1) << k;
720 for (mp_size_t i = 0; i < K; ++i) {
721 mp_size_t coeflen;
722 pfcb[i] = (mp_ptr)(pfcb + K) + i * nlen;
723 if (nbrest > 0) {
724 coeflen = M + (i == K - 1);
725 coeflen = LMMP_MIN(nbrest, coeflen);
726 nbrest -= coeflen;
727 lmmp_fft_extract_coef_(pfcb[i], numb, M * i, coeflen, msr.lenw);
728 if (i > 0)
729 lmmp_fft_shl_coef_(&msr, pfcb + i, i * n >> k);
730 } else {
731 lmmp_zero(pfcb[i], nlen);
732 }
733 }
734 lmmp_fft_(&msr, pfcb, k, n >> (k - 1));
735 }
736
737 lmmp_mul_fermat_recurse_(&msr, pfca, pfcb, K);
738
739 lmmp_ifft_(&msr, pfca, k, n >> (k - 1));
740
741 lmmp_mul_fermat_recombine_(&msr, dst, pfca, K, k, n, M, rn);
742
743 // 处理模 B^rn+1 的溢出
744 if (dst[rn] && !lmmp_zero_q_(dst, rn)) {
745 dst[rn] = 0;
746 lmmp_dec(dst);
747 }
748
749 lmmp_fft_memstack_(&msr, 0);
750}
lmmp_zero_q_
static int lmmp_zero_q_(mp_srcptr p, mp_size_t n)
大数判零函数（内联）
Definition lmmpn.h:1027
LIMB_BYTES
#define LIMB_BYTES
Definition mparam.h:85
fft_memstack::temp_coef
mp_ptr temp_coef
Definition mul_fft.c:60
fft_memstack::maxdepth
mp_ssize_t maxdepth
Definition mul_fft.c:62
lmmp_mul_fermat_recurse_
static void lmmp_mul_fermat_recurse_(fft_memstack *ms, mp_ptr *pc1, mp_ptr *pc2, mp_size_t K0)
费马变换乘法递归函数（核心乘法逻辑）
Definition mul_fft.c:600
lmmp_fft_
static void lmmp_fft_(fft_memstack *ms, mp_ptr *coef, mp_size_t k, mp_size_t w)
Definition mul_fft.c:442
lmmp_fft_memstack_
static void * lmmp_fft_memstack_(fft_memstack *ms, mp_size_t size)
FFT内存栈的分配/释放接口
Definition mul_fft.c:98
lmmp_fft_shl_coef_
static void lmmp_fft_shl_coef_(fft_memstack *ms, mp_ptr *coef, mp_size_t shl)
对模 2^n+1 的系数执行左移操作
Definition mul_fft.c:150
lmmp_mul_fermat_recombine_
static void lmmp_mul_fermat_recombine_(fft_memstack *ms, mp_ptr dst, mp_ptr *pfca, mp_size_t K, mp_size_t k, mp_size_t n, mp_size_t M, mp_size_t rn)
费马变换 模 B^n+1 乘法的结果合并
Definition mul_fft.c:517
lmmp_ifft_
static void lmmp_ifft_(fft_memstack *ms, mp_ptr *coef, mp_size_t k, mp_size_t w)
Definition mul_fft.c:495
lmmp_fft_extract_coef_
static void lmmp_fft_extract_coef_(mp_ptr dst, mp_srcptr numa, mp_size_t bitoffset, mp_size_t bits, mp_size_t lenw)
[dst,lenw+1] = [(bit*)numa+bitoffset, bits]
Definition mul_fft.c:124
fft_memstack::tempdepth
mp_ssize_t tempdepth
Definition mul_fft.c:63
fft_memstack::lenw
mp_size_t lenw
Definition mul_fft.c:61
fft_memstack
Definition mul_fft.c:59

引用了 k, fft_memstack::lenw, LIMB_BITS, LIMB_BYTES, lmmp_debug_assert, lmmp_dec, lmmp_fft_(), lmmp_fft_best_k_(), lmmp_fft_extract_coef_(), lmmp_fft_memstack_(), lmmp_fft_shl_coef_(), lmmp_ifft_(), LMMP_MIN, lmmp_mul_fermat_recombine_(), lmmp_mul_fermat_recurse_(), lmmp_zero, lmmp_zero_q_(), fft_memstack::maxdepth, fft_memstack::temp_coef , 以及 fft_memstack::tempdepth.

被这些函数引用 lmmp_mul_fft_() , 以及 lmmp_mullo_fft_().

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

◆ lmmp_mul_fft_()

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

FFT乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
???<=nb<=na, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

在文件 mul_fft.c 第 1085 行定义.

1085 {
1086 lmmp_param_assert(na > 0 && nb > 0);
1087 lmmp_param_assert(na >= nb);
1088 mp_size_t hn = lmmp_fft_next_size_((na + nb + 1) >> 1);
1089 lmmp_assert(na + nb > hn);
1090 mp_ptr tp = ALLOC_TYPE(hn + 1, mp_limb_t);
1091
1092 mp_srcptr amodm = numa;
1093 mp_size_t nam = na;
1094 if (na > hn) {
1095 /*
1096 Z = B^hb - 1
1097 amodm = a mod Z
1098 */
1099 if (lmmp_add_(dst, numa, hn, numa + hn, na - hn))
1100 lmmp_inc(dst);
1101 amodm = dst;
1102 nam = hn;
1103 }
1104 lmmp_mul_mersenne_(dst, hn, amodm, nam, numb, nb);
1105
1106 mp_srcptr amodp = numa;
1107 mp_size_t nap = na;
1108 if (na > hn) {
1109 /*
1110 Z = B^hp - 1
1111 amodp = a mod Z
1112 */
1113 tp[hn] = 0;
1114 if (lmmp_sub_(tp, numa, hn, numa + hn, na - hn))
1115 lmmp_inc(tp);
1116 amodp = tp;
1117 nap = hn + 1;
1118 }
1119 lmmp_mul_fermat_(tp, hn, amodp, nap, numb, nb);
1120
1121 mp_limb_t cy = lmmp_shr1add_nc_(dst, dst, tp, hn, tp[hn]);
1122 cy <<= LIMB_BITS - 1;
1123 dst[hn - 1] += cy;
1124 if (dst[hn - 1] < cy)
1125 lmmp_inc(dst);
1126
1127 if (na + nb == 2 * hn) {
1128 cy = tp[hn] + lmmp_sub_n_(dst + hn, dst, tp, hn);
1129 // cy==1 means [tp,hn+1]!=0, then [dst,hn]!=0
1130 // cy==2 is impossible since [tp,hn+1] is normalized.
1131 // so the following dec won't overflow.
1132 lmmp_dec_1(dst, cy);
1133 } else {
1134 cy = lmmp_sub_n_(dst + hn, dst, tp, na + nb - hn);
1135 cy = tp[hn] + lmmp_sub_nc_(tp + na + nb - hn, dst + na + nb - hn, tp + na + nb - hn, 2 * hn - (na + nb), cy);
1136 cy = lmmp_sub_1_(dst, dst, na + nb, cy);
1137 }
1138 lmmp_free(tp);
1139}
mp_srcptr
const mp_limb_t * mp_srcptr
Definition lmmp.h:216
lmmp_free
void lmmp_free(void *ptr)
内存释放函数（调用lmmp_heap_free_fn）
Definition memory.c:204
lmmp_add_
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:1058
lmmp_shr1add_nc_
mp_limb_t lmmp_shr1add_nc_(mp_ptr dst, mp_srcptr numa, mp_srcptr numb, mp_size_t n, mp_limb_t c)
带进位加法后右移1位 [dst,n] = ([numa,n] + [numb,n] + c) >> 1
Definition shr.c:79
lmmp_sub_
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:1072
lmmp_mul_mersenne_
void lmmp_mul_mersenne_(mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
梅森数模乘法 [dst,rn] = [numa,na]*[numb,nb] mod B^rn-1
Definition mul_fft.c:752
lmmp_mul_fermat_
void lmmp_mul_fermat_(mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
费马数模乘法 [dst,rn+1]=[numa,na]*[numb,nb] mod B^rn+1
Definition mul_fft.c:677
lmmp_fft_next_size_
mp_size_t lmmp_fft_next_size_(mp_size_t n)
计算FFT运算所需的最小规整化长度（向上取整到2^k的倍数）
Definition mul_fft.c:84
ALLOC_TYPE
#define ALLOC_TYPE(n, type)
Definition tmp_alloc.h:112

引用了 ALLOC_TYPE, LIMB_BITS, lmmp_add_(), lmmp_assert, lmmp_dec_1, lmmp_fft_next_size_(), lmmp_free(), lmmp_inc, lmmp_mul_fermat_(), lmmp_mul_mersenne_(), lmmp_param_assert, lmmp_shr1add_nc_(), lmmp_sub_(), lmmp_sub_1_(), lmmp_sub_n_(), lmmp_sub_nc_() , 以及 tp.

被这些函数引用 lmmp_mul_(), lmmp_mul_n_() , 以及 lmmp_sqr_().

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

◆ lmmp_mul_fft_unbalance_()

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

FFT不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
hnFFT模域参数
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
???<=nb<=na, na>=3*nb, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_mersenne_()

void lmmp_mul_mersenne_ ( mp_ptr  dst,
mp_size_t  rn,
mp_srcptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)

梅森数模乘法 [dst,rn] = [numa,na]*[numb,nb] mod B^rn-1

参数
dst输出结果缓冲区，长度至少为 rn
rn模运算的阶数参数，rn = lmmp_fft_next_size_((na + nb + 1) >> 1)
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
0<=[numa,na]<B^rn, 0<=[numb,nb]<B^rn, rn = lmmp_fft_next_size_((na + nb + 1) >> 1)
返回
无返回值，结果存储在dst中，

在文件 mul_fft.c 第 752 行定义.

752 {
753 int nsqr = numa != numb || na != nb; // 判断是否为平方运算
754 mp_size_t N = rn * LIMB_BITS; // 结果总比特数
755 mp_size_t k = lmmp_fft_best_k_(rn); // 最优FFT层数
756 mp_size_t K = ((mp_size_t)1) << k; // FFT块数（2^k）
757 // 断言：N必须是K的整数倍
758 lmmp_debug_assert(!(N & (K - 1)));
759 mp_size_t M = N >> k; // 每个块的比特数
760 mp_size_t n = 2 * M + k; // 扩展系数长度（梅森数比费马数少2）
761
762 // 规整n：必须是LIMB_BITS和K/2的整数倍
763 n = (n + LIMB_BITS - 1) & (-LIMB_BITS);
764 n = (((n - 1) >> (k - 1)) + 1) << (k - 1);
765
766 // 初始化内存栈
767 fft_memstack msr;
768 msr.maxdepth = -1;
769 msr.tempdepth = -1;
770 msr.lenw = n / LIMB_BITS; // 系数长度（机器字）
771 mp_size_t nlen = msr.lenw + 1; // 系数总长度
772
773 msr.temp_coef = (mp_ptr)lmmp_fft_memstack_(&msr, (((nlen + 1) << (k + nsqr)) + nlen) * LIMB_BYTES);
774
775 mp_ptr *pfca = (mp_ptr*)(msr.temp_coef + nlen), *pfcb = pfca;
776 mp_size_t narest = na * LIMB_BITS, nbrest = nb * LIMB_BITS;
777
778 for (mp_size_t i = 0; i < K; ++i) {
779 mp_size_t coeflen;
780 pfca[i] = (mp_ptr)(pfca + K) + i * nlen;
781 if (narest > 0) {
782 coeflen = LMMP_MIN(narest, M);
783 narest -= coeflen;
784 lmmp_fft_extract_coef_(pfca[i], numa, M * i, coeflen, msr.lenw);
785 } else {
786 lmmp_zero(pfca[i], nlen);
787 }
788 }
789 lmmp_fft_(&msr, pfca, k, n >> (k - 1));
790
791 if (nsqr) {
792 pfcb += (nlen + 1) << k;
793 for (mp_size_t i = 0; i < K; ++i) {
794 mp_size_t coeflen;
795 pfcb[i] = (mp_ptr)(pfcb + K) + i * nlen;
796 if (nbrest > 0) {
797 coeflen = LMMP_MIN(nbrest, M);
798 nbrest -= coeflen;
799 lmmp_fft_extract_coef_(pfcb[i], numb, M * i, coeflen, msr.lenw);
800 } else {
801 lmmp_zero(pfcb[i], nlen);
802 }
803 }
804 lmmp_fft_(&msr, pfcb, k, n >> (k - 1));
805 }
806
807 lmmp_mul_fermat_recurse_(&msr, pfca, pfcb, K);
808
809 lmmp_ifft_(&msr, pfca, k, n >> (k - 1));
810
811 mp_size_t rhead = 0, maxc = 0;
812 for (mp_size_t i = 0; i < K; ++i) {
813 lmmp_fft_shr_coef_(&msr, pfca + i, k);
814 mp_ptr nums = pfca[i];
815
816 if (nums[nlen - 1]) {
817 lmmp_dec(nums);
818 lmmp_debug_assert(nums[nlen - 1] == 1);
819 nums[nlen - 1] = 0;
820 }
821
822 mp_size_t roffset = i * M;
823 mp_size_t shl = roffset & (LIMB_BITS - 1);
824 roffset /= LIMB_BITS;
825
826 if (shl)
827 lmmp_shl_(nums, nums, nlen, shl);
828
829 if (i == 0) {
830 lmmp_copy(dst, nums, nlen);
831 rhead = nlen;
832 } else if (roffset + nlen <= rn) {
833 lmmp_add_(dst + roffset, nums, nlen, dst + roffset, rhead - roffset);
834 rhead = roffset + nlen;
835 } else {
836 maxc += lmmp_add_(dst + roffset, nums, rn - roffset, dst + roffset, rhead - roffset);
837 maxc += lmmp_add_(dst, dst, rn, nums + rn - roffset, nlen + roffset - rn);
838 rhead = rn;
839 }
840 }
841
842 if (!lmmp_add_1_(dst, dst, rn, 1 + maxc))
843 lmmp_dec(dst);
844
845 lmmp_fft_memstack_(&msr, 0);
846}
lmmp_fft_shr_coef_
static void lmmp_fft_shr_coef_(fft_memstack *ms, mp_ptr *coef, mp_size_t shr)
对模 2^n+1 的系数执行右移操作 右移shr位 = 左移(2n - shr)位（mod 2^n+1的循环特性）
Definition mul_fft.c:217

引用了 k, fft_memstack::lenw, LIMB_BITS, LIMB_BYTES, lmmp_add_(), lmmp_add_1_(), lmmp_copy, lmmp_debug_assert, lmmp_dec, lmmp_fft_(), lmmp_fft_best_k_(), lmmp_fft_extract_coef_(), lmmp_fft_memstack_(), lmmp_fft_shr_coef_(), lmmp_ifft_(), LMMP_MIN, lmmp_mul_fermat_recurse_(), lmmp_shl_(), lmmp_zero, fft_memstack::maxdepth, fft_memstack::temp_coef , 以及 fft_memstack::tempdepth.

被这些函数引用 binvert_mulhi_(), lmmp_div_mulinv_(), lmmp_invappr_newton_(), lmmp_invsqrt_newton_(), lmmp_mul_fft_() , 以及 lmmp_mullo_fft_().

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

◆ lmmp_mul_n_()

void lmmp_mul_n_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)

等长大数乘法操作 [dst,2*n] = [numa,n] * [numb,n]

警告
n>0, sep(dst,[numa|numb]) 特殊情况: n==1时dst<=numa+1是允许的 n==2时dst<=numa是允许的
参数
dst乘积结果输出指针（需要 2*n 的 limb 长度）
numa第一个乘数指针
numb第二个乘数指针
nlimb长度

在文件 mul.c 第 99 行定义.

99 {
100 if (n < MUL_TOOM22_THRESHOLD)
101 lmmp_mul_basecase_(dst, numa, n, numb, n);
102 else if (n < MUL_TOOM33_THRESHOLD)
103 lmmp_mul_toom22_(dst, numa, n, numb, n);
104 else if (n < MUL_TOOM44_THRESHOLD)
105 lmmp_mul_toom33_(dst, numa, n, numb, n);
106 else if (n < MUL_FFT_THRESHOLD)
107 lmmp_mul_toom44_(dst, numa, n, numb, n);
108 else
109 lmmp_mul_fft_(dst, numa, n, numb, n);
110}
lmmp_mul_toom22_
void lmmp_mul_toom22_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
Toom-22乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
lmmp_mul_toom44_
void lmmp_mul_toom44_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
Toom-44乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
lmmp_mul_basecase_
void lmmp_mul_basecase_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
基础乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
lmmp_mul_fft_
void lmmp_mul_fft_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
FFT乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
Definition mul_fft.c:1085
lmmp_mul_toom33_
void lmmp_mul_toom33_(mp_ptr dst, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
Toom-33乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]
MUL_FFT_THRESHOLD
#define MUL_FFT_THRESHOLD
Definition mparam.h:54
MUL_TOOM22_THRESHOLD
#define MUL_TOOM22_THRESHOLD
Definition mparam.h:46
MUL_TOOM33_THRESHOLD
#define MUL_TOOM33_THRESHOLD
Definition mparam.h:50
MUL_TOOM44_THRESHOLD
#define MUL_TOOM44_THRESHOLD
Definition mparam.h:52

引用了 lmmp_mul_basecase_(), lmmp_mul_fft_(), lmmp_mul_toom22_(), lmmp_mul_toom33_(), lmmp_mul_toom44_(), MUL_FFT_THRESHOLD, MUL_TOOM22_THRESHOLD, MUL_TOOM33_THRESHOLD , 以及 MUL_TOOM44_THRESHOLD.

被这些函数引用 binvert_mulhi_(), lmmp_div_mulinv_(), lmmp_invappr_newton_(), lmmp_invsqrt_newton_(), lmmp_mul_(), lmmp_mul_fermat_recurse_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom33_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom43_(), lmmp_mul_toom44_(), lmmp_mul_toom52_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_() , 以及 lmmp_mullo_dc_().

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

◆ lmmp_mul_toom22_()

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

Toom-22乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
4/5<=nb/na<=1, nb>=5, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_() , 以及 lmmp_mul_n_().

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

◆ lmmp_mul_toom32_()

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

Toom-32乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
5/9<=nb/na<=4/5, nb>=12, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom33_()

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

Toom-33乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
4/5<=nb/na<=1, nb>=26, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_() , 以及 lmmp_mul_n_().

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

◆ lmmp_mul_toom42_()

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

Toom-42乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
1/3<=nb/na<=5/9, nb>=20, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom42_unbalance_()

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

Toom-42不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
na>=3*nb, nb>=20, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom43_()

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

Toom-43乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
3/5<=nb/na<=4/5, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom44_()

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

Toom-44乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
4/5<=nb/na<=1, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_() , 以及 lmmp_mul_n_().

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

◆ lmmp_mul_toom52_()

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

Toom-52乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
1/3<=nb/na<=9/20, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom53_()

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

Toom-53乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
9/20<=nb/na<=3/5, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom62_()

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

Toom-62乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
1/5<=nb/na<=1/3, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mul_toom62_unbalance_()

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

Toom-62不平衡乘法运算 [dst,na+nb] = [numa,na] * [numb,nb]

参数
dst输出结果缓冲区，长度至少为 na+nb
numa第一个输入操作数，长度为 na
na第一个操作数的 limb 长度
numb第二个输入操作数，长度为 nb
nb第二个操作数的 limb 长度
警告
na>=5*nb, nb>=??, sep(dst,[numa|numb])
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_mul_().

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

◆ lmmp_mulh_()

mp_limb_t lmmp_mulh_ ( mp_limb_t  a,
mp_limb_t  b 
)

计算两个64位无符号整数相乘的高位结果 (a*b)/2^64

参数
a第一个64位无符号整数
b第二个64位无符号整数
返回
乘积的高64位结果

在文件 tiny.c 第 73 行定义.

73 {
74#if (defined(__GNUC__) || defined(__clang__)) && defined(__SIZEOF_INT128__)
75 __uint128_t t = (__uint128_t)a * (__uint128_t)b;
76 return (mp_limb_t)(t >> 64);
77#elif defined(_MSC_VER) && (defined(_M_X64) || defined(_M_ARM64))
78 return __umulh(a, b);
79#else
80 uint64_t ah = a >> 32, bh = b >> 32;
81 a = (uint32_t)a, b = (uint32_t)b;
82 uint64_t r0 = a * b, r1 = a * bh, r2 = ah * b, r3 = ah * bh;
83 r3 += (r1 >> 32) + (r2 >> 32);
84 r1 = (uint32_t)r1, r2 = (uint32_t)r2;
85 r1 += r2;
86 r1 += (r0 >> 32);
87 return r3 + (r1 >> 32);
88#endif
89}

引用了 r0, r1, r2 , 以及 r3.

被这些函数引用 lmmp_from_str_len_(), lmmp_to_str_basecase_() , 以及 lmmp_to_str_len_().

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

◆ lmmp_mullh_()

void lmmp_mullh_ ( mp_limb_t  a,
mp_limb_t  b,
mp_ptr  dst 
)

计算两个64位无符号整数相乘的128位结果 (a*b)

参数
dst输出结果缓冲区，存储乘积结果，长度为 2
a第一个64位无符号整数
b第二个64位无符号整数
警告
dst 必须指向一个长度为 2 的数组
返回
无返回值

被这些函数引用 lmmp_u16_pow_1_(), lmmp_u32_pow_1_() , 以及 lmmp_u64_pow_1_().

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

◆ lmmp_mullo_()

void lmmp_mullo_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)

低位乘法 [dst,n] = [numa,n] * [numb,n] mod B^n

参数
dst输出结果缓冲区，长度至少为 n
numa第一个输入操作数，长度为 n
numb第二个输入操作数，长度为 n
nlimb长度
警告
n>0, sep(dst,[numa|numb]) 特殊情况：当 n >= MULLO_DC_THRESHOLD 时，eqsep(dst,[numa|numb])是允许的
返回
无返回值，结果存储在dst中，[dst,n]=[numa,n] * [numb,n] mod B^n

◆ lmmp_mullo_dc_()

void lmmp_mullo_dc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_ptr  tp,
mp_size_t  n 
)

低位乘法 [dst,n] = [numa,n] * [numb,n] mod B^n

参数
dst输出结果缓冲区，长度至少为 n
numa第一个输入操作数，长度为 n
numb第二个输入操作数，长度为 n
tp临时缓冲区，长度至少为 2*n
nlimb长度
警告
n>0, sep(dst,[numa|numb],tp)
返回
无返回值，结果存储在dst中，[dst,n]=[numa,n] * [numb,n] mod B^n

被这些函数引用 lmmp_mullo_n_().

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

◆ lmmp_mullo_fft_()

void lmmp_mullo_fft_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_ptr  scratch 
)

低位FFT乘法 [dst,n] = [numa,n] * [numb,n] mod B^n

参数
dst输出结果缓冲区，长度至少为 n
numa第一个输入操作数，长度为 n
numb第二个输入操作数，长度为 n
scratch临时缓冲区，长度至少为 2*n
n缓冲区 limb 长度
警告
???<n, sep(scratch,[numa|numb]), eqsep(dst,scratch)
返回
无返回值，结果存储在dst中，[dst,n]=[numa,n] * [numb,n] mod B^n

在文件 mullo.c 第 11 行定义.

11 {
12 lmmp_param_assert(n > 0);
13 mp_size_t hn = lmmp_fft_next_size_((n + n + 1) >> 1);
14 lmmp_assert(n + n > hn);
15 mp_ptr tp = ALLOC_TYPE(hn + 1, mp_limb_t);
16
17 mp_srcptr amodm = numa;
18 mp_size_t nam = n;
19 if (n > hn) {
20 /*
21 Z = B^hb - 1
22 amodm = a mod Z
23 */
24 if (lmmp_add_(scratch, numa, hn, numa + hn, n - hn))
25 lmmp_inc(scratch);
26 amodm = scratch;
27 nam = hn;
28 }
29 lmmp_mul_mersenne_(scratch, hn, amodm, nam, numb, n);
30
31 mp_srcptr amodp = numa;
32 mp_size_t nap = n;
33 if (n > hn) {
34 /*
35 Z = B^hp - 1
36 amodp = a mod Z
37 */
38 tp[hn] = 0;
39 if (lmmp_sub_(tp, numa, hn, numa + hn, n - hn))
40 lmmp_inc(tp);
41 amodp = tp;
42 nap = hn + 1;
43 }
44 lmmp_mul_fermat_(tp, hn, amodp, nap, numb, n);
45
46 mp_limb_t cy = lmmp_shr1add_nc_(scratch, scratch, tp, hn, tp[hn]);
47 cy <<= LIMB_BITS - 1;
48 scratch[hn - 1] += cy;
49 if (scratch[hn - 1] < cy)
50 lmmp_inc(scratch);
51
52 if (n == hn) {
53 cy = tp[hn] + lmmp_sub_n_(scratch + hn, scratch, tp, hn);
54 // cy==1 means [tp,hn+1]!=0, then [dst,hn]!=0
55 // cy==2 is impossible since [tp,hn+1] is normalized.
56 // so the following dec won't overflow.
57 lmmp_dec_1(scratch, cy);
58 } else {
59 mp_size_t n2 = 2 * n;
60 cy = lmmp_sub_n_(scratch + hn, scratch, tp, n2 - hn);
61 cy = tp[hn] + lmmp_sub_nc_(tp + n2 - hn, scratch + n2 - hn, tp + n2 - hn, 2 * hn - n2, cy);
62 cy = lmmp_sub_1_(scratch, scratch, n2, cy);
63 }
64 lmmp_free(tp);
65 lmmp_copy(dst, scratch, n);
66}
scratch
#define scratch
lmmp_mul_fermat_
void lmmp_mul_fermat_(mp_ptr dst, mp_size_t rn, mp_srcptr numa, mp_size_t na, mp_srcptr numb, mp_size_t nb)
费马数模乘法 [dst,rn+1]=[numa,na]*[numb,nb] mod B^rn+1
Definition mul_fft.c:677

引用了 ALLOC_TYPE, LIMB_BITS, lmmp_add_(), lmmp_assert, lmmp_copy, lmmp_dec_1, lmmp_fft_next_size_(), lmmp_free(), lmmp_inc, lmmp_mul_fermat_(), lmmp_mul_mersenne_(), lmmp_param_assert, lmmp_shr1add_nc_(), lmmp_sub_(), lmmp_sub_1_(), lmmp_sub_n_(), lmmp_sub_nc_(), scratch , 以及 tp.

被这些函数引用 lmmp_mullo_(), lmmp_mullo_n_() , 以及 lmmp_sqrlo_n_().

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

◆ lmmp_not_()

void lmmp_not_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

大数按位取反操作 [dst,na] = ~[numa,na] (对每个limb执行按位非操作)

参数
dst结果输出指针
numa源操作数指针
nalimb长度
警告
na>0, eqsep(dst,numa)

被这些函数引用 lmmp_binvert_n_dc_(), lmmp_fft_shl_coef_(), lmmp_inv_basecase_(), lmmp_invappr_newton_() , 以及 lmmp_invsqrt_newton_().

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

◆ lmmp_shl_()

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

参数
dst结果输出指针
numa源操作数指针
nalimb长度
shl左移的位数 (0~63)
警告
na>0, 0<=shl<64, eqsep(dst,numa) 允许dst指针地址大于numa（即支持原地长移位操作）
返回
其最低shl个比特位填充[numa,na]被移出的shl个最高位，其余比特位为0

在文件 shl.c 第 9 行定义.

9 {
10 if (shr == 0) {
11 lmmp_copy(dst, numa, na);
12 return 0;
13 } else {
14 const mp_limb_t rshr = LIMB_BITS - shr;
15 mp_limb_t high_limb, low_limb;
16 mp_limb_t retval;
17 numa += na;
18 dst += na;
19 low_limb = *--numa;
20 retval = low_limb >> rshr;
21 high_limb = (low_limb << shr);
22 while (--na != 0) {
23 low_limb = *--numa;
24 *--dst = high_limb | (low_limb >> rshr);
25 high_limb = (low_limb << shr);
26 }
27 *--dst = high_limb;
28 return retval;
29 }
30}

引用了 LIMB_BITS , 以及 lmmp_copy.

被这些函数引用 lmmp_10pow_1_(), lmmp_12pow_1_(), lmmp_14pow_1_(), lmmp_6pow_1_(), lmmp_arith_seqprod_(), lmmp_bninv_(), lmmp_bninv_appr_newton_(), lmmp_div_(), lmmp_factorial_(), lmmp_fft_shl_coef_(), lmmp_inv_(), lmmp_mul_fermat_recombine_(), lmmp_mul_mersenne_(), lmmp_mul_mersenne_single_(), lmmp_mul_toom43_(), lmmp_mul_toom44_(), lmmp_mul_toom53_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_(), lmmp_mul_toom62_cache_init_(), lmmp_multinomial_(), lmmp_nCr_(), lmmp_nPr_(), lmmp_pow_1_(), lmmp_sqr_toom4_(), lmmp_sqrt_(), lmmp_sqrt_newton_(), lmmp_to_str_(), lmmp_to_str_divide_(), lmmp_toom_eval_dgr3_pm2_(), lmmp_toom_eval_pm2_(), lmmp_toom_interp5_(), lmmp_toom_interp6_(), lmmp_toom_interp7_() , 以及 pow_nPr_().

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

◆ lmmp_shl_c_()

mp_limb_t lmmp_shl_c_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  shl,
mp_limb_t  c 
)

带进位的大数左移操作 [dst,na] = [numa,na]<<shl，dst的低shl位填充c的低shl位

参数
dst结果输出指针
numa源操作数指针
nalimb长度
shl左移的位数 (0~63)
c进位值（其高(64-shl)位必须为0）
警告
na>0, 0<=shl<64, eqsep(dst,numa) c的高(64-shl)位必须为0 允许dst指针地址大于numa（即支持原地长移位操作）
返回
其最低shl个比特位填充[numa,na]被移出的shl个最高位，其余比特位为0

在文件 shl.c 第 32 行定义.

32 {
33 if (shr == 0) {
34 lmmp_copy(dst, numa, na);
35 return 0;
36 } else {
37 const mp_limb_t rshr = LIMB_BITS - shr;
38 mp_limb_t high_limb, low_limb;
39 mp_limb_t retval;
40 numa += na;
41 dst += na;
42 low_limb = *--numa;
43 retval = low_limb >> rshr;
44 high_limb = (low_limb << shr);
45 while (--na != 0) {
46 low_limb = *--numa;
47 *--dst = high_limb | (low_limb >> rshr);
48 high_limb = (low_limb << shr);
49 }
50 c &= ((mp_limb_t)1 << shr) - 1;
51 *--dst = high_limb | c;
52 return retval;
53 }
54}

引用了 LIMB_BITS , 以及 lmmp_copy.

被这些函数引用 lmmp_fft_bfy_().

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

◆ lmmp_shlnot_()

mp_limb_t lmmp_shlnot_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  shl 
)

左移后按位取反操作 [dst,na] = ~([numa,na] << shl)，dst的低shl位填充1

参数
dst结果输出指针
numa源操作数指针
nalimb长度
shl左移的位数 (0~63)
警告
na>0, 0<=shl<64, eqsep(dst,numa)
返回
其最低shl个比特位填充[numa,na]被移出的shl个最高位，其余比特位为0

被这些函数引用 lmmp_fft_shl_coef_() , 以及 lmmp_invsqrt_newton_().

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

◆ lmmp_shr1add_n_()

mp_limb_t lmmp_shr1add_n_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)

加法后右移1位 [dst,n] = ([numa,n] + [numb,n]) >> 1

参数
dst结果输出指针
numa第一个加数指针
numb第二个加数指针
nlimb长度
警告
n>0, eqsep(dst,[numa|numb])
返回
右移操作产生的进位值 [0|1]

在文件 shr.c 第 52 行定义.

52 {
53 mp_size_t i = 0, c = 0, l = 0;
54 mp_limb_t a, b, r;
55
56 a = numa[i];
57 b = numb[i];
58 r = a + c;
59 c = (r < c);
60 r += b;
61 c += (r < b);
62 dst[i] = r >> 1;
63 l = r & 1;
64
65 for (i = 1; i < n; i++) {
66 a = numa[i];
67 b = numb[i];
68 r = a + c;
69 c = (r < c);
70 r += b;
71 c += (r < b);
72 dst[i - 1] |= r << (LIMB_BITS - 1);
73 dst[i] = r >> 1;
74 }
75 dst[n - 1] |= c << (LIMB_BITS - 1);
76 return l;
77}

引用了 LIMB_BITS.

被这些函数引用 lmmp_mul_toom32_(), lmmp_toom_interp5_(), lmmp_toom_interp6_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_shr1add_nc_()

mp_limb_t lmmp_shr1add_nc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_limb_t  c 
)

带进位加法后右移1位 [dst,n] = ([numa,n] + [numb,n] + c) >> 1

参数
dst结果输出指针
numa第一个加数指针
numb第二个加数指针
nlimb长度
c初始进位值 [0|1]
警告
n>0, c=[0|1], eqsep(dst,[numa|numb])
返回
右移操作产生的进位值 [0|1]

在文件 shr.c 第 79 行定义.

79 {
80 mp_size_t i = 0, l = 0;
81 mp_limb_t a, b, r;
82
83 a = numa[i];
84 b = numb[i];
85 r = a + c;
86 c = (r < c);
87 r += b;
88 c += (r < b);
89 dst[i] = r >> 1;
90 l = r & 1;
91
92 for (i = 1; i < n; i++) {
93 a = numa[i];
94 b = numb[i];
95 r = a + c;
96 c = (r < c);
97 r += b;
98 c += (r < b);
99 dst[i - 1] |= r << (LIMB_BITS - 1);
100 dst[i] = r >> 1;
101 }
102 dst[n - 1] |= c << (LIMB_BITS - 1);
103 return l;
104}

引用了 LIMB_BITS.

被这些函数引用 lmmp_mul_fft_(), lmmp_mul_fft_cache_() , 以及 lmmp_mullo_fft_().

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

◆ lmmp_shr1sub_n_()

mp_limb_t lmmp_shr1sub_n_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)

减法后右移1位 [dst,n] = ([numa,n] - [numb,n]) >> 1

参数
dst结果输出指针
numa被减数指针
numb减数指针
n操作数的位数（limb数量）
警告
n>0, eqsep(dst,[numa|numb])
返回
右移操作产生的进位值 (0或1)

在文件 shr.c 第 106 行定义.

106 {
107 mp_size_t i = 0, c = 0, l = 0;
108 mp_limb_t a, b, r;
109
110 a = numa[i];
111 b = numb[i];
112 b += c;
113 c = (b < c);
114 c += (a < b);
115 r = a - b;
116 dst[i] = r >> 1;
117 l = r & 1;
118
119 for (i = 1; i < n; i++) {
120 a = numa[i];
121 b = numb[i];
122 b += c;
123 c = (b < c);
124 c += (a < b);
125 r = a - b;
126 dst[i - 1] |= r << (LIMB_BITS - 1);
127 dst[i] = r >> 1;
128 }
129 dst[n - 1] |= c << (LIMB_BITS - 1);
130 return l;
131}

引用了 LIMB_BITS.

被这些函数引用 lmmp_mul_toom32_(), lmmp_toom_interp5_(), lmmp_toom_interp6_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_shr1sub_nc_()

mp_limb_t lmmp_shr1sub_nc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_limb_t  c 
)

带借位减法后右移1位 [dst,n] = ([numa,n] - [numb,n] - c) >> 1

参数
dst结果输出指针
numa被减数指针
numb减数指针
nlimb长度
c初始借位值 [0|1]
警告
n>0, c=[0|1], eqsep(dst,[numa|numb])
返回
右移操作产生的进位值 [0|1]

在文件 shr.c 第 133 行定义.

133 {
134 mp_size_t i = 0, l = 0;
135 mp_limb_t a, b, r;
136
137 a = numa[i];
138 b = numb[i];
139 b += c;
140 c = (b < c);
141 c += (a < b);
142 r = a - b;
143 dst[i] = r >> 1;
144 l = r & 1;
145
146 for (i = 1; i < n; i++) {
147 a = numa[i];
148 b = numb[i];
149 b += c;
150 c = (b < c);
151 c += (a < b);
152 r = a - b;
153 dst[i - 1] |= r << (LIMB_BITS - 1);
154 dst[i] = r >> 1;
155 }
156 dst[n - 1] |= c << (LIMB_BITS - 1);
157 return l;
158}

引用了 LIMB_BITS.

◆ lmmp_shr_()

mp_limb_t lmmp_shr_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  shr 
)

大数右移操作 [dst,na] = [numa,na] >> shr，dst的高shr位填充0

参数
dst结果输出指针
numa源操作数指针
nalimb长度
shr右移的位数 (0~63)
警告
na>0, 0<=shr<64, eqsep(dst,numa) 允许dst指针地址小于numa（即支持原地长移位操作）
返回
其最高shr个比特位填充[numa,na]被移出的shr个最低位，其余比特位为0

在文件 shr.c 第 9 行定义.

9 {
10 if (shr == 0) {
11 lmmp_copy(dst, numa, na);
12 return 0;
13 } else {
14 mp_limb_t high_limb, low_limb;
15 const mp_size_t rshr = LIMB_BITS - shr;
16 mp_limb_t retval;
17 high_limb = *numa++;
18 retval = (high_limb << rshr);
19 low_limb = high_limb >> shr;
20 while (--na != 0) {
21 high_limb = *numa++;
22 *dst++ = low_limb | (high_limb << rshr);
23 low_limb = high_limb >> shr;
24 }
25 *dst = low_limb;
26 return retval;
27 }
28}

引用了 LIMB_BITS , 以及 lmmp_copy.

被这些函数引用 lmmp_bninv_(), lmmp_div_(), lmmp_fft_extract_coef_(), lmmp_invsqrt_newton_(), lmmp_sqrt_(), lmmp_sqrt_newton_(), lmmp_to_str_divide_(), lmmp_toom_interp6_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_shr_c_()

mp_limb_t lmmp_shr_c_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  shr,
mp_limb_t  c 
)

带进位的大数右移操作 [dst,na] = [numa,na]>>shr，dst的高shr位填充c的高shr位

参数
dst结果输出指针
numa源操作数指针
nalimb长度
shr右移的位数 (0~63)
c进位值（其低(64-shr)位必须为0）
警告
na>0, 0<=shr<64, eqsep(dst,numa) c的低(64-shr)位必须为0 允许dst指针地址小于numa（即支持原地长移位操作）
返回
其最高shr个比特位填充[numa,na]被移出的shr个最低位，其余比特位为0

在文件 shr.c 第 30 行定义.

30 {
31 if (shr == 0) {
32 lmmp_copy(dst, numa, na);
33 return 0;
34 } else {
35 mp_limb_t high_limb, low_limb;
36 const mp_size_t rshr = LIMB_BITS - shr;
37 mp_limb_t retval;
38 high_limb = *numa++;
39 retval = (high_limb << rshr);
40 low_limb = high_limb >> shr;
41 while (--na != 0) {
42 high_limb = *numa++;
43 *dst++ = low_limb | (high_limb << rshr);
44 low_limb = high_limb >> shr;
45 }
46 c &= ~(((mp_limb_t)1 << rshr) - 1);
47 *dst = low_limb | c;
48 return retval;
49 }
50}

引用了 LIMB_BITS , 以及 lmmp_copy.

被这些函数引用 lmmp_ifft_bfy_(), lmmp_inv_() , 以及 lmmp_sqrt_divide_().

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

◆ lmmp_sqr_()

void lmmp_sqr_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

大数平方操作 [dst,2*na] = [numa,na]^2

警告
na>0, sep(dst,numa)
参数
dst平方结果输出指针（需要2*na的limb长度）
numa源操作数指针
nalimb长度

在文件 sqr.c 第 10 行定义.

10 {
11 if (na < MUL_TOOM22_THRESHOLD)
12 lmmp_sqr_basecase_(dst, numa, na);
13 else if (na < MUL_TOOM33_THRESHOLD)
14 lmmp_sqr_toom2_(dst, numa, na);
15 else if (na < MUL_TOOM44_THRESHOLD)
16 lmmp_sqr_toom3_(dst, numa, na);
17 else if (na < MUL_FFT_THRESHOLD)
18 lmmp_sqr_toom4_(dst, numa, na);
19 else
20 lmmp_mul_fft_(dst, numa, na, numa, na);
21}
lmmp_sqr_basecase_
void lmmp_sqr_basecase_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
基础平方运算 [dst,2*na] = [numa,na]^2
lmmp_sqr_toom2_
void lmmp_sqr_toom2_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
Toom-2平方运算 [dst,2*na] = [numa,na]^2
lmmp_sqr_toom3_
void lmmp_sqr_toom3_(mp_ptr dst, mp_srcptr numa, mp_size_t na)
Toom-3平方运算 [dst,2*na] = [numa,na]^2
lmmp_sqr_toom4_
void lmmp_sqr_toom4_(mp_ptr pp, mp_srcptr ap, mp_size_t an)
Toom-4平方运算 [dst,2*na] = [numa,na]^2

引用了 lmmp_mul_fft_(), lmmp_sqr_basecase_(), lmmp_sqr_toom2_(), lmmp_sqr_toom3_(), lmmp_sqr_toom4_(), MUL_FFT_THRESHOLD, MUL_TOOM22_THRESHOLD, MUL_TOOM33_THRESHOLD , 以及 MUL_TOOM44_THRESHOLD.

被这些函数引用 lmmp_3pow_1_(), lmmp_9pow_1_(), lmmp_bninv_appr_newton_(), lmmp_factors_mul_(), lmmp_from_str_(), lmmp_invsqrt_newton_(), lmmp_mul_(), lmmp_mul_fermat_recurse_(), lmmp_pow_(), lmmp_pow_basecase_(), lmmp_pow_win2_(), lmmp_remove_(), lmmp_sqr_signed_(), lmmp_sqr_toom2_(), lmmp_sqr_toom3_(), lmmp_sqr_toom4_(), lmmp_sqrlo_dc_(), lmmp_sqrt_divide_(), lmmp_to_str_(), lmmp_u16_pow_1_(), lmmp_u32_pow_1_(), lmmp_u64_pow_1_() , 以及 lmmp_u8_pow_1_().

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

◆ lmmp_sqr_basecase_()

void lmmp_sqr_basecase_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

基础平方运算 [dst,2*na] = [numa,na]^2

参数
dst输出结果缓冲区，长度至少为2*na
numa输入操作数，长度为na
na输入操作数的 limb 长度
警告
0<na, sep(dst,numa)
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_sqr_(), lmmp_u32_pow_1_() , 以及 lmmp_u64_pow_1_().

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

◆ lmmp_sqr_toom2_()

void lmmp_sqr_toom2_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

Toom-2平方运算 [dst,2*na] = [numa,na]^2

参数
dst输出结果缓冲区，长度至少为 2*na
numa输入操作数，长度为na
na输入操作数的 limb 长度
警告
??<na, sep(dst,numa)
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_sqr_().

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

◆ lmmp_sqr_toom3_()

void lmmp_sqr_toom3_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na 
)

Toom-3平方运算 [dst,2*na] = [numa,na]^2

参数
dst输出结果缓冲区，长度至少为2*na
numa输入操作数，长度为na
na输入操作数的单精度数(limb)长度
警告
??<na, sep(dst,numa)
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_sqr_().

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

◆ lmmp_sqr_toom4_()

void lmmp_sqr_toom4_ ( mp_ptr  pp,
mp_srcptr  ap,
mp_size_t  an 
)

Toom-4平方运算 [dst,2*na] = [numa,na]^2

参数
dst输出结果缓冲区，长度至少为2*na
numa输入操作数，长度为na
na输入操作数的单精度数(limb)长度
警告
??<na, sep(dst,numa)
返回
无返回值，结果存储在dst中

被这些函数引用 lmmp_sqr_().

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

◆ lmmp_sqrlo_dc_()

void lmmp_sqrlo_dc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_ptr  tp,
mp_size_t  n 
)

低位平方 [dst,n] = [numa,n]^2 mod B^n

参数
dst输出结果缓冲区，长度至少为 n
numa第一个输入操作数，长度为 n
tp临时缓冲区，长度至少为 2*n
nlimb长度
警告
n>0, sep(dst,numa,tp)
返回
无返回值，结果存储在dst中，[dst,n]=[numa,n]^2 mod B^n

被这些函数引用 lmmp_sqrlo_n_().

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

◆ lmmp_sqrt_()

void lmmp_sqrt_ ( mp_ptr  dsts,
mp_ptr  dstr,
mp_srcptr  numa,
mp_size_t  na,
mp_size_t  nf 
)

大数平方根和取余操作

注解
如果dstr不为NULL: [dsts,nf+na/2+1], [dstr,nf+na/2+1] = sqrtrem([numa,na]*B^(2*nf)) 也即 [numa,na] × B^(2×nf) = [dsts,nf+na/2+1]^2 + [dstr,nf+na/2+1] 且 0 <= [dstr,nf+na/2+1] < 2 * [dsts,nf+na/2+1] + 1 如果dstr为NULL: [dsts,nf+na/2+1] = [round|floor](sqrt([numa,na]*B^(2*nf)))
警告
na>0, numa[na-1]!=0, eqsep(dsts,numa), eqsep(dstr,numa)
参数
dsts平方根结果输出指针
dstr余数结果输出指针（NULL表示不计算余数）
numa源操作数指针
na操作数的 limb 长度
nf精度因子

在文件 sqrt.c 第 323 行定义.

323 {
324 lmmp_debug_assert(na > 0);
325 lmmp_debug_assert(numa[na - 1] > 0);
326 mp_limb_t high = numa[na - 1];
327 int nsh = lmmp_leading_zeros_(high) / 2;
328 mp_size_t nl = na + 2 * nf;
329 if (nl == 1) {
330 mp_limb_t srt;
331 lmmp_sqrt_1_(&srt, high << nsh * 2);
332 srt >>= nsh;
333 dsts[0] = srt;
334 if (dstr)
335 dstr[0] = high - srt * srt;
336 } else if (!dstr && nf >= 10 * na + SQRT_NEWTON_THRESHOLD) {
337 lmmp_sqrt_newton_(dsts, numa, na, nf);
338 } else {
339 TEMP_DECL;
340 mp_limb_t ns = (nl + 1) / 2;
341 mp_ptr numa2 = TALLOC_TYPE(2 * ns, mp_limb_t);
342 if (nf)
343 lmmp_zero(numa2, 2 * nf);
344 if (nsh)
345 lmmp_shl_(numa2 + 2 * ns - na, numa, na, nsh * 2);
346 else
347 lmmp_copy(numa2 + 2 * ns - na, numa, na);
348 if (nl & 1) {
349 numa2[2 * nf] = 0;
350 nsh += LIMB_BITS / 2;
351 } else {
352 dsts[ns] = 0;
353 }
354 mp_limb_t rh = lmmp_sqrt_divide_(dsts, numa2, ns, dstr ? 0 : nsh);
355 if (nsh) {
356 if (dstr) {
357 mp_limb_t ds = dsts[0] & (((mp_limb_t)1 << nsh) - 1);
358 rh += lmmp_addmul_1_(numa2, dsts, ns, 2 * ds);
359 mp_limb_t b = lmmp_submul_1_(numa2, &ds, 1, ds);
360 if (ns == 1)
361 rh -= b;
362 else
363 rh -= lmmp_sub_1_(numa2 + 1, numa2 + 1, ns - 1, b);
364 }
365 lmmp_shr_(dsts, dsts, ns, nsh);
366 }
367 if (dstr) {
368 numa2[ns] = rh;
369 nsh *= 2;
370 if (nsh >= LIMB_BITS) {
371 nsh -= LIMB_BITS;
372 ++numa2;
373 } else
374 ++ns;
375 if (nsh)
376 lmmp_shr_(dstr, numa2, ns, nsh);
377 else
378 lmmp_copy(dstr, numa2, ns);
379 }
380
381 TEMP_FREE;
382 }
383}
lmmp_addmul_1_
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
SQRT_NEWTON_THRESHOLD
#define SQRT_NEWTON_THRESHOLD
Definition mparam.h:41
lmmp_sqrt_divide_
static mp_limb_t lmmp_sqrt_divide_(mp_ptr dsts, mp_ptr numa, mp_size_t ns, int nsh)
Definition sqrt.c:107
lmmp_sqrt_1_
static mp_limb_t lmmp_sqrt_1_(mp_ptr dsts, mp_limb_t x)
Definition sqrt.c:38
lmmp_sqrt_newton_
static void lmmp_sqrt_newton_(mp_ptr dsts, mp_srcptr numa, mp_size_t na, mp_size_t nf)
Definition sqrt.c:275

引用了 LIMB_BITS, lmmp_addmul_1_(), lmmp_copy, lmmp_debug_assert, lmmp_leading_zeros_(), lmmp_shl_(), lmmp_shr_(), lmmp_sqrt_1_(), lmmp_sqrt_divide_(), lmmp_sqrt_newton_(), lmmp_sub_1_(), lmmp_submul_1_(), lmmp_zero, SQRT_NEWTON_THRESHOLD, TALLOC_TYPE, TEMP_DECL , 以及 TEMP_FREE.

+ 函数调用图:

◆ lmmp_sub_()

static mp_limb_t lmmp_sub_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_srcptr  numb,
mp_size_t  nb 
)
inlinestatic

大数减法静态内联函数 [dst,na]=[numa,na]-[numb,nb]

参数
dst输出结果缓冲区，存储numa - numb
numa被减数，长度为na
na被减数的 limb 长度
numb减数，长度为nb
nb减数的 limb 长度
返回
借位标志（1表示有借位，0表示无借位）
警告
0<nb<=na, eqsep(dst,[numa|numb])

在文件 lmmpn.h 第 1072 行定义.

1072 {
1073 LMMP_AORS_(lmmp_sub_n_, ((dst[nb++] = _x_ - 1), _x_ == 0));
1074}

引用了 LMMP_AORS_ , 以及 lmmp_sub_n_().

被这些函数引用 lmmp_add_signed_(), lmmp_lehmer_mul_(), lmmp_mul_fermat_recombine_(), lmmp_mul_fft_(), lmmp_mul_fft_cache_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom52_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_(), lmmp_mullo_fft_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_sub_1_()

static mp_limb_t lmmp_sub_1_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_size_t  na,
mp_limb_t  x 
)
inlinestatic

大数减单精度数静态内联函数 [dst,na]=[numa,na]-x

参数
dst输出结果缓冲区，存储numa - x
numa被减数，长度为na
na被减数的 limb 长度
x减数（单个 limb ）
返回
借位标志（1表示有借位，0表示无借位）
警告
na>0, eqsep(dst,numa)

在文件 lmmpn.h 第 1122 行定义.

1122{ LMMP_AORS_1_(-, LMMP_SUBCB_); }
LMMP_SUBCB_
#define LMMP_SUBCB_(r, x, y)
Definition lmmpn.h:1100

引用了 LMMP_AORS_1_ , 以及 LMMP_SUBCB_.

被这些函数引用 lmmp_div_divide_n_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_mul_fft_(), lmmp_mul_fft_cache_(), lmmp_mul_toom32_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_init_(), lmmp_mullo_fft_(), lmmp_sqrt_() , 以及 lmmp_sqrt_divide_().

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

◆ lmmp_sub_n_()

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]

参数
dst结果输出指针
numa被减数指针
numb减数指针
nlimb长度
警告
n>0, eqsep(dst,[numa|numb])
返回
运算后的最终借位值 [0|1]

在文件 sub_n.c 第 70 行定义.

70 {
71 mp_size_t i = 0;
72 mp_limb_t cy = 0;
73
74 for (; i + 4 <= n; i += 4) {
75 mp_limb_t a0, b0, r0;
76 mp_limb_t a1, b1, r1;
77 mp_limb_t a2, b2, r2;
78 mp_limb_t a3, b3, r3;
79
80 a0 = numa[i + 0];
81 b0 = numb[i + 0];
82
83 a1 = numa[i + 1];
84 b1 = numb[i + 1];
85
86 a2 = numa[i + 2];
87 b2 = numb[i + 2];
88
89 a3 = numa[i + 3];
90 b3 = numb[i + 3];
91
92 b0 += cy;
93 cy = (b0 < cy);
94 cy += (a0 < b0);
95 r0 = a0 - b0;
96
97 b1 += cy;
98 cy = (b1 < cy);
99 cy += (a1 < b1);
100 r1 = a1 - b1;
101
102 b2 += cy;
103 cy = (b2 < cy);
104 cy += (a2 < b2);
105 r2 = a2 - b2;
106
107 b3 += cy;
108 cy = (b3 < cy);
109 cy += (a3 < b3);
110 r3 = a3 - b3;
111
112 dst[i + 0] = r0;
113 dst[i + 1] = r1;
114 dst[i + 2] = r2;
115 dst[i + 3] = r3;
116 }
117
118 for (; i < n; i++) {
119 mp_limb_t a, b;
120 a = numa[i];
121 b = numb[i];
122 b += cy;
123 cy = (b < cy);
124 cy += (a < b);
125 dst[i] = a - b;
126 }
127
128 return cy;
129}

引用了 a0, a1, a2, a3, b0, b1, b2, b3, r0, r1, r2 , 以及 r3.

被这些函数引用 lmmp_bninv_appr_newton_(), lmmp_div_(), lmmp_div_basecase_(), lmmp_div_divide_n_(), lmmp_div_mulinv_(), lmmp_div_s_(), lmmp_invappr_newton_(), lmmp_mul_fermat_recurse_(), lmmp_mul_fft_(), lmmp_mul_fft_cache_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom52_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_(), lmmp_mullo_fft_(), lmmp_sqr_toom2_(), lmmp_sqrt_divide_(), lmmp_sub_(), lmmp_toom_interp5_(), lmmp_toom_interp6_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_sub_nc_()

mp_limb_t lmmp_sub_nc_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_limb_t  c 
)

带借位的n位减法 [dst,n] = [numa,n] - [numb,n] - c

参数
dst结果输出指针
numa被减数指针
numb减数指针
nlimb长度
c初始借位值 [0|1]
警告
c=[0|1], n>0, eqsep(dst,[numa|numb])
返回
运算后的最终借位值 [0|1]

在文件 sub_n.c 第 9 行定义.

9 {
10 mp_size_t i = 0;
11 mp_limb_t cy = c;
12
13 for (; i + 4 <= n; i += 4) {
14 mp_limb_t a0, b0, r0;
15 mp_limb_t a1, b1, r1;
16 mp_limb_t a2, b2, r2;
17 mp_limb_t a3, b3, r3;
18
19 a0 = numa[i + 0];
20 b0 = numb[i + 0];
21
22 a1 = numa[i + 1];
23 b1 = numb[i + 1];
24
25 a2 = numa[i + 2];
26 b2 = numb[i + 2];
27
28 a3 = numa[i + 3];
29 b3 = numb[i + 3];
30
31 b0 += cy;
32 cy = (b0 < cy);
33 cy += (a0 < b0);
34 r0 = a0 - b0;
35
36 b1 += cy;
37 cy = (b1 < cy);
38 cy += (a1 < b1);
39 r1 = a1 - b1;
40
41 b2 += cy;
42 cy = (b2 < cy);
43 cy += (a2 < b2);
44 r2 = a2 - b2;
45
46 b3 += cy;
47 cy = (b3 < cy);
48 cy += (a3 < b3);
49 r3 = a3 - b3;
50
51 dst[i + 0] = r0;
52 dst[i + 1] = r1;
53 dst[i + 2] = r2;
54 dst[i + 3] = r3;
55 }
56
57 for (; i < n; i++) {
58 mp_limb_t a, b;
59 a = numa[i];
60 b = numb[i];
61 b += cy;
62 cy = (b < cy);
63 cy += (a < b);
64 dst[i] = a - b;
65 }
66
67 return cy;
68}

引用了 a0, a1, a2, a3, b0, b1, b2, b3, r0, r1, r2 , 以及 r3.

被这些函数引用 lmmp_add_n_sub_n_(), lmmp_div_(), lmmp_div_mulinv_(), lmmp_fft_bfy_(), lmmp_ifft_bfy_(), lmmp_invappr_newton_(), lmmp_mul_fft_(), lmmp_mul_fft_cache_() , 以及 lmmp_mullo_fft_().

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

◆ lmmp_submul_1_()

mp_limb_t lmmp_submul_1_ ( mp_ptr  numa,
mp_srcptr  numb,
mp_size_t  n,
mp_limb_t  b 
)

大数乘以单limb并累减操作 [numa,n] -= [numb,n] * b

参数
numa被减数指针（结果也存储在此）
numb乘数指针
nlimb长度
b乘数
警告
n>0, eqsep(numa,numb))
返回
运算后的借位limb值

被这些函数引用 lmmp_div_(), lmmp_div_basecase_(), lmmp_sqrt_() , 以及 lmmp_toom_interp7_().

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

◆ lmmp_subshl1_n_()

mp_limb_t lmmp_subshl1_n_ ( mp_ptr  dst,
mp_srcptr  numa,
mp_srcptr  numb,
mp_size_t  n 
)

减法结合左移1位操作 [dst,n] = [numa,n] - ([numb,n] << 1)

参数
dst结果输出指针
numa被减数指针
numb减数指针（先左移1位）
nlimb长度
警告
n>0, eqsep(dst,[numa|numb])
返回
运算后的借位值 [0|1|2]

在文件 shl.c 第 73 行定义.

73 {
74 mp_size_t i, c = 0, mb = 0;
75
76 for (i = 0; i < n; i++) {
77 mp_limb_t a, b;
78 a = numa[i];
79 b = (numb[i] << 1) + mb;
80 mb = numb[i] >> (LIMB_BITS - 1);
81 b += c;
82 c = (b < c);
83 c += (a < b);
84 dst[i] = a - b;
85 }
86 return c + mb;
87}

引用了 LIMB_BITS.

◆ lmmp_tailing_zeros_()

int lmmp_tailing_zeros_ ( mp_limb_t  x)

计算一个单精度数(limb)中末尾零的个数

参数
x输入的64位无符号整数
返回
末尾零的位数（范围：0~64）

在文件 tiny.c 第 54 行定义.

54 {
55 if (x == 0) return 64;
56#ifdef __GNUC__
57 return __builtin_ctzll(x);
58#elif _MSC_VER
59 unsigned long index;
60 _BitScanForward64(&index, x);
61 return (int)index;
62#else
63 int n = 0;
64 if ((x & 0x00000000FFFFFFFF) == 0) { n += 32; x >>= 32; }
65 if ((x & 0x000000000000FFFF) == 0) { n += 16; x >>= 16; }
66 if ((x & 0x00000000000000FF) == 0) { n += 8; x >>= 8; }
67 if ((x & 0x000000000000000F) == 0) { n += 4; x >>= 4; }
68 if ((x & 0x0000000000000003) == 0) { n += 2; x >>= 2; }
69 if ((x & 0x0000000000000001) == 0) { n += 1; x >>= 1; }
70#endif
71}

被这些函数引用 lmmp_from_str_(), lmmp_gcd_11_(), lmmp_gcd_22_(), lmmp_pow_(), lmmp_pow_1_(), lmmp_to_str_() , 以及 pow_nPr_().

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

◆ lmmp_to_str_()

mp_size_t lmmp_to_str_ ( mp_byte_t *  dst,
mp_srcptr  numa,
mp_size_t  na,
int  base 
)

大数转字符串操作 [numa,na,B] to [dst,return value,base]

警告
na>=0, 2<=base<=256
参数
dst字符串结果输出指针
numa大数源指针
na大数的 limb 长度
base目标字符串的进制基数
返回
转换后的字符串长度

在文件 to_str.c 第 147 行定义.

147 {
148 do {
149 if (na == 0)
150 return 0;
151 } while (numa[--na] == 0);
152 ++na;
153
154 mp_size_t digits;
155 if (LMMP_POW2_Q(base)) {
156 mp_limb_t curlimb = numa[na - 1];
157 int cnt = lmmp_bases_table[base - 2].large_base;
158 int bitsh = lmmp_limb_bits_(curlimb);
159 int mask = (1 << cnt) - 1;
160 mp_size_t bits = bitsh + LIMB_BITS * (na - 1);
161 digits = (bits - 1) / cnt + 1;
162 dst += digits;
163 int bitpos = digits * cnt - LIMB_BITS * (na - 1);
164
165 do {
166 while ((bitpos -= cnt) >= 0) {
167 *--dst = curlimb >> bitpos & mask;
168 }
169 if (--na == 0)
170 break;
171 mp_limb_t prevlimb = curlimb << -bitpos;
172 curlimb = numa[na - 1];
173 bitpos += LIMB_BITS;
174 *--dst = (prevlimb | curlimb >> bitpos) & mask;
175 } while (1);
176 } else if (na < TO_STR_BASEPOW_THRESHOLD) {
177 digits = lmmp_to_str_basecase_(dst, numa, na, base);
178 } else {
179 TEMP_DECL;
180 mp_basepow_t powers[LIMB_BITS];
181 mp_size_t exps[LIMB_BITS];
182 mp_limb_t lbase = lmmp_bases_table[base - 2].large_base;
183 mp_size_t digitspl = lmmp_bases_table[base - 2].digits_in_limb;
184 mp_size_t bexp = (lmmp_to_str_len_(numa, na, base) - 1) / digitspl + 1;
185 mp_size_t tzbit = lmmp_tailing_zeros_(lbase);
186 // numa 的拷贝空间，多一个 limb 预留规整化移位所需
187 mp_size_t alloc_size = na + 1;
188 mp_limb_t cy;
189 mp_ptr tp;
190
191 int cpow = 0;
192
193 do {
194 bexp = (bexp + 1) >> 1;
195 exps[cpow] = bexp;
196 ++cpow;
197
198 // we will calculate lbase^(bexp-1) first, and trim it s. t.
199 // it contains at most 2 tailing 0 limb, then multiply it by lbase,
200 // so we need npow limbs to store lbase^bexp
201 mp_size_t npow = lmmp_from_str_len_(0, (bexp - 1) * digitspl + 1, base) + 1;
202
203 // space needed for quotients in recursive calls,
204 // quotients are smaller than lbase^bexp
205 alloc_size += npow + 1;
206
207 if (tzbit) {
208 mp_size_t tzlimb = tzbit * (bexp - 1) / LIMB_BITS;
209 if (tzlimb >= 2)
210 npow -= tzlimb - 2;
211 }
212
213 // space needed for a trimmed npow-limb lbase^bexp and its inverse
214 alloc_size += npow * 2;
215 } while (bexp > 1);
216
217 tp = BALLOC_TYPE(alloc_size, mp_limb_t);
218
219 for (int i = 0; i < 2; ++i) {
220 tp[0] = lbase;
221 powers[i].p = tp;
222 powers[i].np = 1;
223 tp += i + 1;
224 powers[i].zeros = 0;
225 powers[i].digits = digitspl * (i + 1);
226 powers[i].base = base;
227 }
228
229 mp_ptr p = powers[1].p;
230 mp_size_t zeros = 0, np = 1;
231 bexp = 1;
232 for (int i = 2; i < cpow; ++i) {
233 lmmp_sqr_(tp, p, np);
234 bexp *= 2;
235 np *= 2;
236 np -= tp[np - 1] == 0;
237 if (bexp + 1 < exps[cpow - 1 - i]) {
238 cy = lmmp_mul_1_(tp, tp, np, lbase);
239 tp[np] = cy;
240 np += cy != 0;
241 ++bexp;
242 }
243 zeros *= 2;
244 while (tp[0] == 0) {
245 // at most 2 tailing 0 limb here
246 ++zeros;
247 ++tp;
248 --np;
249 }
250 p = tp;
251 powers[i].p = p;
252 powers[i].np = np;
253 powers[i].zeros = zeros;
254 powers[i].digits = digitspl * (bexp + 1);
255 powers[i].base = base;
256 tp += np + 1;
257 }
258
259 for (int i = 1; i < cpow; ++i) {
260 p = powers[i].p;
261 np = powers[i].np;
262 cy = lmmp_mul_1_(p, p, np, lbase);
263 p[np] = cy;
264 np += cy != 0;
265 if (p[0] == 0) {
266 ++powers[i].zeros;
267 ++p;
268 --np;
269 }
270
271 powers[i].p = p;
272 powers[i].np = np;
273
274 // Note: all powers except powers[0] are normalized
275 // ASSERT: powers[0] will be never used in lmmp_to_str_divide_
276 // i.e. TO_STR_DIVIDE_THRESHOLD >= 3
277 int cnt = lmmp_leading_zeros_(p[np - 1]);
278 if (powers[i].norm_cnt = cnt)
279 lmmp_shl_(p, p, np, cnt);
280
281 if (np < DIV_MULINV_L_THRESHOLD) {
282 // use divs, no need to compute inv
283 powers[i].invp = 0;
284 powers[i].ni = 0;
285 } else {
286 // pre-compute inverse
287 mp_size_t ni = lmmp_div_inv_size_(np + powers[i].zeros, np);
288 lmmp_inv_prediv_(tp, p, np, ni);
289 powers[i].invp = tp;
290 powers[i].ni = ni;
291 tp += ni;
292 }
293 }
294
295 lmmp_copy(tp, numa, na);
296 digits = lmmp_to_str_divide_(dst, tp, na, powers + cpow - 1, tp + na + 1);
297
298 TEMP_FREE;
299 }
300
301 return digits;
302}
mp_basepow_t::ni
mp_size_t ni
Definition base_table.h:50
mp_basepow_t::invp
mp_ptr invp
Definition base_table.h:48
lmmp_from_str_len_
mp_size_t lmmp_from_str_len_(const mp_byte_t *src, mp_size_t len, int base)
计算字符串转大数所需的 limb 缓冲区长度
Definition from_str.c:13
TO_STR_BASEPOW_THRESHOLD
#define TO_STR_BASEPOW_THRESHOLD
Definition mparam.h:70
lmmp_to_str_basecase_
static mp_size_t lmmp_to_str_basecase_(mp_byte_t *dst, mp_srcptr numa, mp_size_t na, int base)
Definition to_str.c:26
lmmp_to_str_len_
mp_size_t lmmp_to_str_len_(mp_srcptr numa, mp_size_t na, int base)
计算大数转换为字符串，字符串需要的缓冲区长度
Definition to_str.c:12
lmmp_to_str_divide_
static mp_size_t lmmp_to_str_divide_(mp_byte_t *dst, mp_ptr numa, mp_size_t na, mp_basepow_t *pow, mp_ptr tpq)
Definition to_str.c:69

引用了 BALLOC_TYPE, mp_basepow_t::base, mp_basepow_t::digits, mp_base_t::digits_in_limb, DIV_MULINV_L_THRESHOLD, mp_basepow_t::invp, mp_base_t::large_base, LIMB_BITS, lmmp_bases_table, lmmp_copy, lmmp_div_inv_size_(), lmmp_from_str_len_(), lmmp_inv_prediv_(), lmmp_leading_zeros_(), lmmp_limb_bits_(), lmmp_mul_1_(), LMMP_POW2_Q, lmmp_shl_(), lmmp_sqr_(), lmmp_tailing_zeros_(), lmmp_to_str_basecase_(), lmmp_to_str_divide_(), lmmp_to_str_len_(), mp_basepow_t::ni, mp_basepow_t::np, mp_basepow_t::p, TEMP_DECL, TEMP_FREE, TO_STR_BASEPOW_THRESHOLD, tp , 以及 mp_basepow_t::zeros.

+ 函数调用图:

◆ lmmp_to_str_len_()

mp_size_t lmmp_to_str_len_ ( mp_srcptr  numa,
mp_size_t  na,
int  base 
)

计算大数转换为字符串，字符串需要的缓冲区长度

参数
numa输入大数，长度为na
na大数的 limb 长度
base目标基数（2~256）
返回
大数在指定基数下的位数
警告
na>=0, 2<=base<=256
注解
将会忽略numa的前导零，
  1. if (numa!=NULL) 返回的长度可能会多分配一个字符空间
  2. if (numa==NULL) 返回na个limb长度的数的最大可能字符长度（最坏情况）

在文件 to_str.c 第 12 行定义.

12 {
13 lmmp_param_assert(base >= 2 && base <= 256);
14 int mslbits = 0;
15 if (numa) {
16 do {
17 if (na == 0)
18 return 1;
19 } while (numa[--na] == 0);
20 mslbits = lmmp_limb_bits_(numa[na]);
21 }
22 return lmmp_mulh_(na * LIMB_BITS + mslbits, lmmp_bases_table[base - 2].inv_lg_base) + 1;
23}

引用了 LIMB_BITS, lmmp_bases_table, lmmp_limb_bits_(), lmmp_mulh_() , 以及 lmmp_param_assert.

被这些函数引用 lmmp_to_str_().

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

◆ lmmp_zero_q_()

static int lmmp_zero_q_ ( mp_srcptr  p,
mp_size_t  n 
)
inlinestatic

大数判零函数（内联）

参数
p指向大数起始位置的指针
n大数的单精度数(limb)长度
返回
1(全零) / 0(非零)
警告
n>0
注解
从最高位开始检查，只要有非零位则返回0

在文件 lmmpn.h 第 1027 行定义.

1027 {
1028 do {
1029 if (p[--n] != 0)
1030 return 0;
1031 } while (n != 0);
1032 return 1;
1033}

被这些函数引用 lmmp_invsqrt_newton_(), lmmp_mul_fermat_(), lmmp_mul_fermat_single_(), lmmp_mul_toom22_(), lmmp_mul_toom32_(), lmmp_mul_toom42_(), lmmp_mul_toom42_cache_init_(), lmmp_mul_toom52_(), lmmp_mul_toom62_(), lmmp_mul_toom62_cache_init_() , 以及 try_div_().

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