问题

一个无向联通图,应该要正边权?，最小割分成两个子图

方法

引理

图G中的两个点s,t,

如果最小割分割了它们，那么直接求s到t的最小割

如果没有，那么合并s-t，得到的新图G1 和原图的最小割一致

因此有方法, 任意选择点a，加入集合A中

每次选择和A距离和最大的点加入集合中，也就是上面的合并操作

注意，这里每次都不是把图里的最大的边去掉，而是集合A距离最大的点加入A中

这样反复操作，一直到剩下3个点A,b,c ，可以得到一个割的值, 除了A以外的两个点b,c在这次分割中一定是被分开了的

我们也能得到一个b和c是被分割的值

重新初始的图

合并 b和c成集合B

重复上面所有操作，直到剩余3个点，B,d,e

重新初始的图

合并b,c成集合B, 合并d和e成集合D

重复上面所有操作，直到剩余3个点，?a,?b,?c,

我们能到到一个分割值，也能让下一轮的操作合并两个点

这样我们最终所有操作中最小的分割就是答案

Proof

对于 每一轮剩余3个点的时候，我们有该轮的初始集合S，和另外两点$p_1,p_2$

因为合并规则一直是与S距离最大的边，所以$p_1,p_2$ 一定是被分割的，

虽然因为边的长度可能在操作过程中相等，也就意味着在每一轮操作开始时，甚至$p_1,p_2$ 是哪两个点，是不确定的,(当然我们可以给所有点提前编号，保持即使值相等的点也有一定的偏序关系，来让结果唯一)

但是可以证明，对于每次计算到剩余3个点时，其中的$p_1,p_2$来说，当前图里面这是这两个点要切割的最小割

我们把这轮的$(p_1,p_2)$ 命作$(s,t)$

对于给定图$G$，和通过上述步骤得到的$s,t$

如果两个点之间没有相连，我们在数值上可以看做有一条相连权重为0的边,在下面的计算上没有影响(因为每次取得是最大的), 连通的判断时权重为0视为不连通

割: 把原图分成两个联通块的边集

$C$是$G$上$s-t$的一个任意的割,即$s,t$在不同连通块

$CP$是按照上述步骤得到的$s-t$的割

$W(边集)$ = 边集的权重和

$W(点集A，点v)$ = $点v$到$点集A$中所有点的边权和

要证:$W(C) \ge W(CP)$

对于上述给定流程$CP$，根据点加入的顺序，命名为$a_1,a_2,a_3…a_{n-1},a_n$, 其中$s = a_{n-1},t = a_n$

对于 $i$,如果$a_{< i-1}$和$a_{i}$ 位于 割C的两侧,那么称为$i$为$active$的(简化后面描述)

• 我们不用讨论$CP$，因为序列就是$CP$产生的，所有操作都是位于$CP$割的同一连通块

点集 $A_i = 集合(a_1…a_{i-1}) $

也就是$a_i$前面的点

边集 $C_i = 集合(边|边的两端均属于集合(a_1…a_i)) ∩ C$

也就是割边中，端点都在前面$a_1..a_i$ 中

我们要证明对于每一个$active$的$i$, $W(A_i, a_i) \le W(C_i)$,

也就是$a_i$和它之前的点的边权和 小于 割$C$中的两端来自前$a_1 … a_i$的边权和

接下来归纳开始 (归纳的目标和总目标一致)

初始: 若$i$是最小使得$active$的，那么$W(A_i,a_i) = W(C_i)$

因为 说明 前面的边,按照割C，都在同一侧，所以$C_i$中所有的边，都一个顶点是$a_i$

对于$i,j,(i < j)$ 都是$active$相邻

$i$和$j$是active, $i$和$j$之间没有其它$active$的,

$W(A_j,a_j) = W(A_i,a_j) + W(A_j - A_i, a_j) $

直接的拆分点集$A_j$

$\le W(C_i) + W(A_j - A_i, a_j) $

根据初始条件和归纳,$W(A_i,a_j)\le W(A_i,a_i) \le W(C_i)$ , 前一半不等式是操作过程中每次选最大距离的连通，后一部分是归纳

$= W(C_j)$

因为$i$和$j$之间没有其它$active$，所以也就是$i$和$j$之间的点不会与小于$j$的点构成属于$C_u$的边，
而对于$W(A_j-A_i,a_j)$的贡献，不会和$C_i$重复因为(不包含$a_j$),

综上

$W(A_n,a_n) \le W(C_n) = W(C)$

左边不等式是归纳的结果，右边是因为所有C中的边且两端的点在所有点中，就是$C$自己

结论:

也就是说，对于图$G$的$s-t$任何的割$C$，都不小于$t$单独划分出来的边权和

