数论类
线性方程组
列主元
/*
* 列主元gauss消去求解a[][] * x[] = b[]
* 返回是否有唯一解，若有解在b[]中
*/
#define fabs(x) ((x) > 0 ? (x) : (-x))
#define eps 1e-10
const int MAXN = 100;
int gaussCpivot(int n, double a[][MAXN], double b[])
{
int i, j, k, row = 0;
double MAXP, temp;
for (k = 0; k < n; k++)
{
for (MAXP = 0, i = k; i < n; i++)
{
if (fabs(a[i][k]) > fabs(MAXP))
{
MAXP = a[row = i][k];
}
}
if (fabs(MAXP) < eps)
{
return 0;
}
if (row != k)
{
for (j = k; j < n; j++)
{
temp = a[k][j];
a[k][j] = a[row][j];
a[row][j] = temp;
temp = b[k];
b[k] = b[row];
b[row] = temp;
}
}
for (j = k + 1; j < n; j++)
{
a[k][j] /= MAXP;
for (i = k + 1; i < n; i++)
{
a[i][j] -= a[i][k] * a[k][j];
}
}
b[k] /= MAXP;
for (i = n - 1; i >= 0; i--)
{
for (j = i + 1; j < n; j++)
{
b[i] -= a[i][j] * b[j];
}
}
}
return 1;
}
全主元
/*
* 全主元gauss消去解a[][] * x[] = b[]
* 返回是否有唯一解，若有解在b[]中
*/
#define fabs(x) ((x) > 0 ? (x) : (-x))
#define eps 1e-10
const int MAXN = 100;
int gaussTpivot(int n, double a[][MAXN], double b[])
{
int i, j, k, row = 0, col = 0, index[MAXN];
double MAXP, temp;
for (i = 0; i < n; i++)
{
index[i] = i;
}
for (k = 0; k < n; k++)
{
for (MAXP = 0, i = k; i < n; i++)
{
for (j = k; j < n; j++)
{
if (fabs(a[i][j] > fabs(MAXP)))
{
MAXP = a[row = i][col = j];
}
}
}
if (fabs(MAXP) < eps)
{
return 0;
}
if (col != k)
{
for (i = 0; i < n; i++)
{
temp = a[i][col];
a[i][col] = a[i][k];
a[i][k] = temp;
}
j = index[col];
index[col] = index[k];
index[k] = j;
}
if (row != k)
{
for (j = k; j < n; j++)
{
temp = a[k][j];
a[k][j] = a[row][j];
a[row][j] = temp;
}
temp = b[k];
b[k] = b[row];
b[row] = temp;
}
for (j = k + 1; j < n; j++)
{
a[k][j] /= MAXP;
for (i = k + 1; i < n; i++)
{
a[i][j] -= a[i][k] * a[k][j];
}
}
b[k] /= MAXP;
for (i = k + 1; i < n; i++)
{
b[i] -= b[k] * a[i][k];
}
}
for (i = n - 1; i >= 0; i--)
{
for (j = i + 1; j < n; j++)
{
b[i] -= a[i][j] * b[j];
}
}
for (k = 0; k < n; k++)
{
a[0][index[k]] = b[k];
}
for (k = 0; k < n; k++)
{
b[k] = a[0][k];
}
return 1;
}
高斯消元（自由变元，一类开关问题，位运算操作）
// 高斯消元法求方程组的解
const int MAXN = 300;
// 有equ个方程，var个变元。增广矩阵行数为equ，列数为var＋1，分别为0到var
int equ, var;
int a[MAXN][MAXN]; // 增广矩阵
int x[MAXN]; // 解集
int free_x[MAXN]; // 用来存储自由变元（多解枚举自由变元可以使用）
int free_num; // 自由变元的个数
// 返回值为－1表示无解，为0是唯一解，否则返回自由变元个数
int Gauss()
{
int max_r, col, k;
free_num = 0;
for (k = 0, col = 0; k < equ && col < var; k++, col++)
{
max_r = k;
for (int i = k + 1; i < equ; i++)
{
if (abs(a[i][col]) > abs(a[max_r][col]))
{
max_r = i;
}
}
if (a[max_r][col] == 0)
{
k--;
free_x[free_num++] = col; // 这是自由变元
continue;
}
if (max_r != k)
{
for (int j = col; j < var + 1; j++)
{
swap(a[k][j], a[max_r][j]);
}
}
for (int i = k + 1; i < equ; i++)
{
if (a[i][col] != 0)
{
for (int j = col; j < var + 1; j++)
{
a[i][j] ^= a[k][j];
}
}
}
}
for (int i = k; i < equ; i++)
{
if (a[i][col] != 0)
{
return -1; // 无解
}
}
if (k < var)
{
return var - k; // 自由变元个数
}
// 唯一解，回代
for (int i = var - 1; i >= 0; i--)
{
x[i] = a[i][var];
for (int j = i + 1; j < var; j++)
{
x[i] ^= (a[i][j] && x[j]);
}
}
return 0;
}
模线性方程（组）
公共部分（扩展GCD）
int extgcd(int a, int b, int &x, int &y)
{
if (b == 0)
{
x = 1;
y = 0;
return a;
}
int d = extgcd(b, a % b, x, y);
int t = x;
x = y;
y = t - a / b * y;
return d;
}
模线性方程
/*
* 模线性方程 a * x = b (% n)
*/
void modeq(int a, int b, int n)
{
int e, i, d, x, y;
d = extgcd(b, a % b, x, y);
if (b % d > 0)
{
cout << "No answer!\n";
}
else
{
e = (x * (b / d)) % n;
for (i = 0; i < d; i++)
{
cout << i + 1 << "-th ans:" << (e + i * (n / d)) % n << '\n';
}
}
return ;
}
模线性方程组（互质）
/*
* 模线性方程组
* a = B[1](% W[1]); a = B[2](% W[2]); ... a = B[k](% W[k]);
* 其中W，B已知，W[i] > 0且W[i]与W[j]互质，求a（中国剩余定理）
*/
int china(int b[], int w[], int k)
{
int i, d, x, y, m, a = 0, n = 1;
for (i = 0; i < k; i++)
{
n *= w[i]; // 注意不能overflow
}
for (i = 0; i < k; i++)
{
m = n / w[i];
d = extgcd(w[i], m, x, y);
a = (a + y * m * b[i]) % n;
}
if (a > 0)
{
return a;
}
else
{
return (a + n);
}
}
模线性方程组（不要求互质）
typedef long long ll;
const int MAXN = 11;
int n, m;
int a[MAXN], b[MAXN];
int main(int argc, const char * argv[])
{
int T;
cin >> T;
while (T--)
{
cin >> n >> m;
for (int i = 0; i < m; i++)
{
cin >> a[i];
}
for (int i = 0; i < m; i++)
{
cin >> b[i];
}
ll ax = a[0], bx = b[0], x, y;
int flag = 0;
for (int i = 1; i < m; i++)
{
ll d = extgcd(ax, a[i], x, y);
if ((b[i] - bx) % d != 0)
{
flag = 1; // 无整数解
break;
}
ll tmp = a[i] / d;
x = x * (b[i] - bx) / d; // 约分
x = (x % tmp + tmp) % tmp;
bx = bx + ax * x;
ax = ax * tmp; // ax = ax * a[i] / d
}
if (flag == 1 || n < bx)
{
puts("0");
}
else
{
ll ans = (n - bx) / ax + 1;
if (bx == 0)
{
ans--;
}
printf("%lld\n", ans);
}
}
return 0;
}
素数相关
随机素数测试
/*
* 随机素数测试（伪素数原理）
* CALL: bool res = miller(n);
* 快速测试n是否满足素数的“必要”条件，出错概率极低
* 对于任意奇数n > 2和正整数s，算法出错概率≤2^(-s)
*/
int witness(int a, int n)
{
int x, d = 1;
int i = ceil(log(n - 1.0) / log(2.0)) - 1;
for (; i >= 0; i--)
{
x = d;
d = (d * d) % n;
if (d == 1 && x != 1 && x != n - 1)
{
return 1;
}
if (((n - 1) & (1 << i)) > 0)
{
d = (d * a) % n;
}
}
return (d == 1 ? 0 : 1);
}
int miller(int n, int s = 50)
{
if (n == 2) // 质数返回1
{
return 1;
}
if (n % 2 == 0) // 偶数返回0
{
return 0;
}
int j, a;
for (j = 0; j < a; j++)
{
a = rand() * (n - 2) / RAND_MAX + 1;
// rand()只能随机产生[0, RAND_MAX)内的整数
// 而且这个RAND_MAX只有32768直接%n的话是永远
// 也产生不了[RAND_MAX, n)之间的数
if (witness(a, n))
{
return 0;
}
}
return 1;
}
大数素数测试
#define MAXL 4
#define M10 1000000000
#define Z10 9
const int zero[MAXL - 1] = {0};
struct bnum
{
int data[MAXL]; // 断成每截9个长度
// 读取字符串并转存
void read()
{
memset(data, 0, sizeof(data));
char buf[32];
scanf("%s", buf);
int len = (int)strlen(buf);
int i = 0, k;
while (len >= Z10)
{
for (k = len - Z10; k < len; ++k)
{
data[i] = data[i] * 10 + buf[k] - '0';
}
++i;
len -= Z10;
}
if (len > 0)
{
for (k = 0; k < len; ++k)
{
data[i] = data[i] * 10 + buf[k] - '0';
}
}
}
bool operator == (const bnum &x)
{
return memcmp(data, x.data, sizeof(data)) == 0;
}
bnum & operator = (const int x)
{
memset(data, 0, sizeof(data));
data[0] = x;
return *this;
}
bnum operator + (const bnum &x)
{
int i, carry = 0;
bnum ans;
for (i = 0; i < MAXL; ++i)
{
ans.data[i] = data[i] + x.data[i] + carry;
carry = ans.data[i] / M10;
ans.data[i] %= M10;
}
return ans;
}
bnum operator - (const bnum &x)
{
int i, carry = 0;
bnum ans;
for (i = 0; i < MAXL; ++i)
{
ans.data[i] = data[i] - x.data[i] - carry;
if (ans.data[i] < 0)
{
ans.data[i] += M10;
carry = 1;
}
else
{
carry = 0;
}
}
return ans;
}
// assume *this < x * 2
bnum operator % (const bnum &x)
{
int i;
for (i = MAXL - 1; i >= 0; --i)
{
if (data[i] < x.data[i])
{
return *this;
}
else if (data[i] > x.data[i])
{
break;
}
}
return ((*this) - x);
}
bnum & div2()
{
int i, carry = 0, tmp;
for (i = MAXL - 1; i >= 0; --i)
{
tmp = data[i] & 1;
data[i] = (data[i] + carry) >> 1;
carry = tmp * M10;
}
return *this;
}
bool is_odd()
{
return (data[0] & 1) == 1;
}
bool is_zero()
{
for (int i = 0; i < MAXL; ++i)
{
if (data[i])
{
return false;
}
}
return true;
}
};
void mulmod(bnum &a0, bnum &b0, bnum &p, bnum &ans)
{
bnum tmp = a0, b = b0;
ans = 0;
while (!b.is_zero())
{
if (b.is_odd())
{
ans = (ans + tmp) % p;
}
tmp = (tmp + tmp) % p;
b.div2();
}
}
void powmod(bnum &a0, bnum &b0, bnum &p, bnum &ans)
{
bnum tmp = a0, b = b0;
ans = 1;
while (!b.is_zero())
{
if (b.is_odd())
{
mulmod(ans, tmp, p, ans);
}
mulmod(tmp, tmp, p, tmp);
b.div2();
}
}
bool MillerRabinTest(bnum &p, int iter)
{
int i, small = 0, j, d = 0;
for (i = 1; i < MAXL; ++i)
{
if (p.data[i])
{
break;
}
}
if (i == MAXL)
{
// small integer test
if (p.data[0] < 2)
{
return false;
}
if (p.data[0] == 2)
{
return true;
}
small = 1;
}
if (!p.is_odd())
{
return false; // even number
}
bnum a, s, m, one, pd1;
one = 1;
s = pd1 = p - one;
while (!s.is_odd())
{
s.div2();
++d;
}
for (i = 0; i < iter; ++i)
{
a = rand();
if (small)
{
a.data[0] = a.data[0] % (p.data[0] - 1) + 1;
}
else
{
a.data[1] = a.data[0] / M10;
a.data[0] %= M10;
}
if (a == one)
{
continue;
}
powmod(a, s, p, m);
for (j = 0; j < d && !(m == one) && !(m == pd1); ++j)
{
mulmod(m, m, p, m);
}
if (!(m == pd1) && j > 0)
{
return false;
}
}
return true;
}
int main()
{
bnum x;
x.read();
puts(MillerRabinTest(x, 5) ? "Yes" : "No");
return 0;
}
组合数
组合数C(a, b) (预处理)
typedef long long ll;
const ll MOD = 1e9 + 7; // 必须为质数才管用
const ll MAXN = 1e5 + 3;
ll fac[MAXN]; // 阶乘
ll inv[MAXN]; // 阶乘的逆元
ll QPow(ll x, ll n)
{
ll ret = 1;
ll tmp = x % MOD;
while (n)
{
if (n & 1)
{
ret = (ret * tmp) % MOD;
}
tmp = tmp * tmp % MOD;
n >>= 1;
}
return ret;
}
void init()
{
fac[0] = 1;
for (int i = 1; i < MAXN; i++)
{
fac[i] = fac[i - 1] * i % MOD;
}
inv[MAXN - 1] = QPow(fac[MAXN - 1], MOD - 2);
for (int i = MAXN - 2; i >= 0; i--)
{
inv[i] = inv[i + 1] * (i + 1) % MOD;
}
}
ll C(ll a, ll b)
{
if (b > a)
{
return 0;
}
if (b == 0)
{
return 1;
}
return fac[a] * inv[b] % MOD * inv[a - b] % MOD;
}
集合划分问题
/*
* n元集合分划为k类的方案数记为S(n, k),称为第二类Stirling数。
* 如{A,B,C}可以划分{{A}, {B}, {C}}, {{A, B}, {C}}, {{B, C}, {A}}, {{A, C}, {B}}, {{A, B, C}}。
* 即一个集合可以划分为不同集合(1...n个)的种类数
* CALL: compute(N); 每当输入一个n,输出B[n]
*/
const int N = 2001;
int data[N][N], B[N];
void NGetM(int m, int n) // m 个数 n 个集合
{
// data[i][j]: i个数分成j个集合
int min, i, j;
data[0][0] = 1;
for (i = 1; i <= m; i++)
{
data[i][0] = 0;
}
for (i = 0; i <= m; i++)
{
data[i][i + 1] = 0;
}
for (i = 1; i <= m; i++)
{
if (i < n)
{
min = i;
}
else
{
min = n;
}
for (j = 1; j <= min; j++)
{
data[i][j] = (j * data[i - 1][j] + data[i - 1][j - 1]);
}
}
return ;
}
void compute(int m)
{
// b[i]: Bell数
NGetM(m, m);
memset(B, 0, sizeof(B));
int i, j;
for (i=1; i <= m; i++)
{
for (j = 0; j <= i; j++)
{
B[i] += data[i][j];
}
}
return ;
}
卢卡斯定理（从(1, 1)到(n, m)的走法，机器人走方格问题）
#define MOD 1000000007
typedef long long LL;
LL quickPower(LL a, LL b)
{
LL ans = 1;
a %= MOD;
while (b)
{
if (b & 1)
{
ans = ans * a % MOD;
}
b >>= 1;
a = a * a % MOD;
}
return ans;
}
LL c(LL n, LL m)
{
if (m > n)
{
return 0;
}
LL ans = 1;
for (int i = 1; i <= m; i++)
{
LL a = (n + i - m) % MOD;
LL b = i % MOD;
ans = ans * (a * quickPower(b, MOD - 2) % MOD) % MOD;
}
return ans;
}
LL lucas(LL n, LL m)
{
if (m == 0)
{
return 1;
}
return c(n % MOD, m % MOD) * lucas(n / MOD, m / MOD) % MOD;
}
int main(int argc, const char * argv[])
{
LL n, m;
while (~scanf("%lld %lld", &n, &m))
{
LL max, min;
max = n + m - 3;
min = m - 1;
printf("%lld\n", lucas(max - 1, min - 1));
}
return 0;
}
Polya计数
/*
* c种颜色的珠子，组成长为s的项链，项链没有方向和起始位置
*/
int gcd(int a, int b)
{
return b ? gcd(b, a % b) : a;
}
int main(int argc, const char * argv[])
{
int c, s;
while (cin >> c >> s)
{
int k;
long long p[64];
p[0] = 1; // power of c
for (k = 0; k < s; k++)
{
p[k + 1] = p[k] * c;
}
// reflection part
long long count = s & 1 ? s * p[s / 2 + 1] : (s / 2) * (p[s / 2] + p[s / 2 + 1]);
// rotation part
for (k = 1 ; k <= s ; k++)
{
count += p[gcd(k, s)];
count /= 2 * s;
}
cout << count << '\n';
}
return 0;
}
解周期性方程
/*
* 周期性方程定义(n = 5)
* |a_1 b_1 c_1 d_1 e_1| = x_1 --- 1
* |e_2 a_2 b_2 c_2 d_2| = x_2 --- 2
* |d_2 e_2 a_2 b_2 c_2| = x_3 --- 3
* |c_4 d_2 e_2 a_4 b_4| = x_4 --- 4
* |b_5 c_5 d_5 e_5 a_5| = x_5 --- 5
* 输入： a[], b[], c[], x[]
* 输出： 求解结果x在x[]中
*/
const int MAXN = 1000;
int a[MAXN];
int b[MAXN];
int c[MAXN];
int x[MAXN];
void run()
{
c[0] /= b[0];
a[0] /= b[0];
x[0] /= b[0];
for (int i = 1; i < MAXN - 1; i++)
{
double temp = b[i] - a[i] * c[i - 1];
c[i] /= temp;
x[i] = (x[i] - a[i] * x[i - 1]) / temp;
a[i] = -a[i] * a[i - 1] / temp;
}
a[MAXN - 2] = -a[MAXN - 2] - c[MAXN - 2];
for (int i = MAXN - 3; i >= 0; i--)
{
a[i] = -a[i] - c[i] * a[i + 1];
x[i] -= c[i] * x[i + 1];
}
x[MAXN - 1] -= (c[MAXN - 1] * x[0] + a[MAXN - 1] * x[MAXN - 2]);
x[MAXN - 1] /= (c[MAXN - 1] * a[0] + a[MAXN - 1] * a[MAXN - 2] + b[MAXN - 1]);
for (int i = MAXN - 2; i >= 0; i --)
{
x[i] += a[i] * x[MAXN - 1];
}
return ;
}
阶乘最后非零位
/*
* 阶乘最后非零位 复杂度O(nlongn)
* 返回改为，n以字符串方式传入
*/
#define MAXN 10000
const int mod[20] = {1, 1, 2, 6, 4, 2, 2, 4, 2, 8, 4, 4, 8, 4, 6, 8, 8, 6, 8, 2};
int lastDigit(char *buf)
{
int len = (int)strlen(buf);
int a[MAXN], i, c, ret = 1;
if (len == 1)
{
return mod[buf[0] - '0'];
}
for (i = 0; i < len; i++)
{
a[i] = buf[len - 1 - i] - '0';
}
for (; len; len -= !a[len - 1])
{
ret = ret * mod[a[1] % 2 * 10 + a[0]] % 5;
for (c = 0, i = len - 1; i >= 0; i--)
{
c = c * 10 + a[i];
a[i] = c / 5;
c %= 5;
}
}
return ret + ret % 2 * 5;
}
n的阶乘的长度
#define PI 3.1415926
int main()
{
int n, a;
while (~scanf(“%d", &n))
{
a = (int)((0.5 * log(2 * PI * n) + n * log(n) - n) / log(10));
printf("%d\n", a + 1);
}
return 0;
}
求阶乘逆元
typedef long long ll;
const ll MOD = 1e9 + 7; // 必须为质数才管用
const ll MAXN = 1e5 + 3;
ll fac[MAXN]; // 阶乘
ll inv[MAXN]; // 阶乘的逆元
ll QPow(ll x, ll n)
{
ll ret = 1;
ll tmp = x % MOD;
while (n)
{
if (n & 1)
{
ret = (ret * tmp) % MOD;
}
tmp = tmp * tmp % MOD;
n >>= 1;
}
return ret;
}
void init()
{
fac[0] = 1;
for (int i = 1; i < MAXN; i++)
{
fac[i] = fac[i - 1] * i % MOD;
}
inv[MAXN - 1] = QPow(fac[MAXN - 1], MOD - 2);
for (int i = MAXN - 2; i >= 0; i--)
{
inv[i] = inv[i + 1] * (i + 1) % MOD;
}
}
整数划分
整数划分（五边形定理）
P
(
n
)
=
∑
P
(
n
−
k
(
3
k
−
1
)
/
2
+
P
(
n
−
k
(
3
k
+
1
)
/
2
,
k
≥
1
n
<
0
时
，
P
(
n
)
=
0
,
n
=
0
时
,
P
(
n
)
=
1
即
可
P(n) = ∑{P(n - k(3k - 1) / 2 + P(n - k(3k + 1) / 2 , k ≥ 1} n < 0时，P(n) = 0, n = 0时, P(n) = 1即可
P(n)=∑P(n−k(3k−1)/2+P(n−k(3k+1)/2,k≥1n<0时，P(n)=0,n=0时,P(n)=1即可
// 划分元素可重复任意次
#define f(x) (((x) * (3 * (x) - 1)) >> 1)
#define g(x) (((x) * (3 * (x) + 1)) >> 1)
const int MAXN = 1e5 + 10;
const int MOD = 1e9 + 7;
int n, ans[MAXN];
int main()
{
scanf("%d", &n);
ans[0] = 1;
for (int i = 1; i <= n; ++i)
{
for (int j = 1; f(j) <= i; ++j)
{
if (j & 1)
{
ans[i] = (ans[i] + ans[i - f(j)]) % MOD;
}
else
{
ans[i] = (ans[i] - ans[i - f(j)] + MOD) % MOD;
}
}
for (int j = 1; g(j) <= i; ++j)
{
if (j & 1)
{
ans[i] = (ans[i] + ans[i - g(j)]) % MOD;
}
else
{
ans[i] = (ans[i] - ans[i - g(j)] + MOD) % MOD;
}
}
}
printf("%d\n", ans[n]);
return 0;
}
整数划分（五边形定理拓展）
F
(
n
,
k
)
=
P
(
n
)
−
划
分
元
素
重
复
次
数
≥
k
次
的
情
况
F(n, k) = P(n) - 划分元素重复次数≥k次的情况
F(n,k)=P(n)−划分元素重复次数≥k次的情况。
// 问一个数n能被拆分成多少种情况
// 且要求拆分元素重复次数不能≥k
const int MOD = 1e9 + 7;
const int MAXN = 1e5 + 10;
int ans[MAXN];
// 此函数求ans[]效率比上一个代码段中求ans[]效率高很多
void init()
{
memset(ans, 0, sizeof(ans));
ans[0] = 1;
for (int i = 1; i < MAXN; ++i)
{
ans[i] = 0;
for (int j = 1; ; j++)
{
int tmp = (3 * j - 1) * j / 2;
if (tmp > i)
{
break;
}
int tmp_ = ans[i - tmp];
if (tmp + j <= i)
{
tmp_ = (tmp_ + ans[i - tmp - j]) % MOD;
}
if (j & 1)
{
ans[i] = (ans[i] + tmp_) % MOD;
}
else
{
ans[i] = (ans[i] - tmp_ + MOD) % MOD;
}
}
}
return ;
}
int solve(int n, int k)
{
int res = ans[n];
for (int i = 1; ; i++)
{
int tmp = k * i * (3 * i - 1) / 2;
if (tmp > n)
{
break;
}
int tmp_ = ans[n - tmp];
if (tmp + i * k <= n)
{
tmp_ = (tmp_ + ans[n - tmp - i * k]) % MOD;
}
if (i & 1)
{
res = (res - tmp_ + MOD) % MOD;
}
else
{
res = (res + tmp_) % MOD;
}
}
return res;
}
int main(int argc, const char * argv[])
{
init();
int T, n, k;
cin >> T;
while (T--)
{
cin >> n >> k;
cout << solve(n, k) << '\n';
}
return 0;
}
A^B约数之和对MOD取模
/*
* 求A^B的约数之和对MOD取模
* 需要素数筛选和合数分解的算法，需要先调用getPrime();
* 1+p+p^2+p^3+...+p^n
*/
const int MOD = 1000000;
long long pow_m(long long a, long long n)
{
long long ret = 1;
long long tmp = a % MOD;
while(n)
{
if (n & 1)
{
ret = (ret * tmp) % MOD;
}
tmp = tmp * tmp % MOD;
n >>= 1;
}
return ret;
}
// 计算1+p＋p^2+...+p^n
long long sum(long long p, long long n)
{
if (p == 0)
{
return 0;
}
if (n == 0)
{
return 1;
}
if (n & 1)
{
return ((1 + pow_m(p, n / 2 + 1)) % MOD * sum(p, n / 2) % MOD) % MOD;
}
else
{
return ((1 + pow_m(p, n / 2 + 1)) % MOD * sum(p, n / 2 - 1) + pow_m(p, n / 2) % MOD) % MOD;
}
}
// 返回A^B的约数之和%MOD
long long solve(long long A, long long B)
{
getFactors(A);
long long ans = 1;
for (int i = 0; i < fatCnt; i++)
{
ans *= sum(factor[i][0], B * factor[i][1]) % MOD;
ans %= MOD;
}
return ans;
}
BSGS
/*
* baby_step giant _step
* a^x = b(mod n) n不要求是素数
* 求解上式0 ≤ x < n的解
*/
#define MOD 76543
int hs[MOD];
int head[MOD];
int _next[MOD];
int id[MOD];
int top;
void insert(int x, int y)
{
int k = x % MOD;
hs[top] = x;
id[top] = y;
_next[top] = head[k];
head[k] = top++;
return ;
}
int find(int x)
{
int k = x % MOD;
for (int i = head[k]; i != -1; i = _next[i])
{
if (hs[i] == x)
{
return id[i];
}
}
return -1;
}
long long BSGS(int a, int b, int n)
{
memset(head, -1, sizeof(head));
top = 1;
if (b == 1)
{
return 0;
}
int m = (int)sqrt(n * 1.0), j;
long long x = 1, p = 1;
for (int i = 0; i < m; i++, p = p * a % n)
{
insert(p * b % n, i);
}
for (long long i = m; ; i++)
{
if ((j = find(x = x * p % n)) != -1)
{
return i - j;
}
if (i > n)
{
break;
}
}
return -1;
}
多项式求根（牛顿法）
/*
* 牛顿法解多项式的根
* 输入:多项式系数c[],多项式度数n,求在[a,b]间的根
* 输出:根 要求保证[a,b]间有根
*/
double fabs(double x)
{
return (x < 0) ? -x : x;
}
double f(int m, double c[], double x)
{
int i;
double p = c[m];
for (i = m; i > 0; i--)
{
p = p * x + c[i - 1];
}
return p;
}
int newton(double x0, double *r, double c[], double cp[], int n, double a, double b, double eps)
{
int MAX_ITERATION = 1000;
int i = 1;
double x1, x2, fp, eps2 = eps / 10.0;
x1 = x0;
while (i < MAX_ITERATION)
{
x2 = f(n, c, x1);
fp = f(n - 1, cp, x1);
if ((fabs(fp) < 0.000000001) && (fabs(x2) > 1.0))
{
return 0;
}
x2 = x1 - x2 / fp;
if (fabs(x1 - x2) < eps2)
{
if (x2 < a || x2 > b)
{
return 0;
}
*r = x2;
return 1;
}
x1 = x2;
i++;
}
return 0;
}
double Polynomial_Root(double c[], int n, double a, double b, double eps)
{
double *cp;
int i;
double root;
cp = (double *)calloc(n, sizeof(double));
for (i = n - 1; i >= 0; i--)
{
cp[i] = (i + 1) * c[i + 1];
}
if (a > b)
{
root = a;
a = b;
b = root;
}
if ((!newton(a, &root, c, cp, n, a, b, eps)) && (!newton(b, &root, c, cp, n, a, b, eps)))
{
newton((a + b) * 0.5, &root, c, cp, n, a, b, eps);
}
free(cp);
if (fabs(root) < eps)
{
return fabs(root);
}
else
return root;
}
星期问题
基姆拉尔森公式：
W
=
(
D
+
2
∗
M
+
3
∗
(
M
+
1
)
5
+
Y
+
Y
4
−
Y
100
+
Y
400
)
M
o
d
7
W = (D + 2 * M + 3 * (M + 1) \ 5 + Y + Y \ 4 - Y \ 100 + Y \ 400) Mod 7
W=(D+2∗M+3∗(M+1) 5+Y+Y 4−Y 100+Y 400)Mod7
基姆拉尔森公式的计算结果是0，1，2，3，4，5，6 七种可能；
结果的对应关系：
0：星期一
1：星期二
2：星期三
3：星期四
4：星期五
5：星期六
6：星期日
/*
* 已知1752年9月3日是Sunday，并且日期控制在1700年2月28日后
*/
char name[][15] = { "monday", "tuesday", "wednesday", "thursday", "friday", "saturday", "sunday"};
int main()
{
int d, m, y, a;
printf("Day: ");
scanf("%d", &d);
printf("Month: ");
scanf("%d", &m);
printf("Year: ");
scanf("%d", &y);
// 1月2月当作前一年的13,14月
if (m == 1 || m == 2)
{
m += 12;
y--;
}
// 判断是否在1752年9月3日之前,实际上合并在一起倒更加省事
if ((y < 1752) || (y == 1752 && m < 9) || (y == 1752 && m == 9 && d < 3))
{
// 因为日期控制在1700年2月28日后，所以不用考虑整百年是否是闰年
a = (d + 2 * m + 3 * (m + 1) / 5 + y + y / 4 + 5) % 7;
}
else
{
// 这里需要考虑整百年是否是闰年的情况
a = (d + 2 * m + 3 * (m + 1) / 5 + y + y / 4 - y / 100 + y / 400) % 7; // 实际上这个可以当做公式背下来
}
printf("it's a %s\n", name[a]);
return 0;
}
1/n循环节长度
/*
* 求1/i的循环节长度的最大值，i<=n
*/
const int MAXN = 1005;
int res[MAXN]; // 循环节长度
int main()
{
memset(res, 0, sizeof(res));
int i, temp, j, n;
for (temp = 1; temp <= 1000; temp++)
{
i = temp;
while (i % 2 == 0)
{
i /= 2;
}
while (i % 5 == 0)
{
i /= 5;
}
n = 1;
for (j = 1; j <= i; j++)
{
n *= 10;
n %= i;
if (n == 1)
{
res[temp] = j;
break;
}
}
}
int max_re;
while (cin >> n)
{
max_re = 1;
for (i = 1; i <= n; i++)
{
if (res[i] > res[max_re])
{
max_re = i;
}
}
cout << max_re << endl;
}
return 0;
}
母函数
/*
* 母函数
* c1是保存各项质量砝码可以组合的数目
* c2是中间量，保存每一次的情况
*/
const int MAXN = 1e4 + 10;
int n;
int c1[MAXN];
int c2[MAXN];
int main()
{
while (cin >> n)
{
for (int i = 0; i <= n; ++i)
{
c1[i] = 1;
c2[i] = 0;
}
for (int i = 2; i <= n; ++i)
{
for (int j = 0; j <= n; ++j)
{
for (int k = 0; k + j <= n; k += i)
{
c2[j + k] += c1[j];
}
}
for (int j = 0; j <= n; ++j)
{
c1[j] = c2[j];
c2[j] = 0;
}
}
cout << c1[n] << endl;
}
return 0;
}
扫一扫
1703
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
