部分知识合集
数学
数论
组合优化
线性规划的单纯形法
具体内容找个时间专门写翻译:先咕着,这里给出 UOJ 的 97pts 的核心代码:
#define rep( i, a, b ) for( int i = (a) ; i <= (b) ; i ++ )
#define per( i, a, b ) for( int i = (a) ; i >= (b) ; i -- )
typedef int DOUBLE;
typedef vector<int> VI;
typedef vector<double> VD;
typedef vector<VD> VVD;
const long double EPS = 1e-9;
struct LP
{
VVD A; VI bId, nId;
int n, m;
// 其中 n 是非基变量个数, m 是基变量个数
LP() { n = m = 0; }
LP( int N, int M ):
A( M + 2, VD( N + 2 ) ), bId( M + 1 ), nId( N + 1 ), n( N ), m( M )
{
rep( i, 1, n ) nId[i] = i;
rep( j, 1, m ) bId[j] = n + j;
}
void Pivot( const int r, const int s )
// 将 x_{n+r} 换做 x_{s}
{
DOUBLE t = A[r][s];
rep( i, 1, m + 1 ) if( i ^ r )
rep( j, 1, n + 1 ) if( j ^ s )
A[i][j] -= A[r][j] / t * A[i][s];
rep( i, 1, m + 1 ) if( i ^ r ) A[i][s] /= - t;
rep( j, 1, n + 1 ) if( j ^ s ) A[r][j] /= t;
A[r][s] = ( DOUBLE ) ( 1 ) / t, swap( bId[r], nId[s] );
}
DOUBLE Result() const { return A[m + 1][n + 1]; }
};
namespace LPSolver
{
int Dcmp( const DOUBLE &a, const DOUBLE &b = 0 )
{ return ( DOUBLE ) fabs( a - b ) < EPS ? 0 : ( a < b ? -1 : 1 ); }
bool Simplex( LP &lp )
{
for( int r, s, t ; ; )
{
r = s = -1;
rep( j, 1, lp.n )
if( Dcmp( lp.A[lp.m + 1][j] ) < 0 && ( s == -1 ||
( t = Dcmp( lp.A[lp.m + 1][j], lp.A[lp.m + 1][s] ) ) < 0 ||
( t == 0 && lp.nId[j] < lp.nId[s] ) ) )
s = j;
if( s == -1 ) break;
rep( i, 1, lp.m )
if( Dcmp( lp.A[i][s] ) > 0 && ( r == -1 ||
( t = Dcmp( lp.A[i][lp.n + 1] / lp.A[i][s], lp.A[r][lp.n + 1] / lp.A[r][s] ) ) < 0 ||
( t == 0 && lp.bId[i] < lp.bId[r] ) ) )
r = i;
if( r == -1 ) return false;
lp.Pivot( r, s );
}
return true;
}
bool Init( LP &lp )
{
VD tmp( lp.n + 1 );
bool flg = true;
rep( i, 1, lp.m )
if( Dcmp( lp.A[i][lp.n + 1] ) < 0 )
{ flg = false; break; }
if( flg ) return true;
LP aux( lp.n + 1, lp.m );
rep( i, 1, aux.m )
{
rep( j, 1, lp.n )
aux.A[i][j] = lp.A[i][j];
aux.A[i][lp.n + 1] = - 1;
aux.A[i][lp.n + 2] = lp.A[i][lp.n + 1];
}
rep( i, 1, lp.n ) aux.nId[i] = lp.nId[i]; aux.nId[lp.n + 1] = 0;
rep( i, 1, lp.m ) aux.bId[i] = lp.bId[i];
aux.A[lp.m + 1][lp.n + 1] = 1;
int idx = -1;
rep( i, 1, lp.m )
if( idx == -1 || Dcmp( aux.A[i][aux.n + 1], aux.A[idx][aux.n + 1] ) < 0 )
idx = i;
aux.Pivot( idx, aux.n );
if( ! Simplex( aux ) ) return false;
if( Dcmp( aux.A[aux.m + 1][aux.n + 1] ) < 0 ) return false;
rep( i, 1, aux.m ) if( ! aux.bId[i] ) { aux.Pivot( 1, i ); break; }
int tot;
rep( i, 1, lp.m )
{
int tot = 0;
rep( j, 1, aux.n )
if( aux.nId[j] )
lp.A[i][++ tot] = aux.A[i][j];
lp.A[i][++ tot] = aux.A[i][aux.n + 1];
}
tot = 0;
rep( j, 1, aux.n )
if( aux.nId[j] )
lp.nId[++ tot] = aux.nId[j];
rep( i, 1, aux.m ) lp.bId[i] = aux.bId[i];
rep( j, 1, lp.n ) tmp[j] = lp.A[lp.m + 1][j], lp.A[lp.m + 1][j] = 0;
rep( i, 1, lp.m )
if( lp.bId[i] <= lp.n )
rep( j, 1, lp.n + 1 )
lp.A[lp.m + 1][j] -= tmp[lp.bId[i]] * lp.A[i][j];
rep( i, 1, lp.n )
if( lp.nId[i] <= lp.n )
lp.A[lp.m + 1][i] += tmp[lp.nId[i]];
return true;
}
int Solve( LP &lp )
{
if( ! Init( lp ) ) return 1;
if( ! Simplex( lp ) ) return 2;
return 0;
}
}
计算几何
点和线
-
求两条直线的交点:
解析几何的方法还是略显复杂,向量的方法比较给力。
设 \(l_1=\vec p_1+n\vec d_1,l_2=\vec p_2+m\vec d_2\),则不难列出方程 \(\vec p_1+n\vec d_1=\vec p_2+m\vec d_2\),大量跳步之后可以整理得到:
\[\begin{bmatrix} x_{d_1}&x_{d_2}\\ y_{d_1}&y_{d_2} \end{bmatrix} \begin{bmatrix} n\\m \end{bmatrix}= \begin{bmatrix} x_{p_2}-x_{p_1}\\ y_{p_2}-y_{p_1} \end{bmatrix} \]使用克莱姆法则直接解得 \(n=\frac{(\vec p_2-\vec p_1)\times \vec d_2}{\vec d_1\times \vec d_2}\),代入算出点即可。
三角形
-
单纯的面积公式:
-
海伦公式 \(S=\sqrt{p(p-a)(p-b)(p-c)}, p=\frac{a+b+c}{2}\)。
其中 \(a,b,c\) 为三角形三边长,\(p\) 为半周长,\(S\) 为面积。
一种用不到三角函数的简单证明方法就是画高线直接算。
-
经典的 \(S=\frac{1}{2}ab\sin C\),不多做解释。
-
-
三角形的若干心:
-
重心,三条中线的交点,有如下性质:
-
重心为任意中线的一个三等分点,证明请使用面积法。
-
重心为三角形内到三个顶点平方和最小的点,证明不知道。
-
重心为三角形内到三边距离最大的点,证明不知道。
-
-
内心,三条角平分线的交点,也是三角形内切圆的圆心,因此它满足到三条边的距离相等。
此外,顺便说到三角形的内接圆半径,为 \(\frac{2S}{C}\),其中 \(S\) 为面积,\(C\) 为周长。这个可以直接算面积证明,
但是很容易搞忘。 -
外心,三条中垂线的交点,也是三角形外接圆的圆心,因此它满足到三个顶点的距离相等。
此外,顺便说到正弦定理。设 \(R\) 为外接圆半径,则有 \(\frac{a}{\sin A}=\frac{b}{\sin B}=\frac{c}{\sin C}=2R\)。证明可以画外接圆,连接一端点和圆心并延长,之后看图说话。
相应的,根据外接圆半径可以导出 \(S=\frac{abc}{4R}\)。
-
垂心,三条高线的交点。OI 中还没有什么亮眼的性质。
-
旁心,任意两条外角平分线的交点都可以叫旁心,因此有三个。旁心对应的是三角形的三个“切一条边并且切另外两条边的延长线”的圆的圆心。
-
多边形
-
Pick 定理:
网格图上的简单(开)多边形,如果在它的内部有 \(a\) 个整点,边上有 \(b\) 个整点,则它的面积 \(S=a+\frac{b}{2}-1\),证明还不会。
图论
最小树形图
问题直接顾名思义,略去。该问题常用朱刘算法,下面这张图很好地展示了它的流程,虽然确实丑得有点离谱了:
本质上朱刘算法就是一个反悔贪心:
-
初始,贪心选择最小入边。
-
没有出现环,我们必然找到了树形图。
-
否则,出现了环,我们必须要拆掉环上的一条边,这个过程通过修改环上所有点的入边边权(减去终点的最小边权)实现,这样在之后选中某条指向环的边的时候,我们就相当于用新边替换了旧边,使得树形图合法;可以发现这个过程一定会发生。
-
最后,收缩环,迭代。
以下是这个算法的 \(O(|V||E|)\) 实现:
int inE[MAXN], pre[MAXN], id[MAXN], vis[MAXN];
int fr[MAXM], to[MAXM], wei[MAXM];
int N, M; // N: 点数; M: 边数
LL ZhuLiu( int rt )
{
LL ret = 0;
int n = N, cnt = 0;
while( true )
{
cnt = 0;
rep( i, 1, n )
vis[i] = id[i] = 0, inE[i] = INF, pre[i] = -1;
rep( i, 1, M )
if( fr[i] ^ to[i] && wei[i] < inE[to[i]] )
inE[to[i]] = wei[i], pre[to[i]] = i;
// 寻找最小入边;注意我们必须忽略 " 自环 " ,因为这里的自环其实也有可能是缩环之前的环边。
rep( i, 1, n )
{
if( i == rt ) continue;
if( inE[i] == INF ) return -1;
else ret += inE[i];
}
rep( i, 1, n )
{
if( i == rt ) continue; int cur = i;
for( ; vis[cur] ^ i && ! id[cur] && cur ^ rt ; cur = fr[pre[cur]] ) vis[cur] = i;
if( ! id[cur] && cur ^ rt )
{
id[cur] = ++ cnt;
for( int u = fr[pre[cur]] ; u ^ cur ; u = fr[pre[u]] )
id[u] = cnt;
//注意我们可能并非一开始就在环上,而是走着走着掉进了一个环里
}
}
//找环,缩环
if( cnt == 0 ) break;
rep( i, 1, n ) if( ! id[i] ) id[i] = ++ cnt;
rep( i, 1, M )
{
int u = fr[i], v = to[i];
fr[i] = id[u], to[i] = id[v], wei[i] -= inE[v];
//保证忽略环边,因此可以对所有边进行修改
//对于不指向环的边,最小边权已经加入答案,且之后的最小边权必然为 0
}
// 修改边权,重新构图
rt = id[rt], n = cnt;
}
return ret;
}
朱刘算法还有其更加高效的实现。具体而言,我们考虑朱刘算法的过程:
- 对于一个点,找出最小入边
- 如果出现了环,合并环上的所有入边,接着回到 1
假如我们想用数据结构维护,不妨分析上述材料:
-
" 最小入边 " 体现了我们需要使用一种可以快速查询最小值的数据结构,堆应该是不错的选择。
-
" 环 " 说明了我们需要快速判断连通性,并查集可以很方便地解决。
-
" 合并入边 " 说明了我们需要一种可以快速合并的数据结构。
-
总之,我们可以使用可并堆维护每个点的出边集合,并查集判断连通性。
如果出现了环,我们就暴力把环上的出边收缩起来(总收缩次数是 \(O(n)\) 的)。在算法过程中我们需要始终维护连通性。
具体可以参考代码:
struct DSU
{
int fa[MAXN];
void MakeSet( const int siz ) { rep( i, 1, siz ) fa[i] = i; }
int FindSet( const int u ) { return fa[u] = ( fa[u] == u ? u : FindSet( fa[u] ) ); }
void UnionSet( const int u, const int v ) { fa[FindSet( u )] = FindSet( v ); }
}cir, link;
int fr[MAXM], to[MAXM];
int lch[MAXM], rch[MAXM], dist[MAXM], w[MAXM], tag[MAXM];
int rt[MAXN], pre[MAXN];
int N, M, R, ans = 0;
void Add( const int x, const int v ) { if( x ) w[x] += v, tag[x] += v; }
void Normalize( const int x ) { Add( lch[x], tag[x] ), Add( rch[x], tag[x] ), tag[x] = 0; }
int Merge( int u, int v )
{
if( ! u || ! v ) return u | v;
if( w[u] > w[v] ) swapp( u, v );
Normalize( u ); rch[u] = Merge( rch[u], v );
if( dist[lch[u]] < dist[rch[u]] ) swapp( lch[u], rch[u] );
dist[u] = dist[rch[u]] + 1; return u;
}
int Pop( const int u )
{
Normalize( u );
return Merge( lch[u], rch[u] );
}
bool Get( const int u )
{
for( ; cir.FindSet( fr[rt[u]] ) == cir.FindSet( u ) ;
rt[u] = Pop( rt[u] ) );
//类似于朴素朱刘,这里不能取到环内边
//一定注意 判连通性 和 判环内 需要两个并查集,因为在环的最后一条边连上之前环其实已经连通了
if( ! rt[u] ) return false;
ans += w[pre[u] = rt[u]]; return true;
}
int main()
{
rep( i, 1, M )
rt[to[i]] = Merge( rt[to[i]], i );
cir.MakeSet( N ), link.MakeSet( N );
rep( i, 1, N )
if( i ^ R ) //每个点可以顺序寻找入边
{
if( ! Get( i ) ) return puts( "-1" ), 0;
while( link.FindSet( fr[pre[i]] ) == link.FindSet( i ) )
//假如出现环,则需要缩环
{
Add( rt[i], - w[pre[i]] );
for( int u = cir.FindSet( fr[pre[i]] ) ; u ^ i ; u = cir.FindSet( fr[pre[u]] ) )
{
Add( rt[u], - w[pre[u]] );
rt[i] = Merge( rt[i], rt[u] );
cir.UnionSet( u, i );
}
//暴力打标记合并,利用环的并查集跳转
if( ! Get( i ) ) return puts( "-1" ), 0;
}
link.UnionSet( i, fr[pre[i]] );
}
write( ans ), putchar( '\n' );
return 0;
}
一般图最大匹配
常用带花树算法。主要思想是处理一般图上特殊的奇环,通过 balabala (这只鸽子决定期末考完了再补)从而说明可以将奇环缩成一个点,然后就可以像二分图匹配一样直接跑增广路增广。
具体可以参考代码:
int q[MAXN], h, t;
int head[MAXN], fa[MAXN], pre[MAXN], match[MAXN], vis[MAXN], travel[MAXN];
int N, M, cnt, tims = 0;
int LCA( int u, int v )
{
u = fa[u], v = fa[v], tims ++;
for( ; travel[u] ^ tims ; u ^= v ^= u ^= v )
if( u ) travel[u] = tims, u = fa[pre[match[u]]];
return u;
}
void Blossom( int u, int v, int lca )
{
while( fa[u] ^ lca )
{
pre[u] = v, v = match[u];
if( vis[v] == 1 ) vis[v] = 2, q[++ t] = v;
fa[u] = fa[v] = lca, u = pre[v];
}
}
bool Augment( const int sta )
{
int u, v; h = 1, t = 0;
rep( i, 1, N ) fa[i] = i, pre[i] = vis[i] = 0;
vis[q[++ t] = sta] = 2;
while( h <= t )
{
u = q[h ++];
for( int i = head[u] ; i ; i = Graph[i].nxt )
if( ! vis[v = Graph[i].to] )
{
vis[v] = 1, pre[v] = u;
if( ! match[v] )
for( int cur = v ; ; )
{
int nxt = match[pre[cur]];
match[cur] = pre[cur], match[pre[cur]] = cur;
if( pre[cur] == sta ) return true;
cur = nxt;
}
else vis[q[++ t] = match[v]] = 2;
}
else if( vis[v] == 2 && fa[u] ^ fa[v] )
{
int lca = LCA( u, v );
Blossom( u, v, lca ), Blossom( v, u, lca );
}
}
return false;
}
int main()
{
int ans = 0;
rep( i, 1, N )
ans += ! match[i] && Augment( i );
}
数据结构
字符串
扩展 KMP
扩展 KMP 有两种情况:
- 给定字符串 \(S\),求 \(S\) 的每个后缀和 \(S\) 的 LCP 长度。
- 给定字符串 \(S,T\),求 \(T\) 的每个后缀和 \(S\) 的 LCP 长度。
当然,两种情况的核心算法都是一样的。
我们先看一下情况 1。一种暴力的算法是枚举后缀并且暴力扫描,这个算法最坏当然是 \(O(|S|^2)\) 的。
优化的想法比较朴素,就是尝试充分利用已知信息来优化。设 \(z_k\) 为后缀 \(S[k:|S|]\) 的答案。考虑在求 \(S[k:|S|]\) 的时候,我们已经计算好了 \(z_1,z_2,\dots,z_{k-1}\),则已有的信息是——对于 \(1\le i<k\),有 \(S[1,z_i]=S[i,i+z_i-1]\)。
- 如果没有一个 \([i,i+z_i)\) 覆盖到了 \(k\),我们似乎看不出什么比较好的方法来优化,干脆直接暴力解决;
- 否则,我们不妨设 \(k\in [i_0,i_0+z_{i_0})\),这个信息暗示我们 \(S[k,i_0+z_{i_0}-1]=S[k-i_0+1,z_{i_0}]\),而这样我们就可以知道 \(z_k\ge \min\{z_{k-i_0+1},i_0+z_{i_0}-k\}\)(\(i_0+z_{i_0}\) 及之后的字符串长什么样还是未知的)。这是一个确定的下界,我们可以在这个下界的基础上再继续扩展。
显然,如果 \(i+z_i<j+z_j\),我们就只需要考虑 \([j,j+z_j)\) 而没有必要考虑 \([i,i+z_i)\)。于是我们只需要保留前缀中 \(i+z_i\) 的一个区间即可。
设 \(r = \max\{i+z_i|1\le i\le k\}\),则这个优化过后的算法可以保证每次暴力扫描时都一定会导致 \(r\) 增大。由于 \(r\le |S|\),所以该算法的复杂度为 \(O(|S|)\)。
再看一下情况 2。实际上这和情况 1 非常类似,我们也可以用同样的思路。
假设我们先已对于 \(S\) 执行了情况 1 的算法,得到了 \(z\)。
设 \(z'_{k}\) 为 \(T[k,|T|]\) 与 \(S\) 的 LCP 的长度。则在计算 \(z'_k\) 时,已知的信息是对于 \(1\le i<k\),都有 \(S[1,z'_i]=T[i,i+z'_i-1]\)(以及 \(S\) 的 \(z\) 信息)。
- 如果没有一个 \([i,i+z'_i)\) 覆盖到了 \(k\),类似略过。
- 否则,我们不妨设 \(k\in [i_0,i_0+z'_{i_0})\),得到 \(T[k,i_0+z'_{i_0}-1]=S[k-i_0+1,z'_{i_0}]\),这样我们就知道 \(z'_k\ge \min\{z_{k-i_0+1},i_0+z'_{i_0}-k\}\)。
之后的过程都类似,于是我们得到了一个 \(O(|T|)\) 的优秀算法。
洛谷模板题,参考代码:
template<typename _T>
_T MIN( const _T a, const _T b ) {
return a < b ? a : b;
}
int z[MAXN], p[MAXN];
char S[MAXN], T[MAXN];
int N, M;
int main() {
scanf( "%s%s", S, T );
N = strlen( S ), M = strlen( T );
z[0] = M;
for( int i = 1, l = 0, r = 0 ; i < M ; i ++ ) {
if( i < r ) z[i] = MIN( r - i, z[i - l] );
for( int &k = z[i] ; k + i < M && T[k] == T[i + k] ; k ++ );
if( i + z[i] > r ) l = i, r = i + z[i];
}
for( int i = 0, l = 0, r = 0 ; i < N ; i ++ ) {
if( i < r ) p[i] = MIN( r - i, z[i - l] );
for( int &k = p[i] ; k < M && i + k < N && S[i + k] == T[k] ; k ++ );
if( i + p[i] > r ) l = i, r = i + p[i];
}
return 0;
}
其他
模拟退火
具体实现细节就不写了,提醒几点在下面:
-
退火其实还是类似于爬山,只不过是增加了跳出的概率。因此一次退火的过程中需要使用过程解进行迭代,如果用最优解迭代其实是没啥用的。
-
概率函数: \(e^{\frac{\Delta}{T}}\) ,注意 \(\Delta\) 的正负。个人认为一个简单的判别方法是:作为一个概率值,函数值应该落在 \([0,1]\) 里面,所以应该有 \(\Delta<0\) 。
-
生成的新解需要在当前解的周围搜出一个点,注意新解的范围不要太小。
-
毕竟是个概率算法,正式场合少用,
除非你很会调参。 -
有的时候模拟退火没有爬山好用,尤其是 OI 的传统题,除非你很会调参