所以这是一个最小分割

回到操作过程

也就是每一轮，都能产生$s-t$分割的最小割，之后的轮次，$s$和$t$就在同一个“连通块”里了

这里有一点 实际上再思考仔细一点，这个过程中其实并不满足题意(分成两个)，因为我们很可能把原图分割成了好几个连通块，也就是不止两个。但是因为如果能分割成多个，反向操作，每次最多把两个合成一个，必定最小割是分成两个的。所以得到的答案也一定是割, 保证了正确性

也就是你合并的s和t可能，它们还并不连通，只不过你从划分上认为属于一块了

综上，算法正确性证毕

Ref

https://basics.sjtu.edu.cn/~dominik/teaching/2016-cs214/presentation-slides/2016-12-06-StoerWagner-BigNews.pdf

https://en.wikipedia.org/wiki/Stoer%E2%80%93Wagner_algorithm

题目

原题链接

直角三角形$(a,b,c)$,其中$c$为斜边

perfect定义

1. 边互质$gcd(a,b) = 1$
2. 斜边是平方数$c = x^2$

super-perfect

1. perfect
2. 直角三角形面积是6和28的倍数

求

$c \le 10^{16}$ 内,有多少perfect的不是super-perfect

题解

前面的PE我们学习到，$gcd(a,b)=1$ 的直角三角形可以唯一表示成

$(a,b,c) = (m^2-n^2,2mn,m^2+n^2)$

其中两个直角边可以对换

其中 $gcd(m,n) = 1, (m+n) = 1 (\bmod 2)$

$lcm(6,28) = 84$

简化成

$m^2+n^2 = c = x^2 \le 10^{16} $

求 $mn(m^2-n^2) != 0 (\bmod 84)$ 的个数

注意到 这里又是一个直角三角形公式

$ordered(m,n,x) = ordered(2uv,u^2-v^2,u^2+v^2)$

即

$u > v, gcd(u,v) = 1, u+v = 1 (\bmod 2)$

$u^2+v^2 \le 10^{8}$

$ordered(m,n) = ordered(2uv,u^2-v^2) $

求 $mn(m^2-n^2) != 0 (\bmod 84)$ 的个数

代码

直接翻译很好写了

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28

perfect_but_not_super_perfect = 0

def gcd(a,b):
    if b == 0:
        return a
    return gcd(b,a%b)


for v in range(1,10**4):
    if v % 100 == 0:
        print("b:", v)
    for u in range(v+1,10**4,2):
        if gcd(u,v) != 1:
            continue

        if u**2+v**2 > 10**8:
            break;

        m,n = u**2-v**2,2*u*v

        if n > m:
            m,n=n,m

        if (m*n*(m**2-n**2)) % 84 != 0:
            perfect_but_not_super_perfect += 1


print("ans",perfect_but_not_super_perfect )

时间

1
2
3

real    0m30.195s
user    0m30.195s
sys     0m0.000s

不要到此为止 上数学

就完了吗

它要是个长长的数字可能也没这篇文章了，因为PE的 one-minute rule 已经满足了

但是它的答案是0那就激起了推导的兴趣

命题

$u > v, gcd(u,v) = 1, u+v = 1 (\bmod 2)$

$(m,n) = (2uv,u^2-v^2) $

则 $mn(m^2-n^2) = 0 (\bmod 84)$

这里我们没有范围限制，考虑模运算的性质，没有了大小顺序正负

代入

$2uv(u^2-v^2)(6u^2v^2-u^4-v^4) = 0 (\bmod 84)$

$uv(u^2-v^2)(6u^2v^2-u^4-v^4) = 0 (\bmod 42)$

$42 = 2 \cdot 3 \cdot 7$

讨论

2 显然，因为$u+v=1(\bmod 2)$ , u 和 v一奇一偶

3 如果u和v其中有3的倍数，那么是3的倍数，如果均不是3的倍数,因为 $1^2=2^2=1 (\bmod 3)$ , 所以 $u^2-v^2 = 0 (\bmod 3)$

7 如果u和v其中有7的倍数，那么是7的倍数，否则都不是7的倍数，还是写一点代码

1
2
3

for u in range(1,7):
    for v in range(1,7):
        print(u,v,((u**2-v**2)*(6*u*u*v*v-u**4-v**4))%7)

输出都是0 得证

题目

1<= n <=50’000’000

求让 $2n^2-1$ 为质数的有多少

尝试

裸暴力

显然我们枚举n,然后每个计算到 $\sqrt{n}$, 时间复杂度$O(n^{1.5})$

在n取这么大,要大概一个小时, 显然不满足pe的期望1min

优化1

显然如果 $ 2n^2-1 = 0 (\bmod b)$

那么有 $ 2(b+n)^2-1 = 0 (\bmod b)$, 且这是充要的

这有两个好的性质,并不要求$2n^2-1$是质数或者合数

对于给定b总能找到一个小于b的n

看上去很美好,而实际操作上,也能降低一定的复杂度,不过这个方法空间要求很大,而上面的暴力是O(1)空间

其2每次进行筛的时候也不算快,也有不少重复筛的

优化2

我们存不下$2n^2-1$的质数,但是可以预运算$\sqrt{2n^2-1}$的质数 进行一定的提速

如何判断质数

Level 1

裸暴力,从2到开根

Level 2

众所周知 费马小定理$a^{p-1} = 1 (\bmod p), gcd(a,p) = 1$

当然这个表达式 p是质数一定成立,但是成立并不一定是质数

Level 3

尝试多个a的费马小定理,

依然是”概率”运气的判断,

依然可能表达式成立,但p不是质数

Level 4 Miller-Rabin

把费马小定理的$p-1$拆解成$p-1 = 2^k (2m+1), m\ge 0$

也就是 2幂次和奇数的乘积

如果p是质数,我们有 $2^{2^k(2m+1)} = 1 (\bmod p)$

也就是 $2^{2^{k-1}(2m+1)} = \pm 1 (\bmod p)$

如果 $2^{2^{k-1}(2m+1)} = 1 (\bmod p)$

那么 $2^{2^{k-2}(2m+1)} = \pm 1 (\bmod p)$

同样,我们可以尝试多个a

这样的话”运气”能好很多,本身还是概率的做法

Level 5

http://miller-rabin.appspot.com/

我们用别人的答案!

Jim Sinclair证明了如果a的取值把下面都测试一遍,那么在 $p<2^{64}$时,都是能判断质数的!

2, 325, 9375, 28178, 450775, 9780504, 1795265022

还是Miller-Rabin,不过我们使用了已有成果,(好奇是证明,还是暴力枚举,毕竟是2011-04-20的事情了

Level 6 AKS

上面要么慢,要么基于概率(或有限范围的历史答案)

$(x-a)^n = x^n-a (\bmod n), gcd(\forall a,n) = 1$,

n是质数时,根据2项式定理显然,

如果n 合数, $(X-a)^n - X^n-a $ 看成 $X$为变量, 其余部分为系数

那么 $X^i$的系数为$C(n,i) a^{n-i}$

把合数表示成其构成的质数相乘$n = p_1^m \cdots$, 那么2项式中$C(n,p_1) \neq 0 (\bmod p_1^m )$,因为$C(n,p_1) = \frac{n!}{p!(n-p)!}= \frac{(n-p+1)…(n)}{1…p}$ 分子只有n是p的倍数,分母只有p是n的倍数

所以$C(n,p_1)a^{n-p_1} \neq 0 (\bmod p_1^m)$ (因为a和n互质,不会包含因子$p_1$) 也说明 $C(n,p_1)a^{n-p_1} \neq 0 (\bmod n)$

注意我们$(X-a)^n - X^n-a $ 看成 $X$为变量, 其余部分为系数,后就是一个 f(X)的表达式, 有系数非零

例如 n = 4

$(x+a)^4 = x^4+C(4,1)ax^3+C(4,2)a^2x^2+C(4,3)a^3x+a^4 = x^4 + C(4,2)a^2x^2+a^4 (\bmod 4)$

$(x+a)^4 -(x^4+a) = 6a^2x^2+a^4-a (\bmod 4)$

想存在$x,gcd(a,4) = 1$,使得 $6a^2x^2+a^4-a \neq (\bmod 4)$, 我们的确能找到 a=1,x=1,有些地方说,

既然是关于x有系数的表达式,或者项数差异,说明它在模意义下非恒为零??? 这是模运算的什么定理吗

小证? 写成 系数乘上 x的从1取到n的范得Vandermonde 矩阵, 系数和结果都对n取mod,所以只能全是零才有解? 甚至如果希望陪结果,可以反向推?脑补证明的感觉可能有误. 比如如果要按mod算的话,那vandermonde矩阵的n列就全是0和部分 位置XD

但是要中间的系数全是n的倍数上面已经证明

有的地方说主要是等式而不是值?

TODO

其它

$(n-1)! = -1 (\bmod p)$ 当且仅当 p为质数,陶哲轩的书上也有这个

效率更差, 在一些特定场景有用

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef uint64_t ull;
#define rep(i,a,n) for (ll i=a;i<n;i++)

// overflow
ll quick_p(__int128_t b, ll p,const ll mod){
  __int128_t r = 1;
  while(p){
    if(p%2)(r*=b)%=mod;
    (b*=b)%=mod;
    p/=2;
  }
  return r%mod;
}

ll mr(ll base,ll v){
  if(base > v)return true;
  ll startp = v-1;
  while(startp%2 == 0)startp>>=1;
  ll p = startp;
  __int128_t r = quick_p(base,p,v);
  while(p != v-1){
    if(r == v-1)return true;
    if(r == 1)return p == startp;
    p*=2;
    // overflow
    (r*=r)%=v;
  }
  return false;
}

bool is_prime_64(ll v){
  if(v < 2)return false;
  if(v < 4)return true;
  // %6 = 1 or 5
  if((v % 6) % 4 != 1)return false;
  ll test_g[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022};
  rep(i,0,7){
    if(!(mr(test_g[i],v)))return false;
  }
  return true;
}

bool isp(ll v){
  rep(i,2,v){
    if(i*i>v)return true;
    if(v%i == 0)return false;
  }
  return true;
}

int main(){
  ll ans = 0;
  rep(i,2,50'000'000+1){
    if(i % 500'000 == 0){
      printf("i:%lld\n",i);
    }
    ans += is_prime_64(2*i*i-1);
  }
  printf("ans:%lld\n",ans);
  return 0;
}

Ref

https://brilliant.org/wiki/prime-testing/

https://primes.utm.edu/prove/prove2_3.html

2002 primes is P

http://www.cs.tau.ac.il/~amnon/Classes/2019-Derandomization/Lectures/Lecture7-AKS-All.pdf

https://en.wikipedia.org/wiki/AKS_primality_test

系数与完整表达式的疑问也有人问过

http://blog.sciencenet.cn/blog-3224443-1115018.html

题目

每次新增一个塑料袋

可以把偶数个之前的塑料袋放入这个新增的塑料袋

求n次后 的状态数

f(4)=5

f(8)=1385

求f(24680)%1020202009

打表

显然 OEIS A000111

DP(n^3)

dp[i][j] = 第i步,剩余j个袋子

dp[i][j] = dp[i-1][j-1+d] * c(j-1+d,d), d = 0..i-j, step = 2,(i-j)%2=0

试图优化,内部运算也没有摆脱O(n^3), n在2w 搞不了

DP(n^2)

f(n+1) =

0个放袋子n+1里c(n,0) * f(0) * f(n)

2个放袋子n+1里c(n,2) * f(2) * f(n-2) ,我们 不关心 这2个 和 这n-2个 内部的关系,而只关心谁在n+1袋子里,谁在袋子外,有了c(n,2), 于是 这2个和n-2互不干扰,各自的结构只与其内部相关,所以方案数也就是f(2)和f(n-2)

2k个放袋子n+1里c(n,2k) * f(2k) * f(n-2k) ,我们 不关心 这2k个 和 这n-2k个 内部的关系,而只关心谁在n+1袋子里,谁在袋子外,有了c(n,2k), 于是 这2k个和n-2k互不干扰,各自的结构只与其内部相关,所以方案数也就是f(2k)和f(n-2k)

f(n+1) = sum( c(n,2k) * f(2k) * f(n-2k) ), k = 0..floor(n/2)

令g(n) = f(n)/(n!), 其实上面dp我也用了这个技巧,但是还是有剩余的阶乘,无法简化到递推

g(n+1) * (n+1) = sum( g(2k) * g(n-2k)), k = 0..floor(n/2)

Ref

https://hiragn.hatenablog.com/entry/2020/12/01/031135

https://gist.github.com/hrgnz/10c4c3afb5a332c8ad3428327e6459c0#file-pe709a-nb

完美洗牌

AAAAAAAAAAAAAAAAABBBBBBBBBBBBBBBBB

变为

ABABABABABABABABABABABABABABABABAB

次数

52张需要8次变回原型

需要8次的方案的张数和为412

问

需要60次的方案的张数和是多少

解

首先 估计已经到2的60次方左右的搜索范围,暴力不可取

递推

f(i,n) = (i%n)*2+int(i>=n)

然后我们分开写

f(i,n) = 2i (i < n)

f(i,n) = 2(i-n)+1 (i >= n)

变形

f(i,n) = 2i - (2n-1) (i >= n)

合并

f(i,n) = (2i)%(2n-1) 这里就很神奇了

相当于每个值的变化都是在乘2,又膜运算有合并性, 所以任何值a的k次后的值为 (a * 2^k) % (2n-1)

所以要所有都会到原位其实就是 任意 a, a = (a * 2^k) % (2n-1)

a取1,也是最大循环节

题目变成 $(2^60-1)%(2n-1) = 0$ 且 $(2^(<60)-1)%(2n-1) != 0$

质因数分解一下就好

总结

关键在 int(i>=n)的那个表达式转换成纯的模运算表达式,后面都好做了