Loading docs/math/linear-programming.md +191 −2 Original line number Diff line number Diff line Loading @@ -25,13 +25,13 @@ 我们可以使用数学语言来描述它： 设 $x_1$ 为花费在修路广告上的钱（千美元） 设 $x_1$ 为花费在修路广告上的钱（千美元） 设 $x_2$ 为花费在枪支管制广告上的钱（千美元） 设 $x_3$ 为花费在农场补贴广告上的钱（千美元） 设 $x_4$ 为花费在汽油税广告上的钱（千美元） 设 $x_4$ 为花费在汽油税广告上的钱（千美元） 那么我们可以将在市区获得至少 50000 张市区选票表述为 Loading Loading @@ -106,3 +106,192 @@ 显而易见的，打了蓝色斜线的区域为三块区域的交集，这就是这个线性规划的所有可行解。因为题目中说明，需要最小化 x 和 y，观察图像可知，点（2.3）为可行解中 x 和 y 最小的一个。因此， $x_{min}=2,y_{min}=3$ 把求解线性规划的图解法总结起来，就是先在坐标系中作出所有的约束条件，然后作出需要求极值的线性函数的定义域。定义域与约束条件的交集就是这个线性规划的解集，而所需求的极值由观察可以轻易得出。 ### 单纯形法 > $$ > \forall k,\sum_{i}a_{k,i}x_i\le b_{k}\\ > \forall i,x_i\ge 0 > $$ > > 求 > > $$ > \max(\sum_{i}c_{i}x_i) > $$ 我们考虑一个下标从 $0$ 开始的矩阵 $A_{0\sim m,0\sim n+m}$: $$ A= \begin{bmatrix} 0 & c_1 & \ldots & c_n & 0 & 0 & \ldots & 0 \\ b_1 & a_{1,1} & \ldots & a_{1,n} & 1 & 0 & \dots & 0 \\ b_2 & a_{2,1} & \ldots & a_{2,n} & 0 & 1 & \dots & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \ddots & \vdots \\ b_m & a_{m,1} & \ldots & a_{m,n} & 0 & 0 & \dots & 1 \end{bmatrix} $$ 其中 $A_{0,0}$ 为这个状态和初始时答案之差，新增的 $m$ 行相当于新建 $m$ 个变量，编号为 $x_{n+1\sim n+m}$，要求它们不小于 $0$ 于是 $$ \forall i\in[1,m]\\ \sum_{1\le j\le n+m}a_{i,j}x_j=b_i $$ 显然，要先求出一套可行的方案 如果 $$ x_i=\begin{cases}0 & i\in[1,n]\\b_{i-n} & i\in[n+1,n+m]\end{cases} $$ 那么显然满足，但是由于 $b_i$ 不一定不小于 $0$，所以这样可能导致 $x_{i}$ 小于 $0$ 然后，介绍一个魔法，通过矩阵变化改变 $b$ 的值来求合法解/最优解 #### 转轴 ~~之所以叫这个名称，是因为可行域是一个超立方体，这个魔法是经过一条轴，但是我不讲~~ 这个魔法是指： 选择两个数 $u,v$，考虑将 $x_v$ 变为非 $0$ 数而 $x_{u+n}$ 变为 $0$，也就是让 $x_v$ 来代替 $x_{u+n}$ 作为等于 $b_u$ 的位置。 如果要实现，就要运用矩阵变换。 我们可以先让第 $u$ 行全部除一个 $A_{u,v}$，然后把其它所有行 $i$ 的位置 $v$ 变为 $0$（即使 $A_{i,v}=0$），可以直接把第 $i$ 行减去 $A_{i,v}$ 倍第 $u$ 行 然后就可以了！ 因为（不算 $u+n,v$ 列的话）后面 $m$ 列只有 $A_{i,i+n}=1$，前面 $n$ 列都是 $0$，所以矩阵变化不会影响其它的 $x$。 当然，最好把第 $u+n,v$ 行换一下，这样会方便很多 （还有，因为矩阵后 $m$ 列值的特殊性，一般不存这 $m$ 列） 顺便一提，这对答案的贡献为： $$ \dfrac{A_{0,v}A_{u,0}}{A_{u,v}} $$ ```cpp int id[N<<1]; double a[N][N]; void pivot(int u,int v) { double k=a[u][v]; swap(id[u+n],id[v]); a[u][v]=1.;//已经知道第v列要被消掉了，就可以直接变为u+n列了（尽管要用一遍，就临时存一下） for(int i=0;i<=n;i++) a[u][i]/=k; for(int i=0;i<=m;i++) if(i!=u&&fabs(a[i][v])>eps) { k=a[i][v]; a[i][v]=0;//同上 for(int j=0;j<=n;j++) a[i][j]-=k*a[u][i]; } } ``` 知道了基本的变（魔）化（法），就来看怎么求一个合法解 #### 初始化 引用一下 Candy? 的博客 > 算法导论上有一个辅助线性规划的做法 > > 但我发现好多人都用了**随机初始化**的黑科技 > > 在所有 $a_{i,0}<0$ 的约束中随机选一个作为 $u$，再随机选一个 $a_{u,v}<0$ 作为 $v$，然后 `pivot(u,v)` 后 $a_{i,0}$ 就变正了... 但这可能出现循环，然后不断损精度，然后 WA/TLE，总之很玄学就对了 但是有「优先前面/后面」的选择倾向再加以小随机效果似乎不错 ```cpp int init() { while (1) { int u=0,v=0; for(int i=1;i<=m;i++) if(a[i][0]<-eps&&(!u||rand()%3==0)) u=i; if(!u)//说明这个解是可行解 return 1; for(int i=1;i<=n;i++) if(a[u][i]<-eps&&(!v||rand()%3==0)) v=i; if(!v) {//说明无解 printf("Infeasible\n"); return 0; } pivot(u,v); } } ``` #### 求解 每次找到一个 $v$ 使得 $A_{0,v}>0$，再找到一个 $u$ 使得 $A_{u,v}>0$，然后 `pivot(u,v)` 即可。 ##### 退化 在旋转的过程中，可能会存在保持答案不变的情况，这种现象称为退化。 ##### Bland 规则 在选择替入变量和替出变量的时候，我们总是选择满足条件的下标最小值。 - 替入变量 $v$ ：目标条件中，系数为正数（如果求最小值就是负数）的第一个作为替入变量 - 替出变量 $u$ ：对所有的约束条件中，选择对 $u$ 约束最紧的第一个（也就是 $\dfrac{A_{u,0}}{A_{u,v}}$ 最小的那个） ```cpp int simplex() { while (1) { int u=0,v=0; double mi=inf; for(int i=1;i<=n;i++) if(a[0][i]>eps) { v=i; return; } if(!v) return 1; for(int i=1;i<=m;i++) if(a[i][v]>eps&&a[i][0]/a[i][v]<mi) {//bland法则 mi=a[i][0]/a[i][v]; u=i; } if(!u) { //没有约束 printf("Unbounded\n"); return 0; } pivot(u,v); } } ``` #### 其它说明 ##### 为什么这能对 每次矩阵变化后的子问题是等价于初始问题的（易证） 并且这类问题有一个特殊性：它的可行域是凸的 意味着不会有除最优解外的局部最优解出现 ##### 这能跑多大 ~~心有多大，梦想有多大，它就能跑多大~~ 对于这种玄学算法，鬼才知道它能跑多大。 不过一般来说， $nm\le10^6$ 是~~不~~可行的 Loading
docs/math/linear-programming.md +191 −2 Original line number Diff line number Diff line Loading @@ -25,13 +25,13 @@ 我们可以使用数学语言来描述它： 设 $x_1$ 为花费在修路广告上的钱（千美元） 设 $x_1$ 为花费在修路广告上的钱（千美元） 设 $x_2$ 为花费在枪支管制广告上的钱（千美元） 设 $x_3$ 为花费在农场补贴广告上的钱（千美元） 设 $x_4$ 为花费在汽油税广告上的钱（千美元） 设 $x_4$ 为花费在汽油税广告上的钱（千美元） 那么我们可以将在市区获得至少 50000 张市区选票表述为 Loading Loading @@ -106,3 +106,192 @@ 显而易见的，打了蓝色斜线的区域为三块区域的交集，这就是这个线性规划的所有可行解。因为题目中说明，需要最小化 x 和 y，观察图像可知，点（2.3）为可行解中 x 和 y 最小的一个。因此， $x_{min}=2,y_{min}=3$ 把求解线性规划的图解法总结起来，就是先在坐标系中作出所有的约束条件，然后作出需要求极值的线性函数的定义域。定义域与约束条件的交集就是这个线性规划的解集，而所需求的极值由观察可以轻易得出。 ### 单纯形法 > $$ > \forall k,\sum_{i}a_{k,i}x_i\le b_{k}\\ > \forall i,x_i\ge 0 > $$ > > 求 > > $$ > \max(\sum_{i}c_{i}x_i) > $$ 我们考虑一个下标从 $0$ 开始的矩阵 $A_{0\sim m,0\sim n+m}$: $$ A= \begin{bmatrix} 0 & c_1 & \ldots & c_n & 0 & 0 & \ldots & 0 \\ b_1 & a_{1,1} & \ldots & a_{1,n} & 1 & 0 & \dots & 0 \\ b_2 & a_{2,1} & \ldots & a_{2,n} & 0 & 1 & \dots & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots & \vdots & \ddots & \vdots \\ b_m & a_{m,1} & \ldots & a_{m,n} & 0 & 0 & \dots & 1 \end{bmatrix} $$ 其中 $A_{0,0}$ 为这个状态和初始时答案之差，新增的 $m$ 行相当于新建 $m$ 个变量，编号为 $x_{n+1\sim n+m}$，要求它们不小于 $0$ 于是 $$ \forall i\in[1,m]\\ \sum_{1\le j\le n+m}a_{i,j}x_j=b_i $$ 显然，要先求出一套可行的方案 如果 $$ x_i=\begin{cases}0 & i\in[1,n]\\b_{i-n} & i\in[n+1,n+m]\end{cases} $$ 那么显然满足，但是由于 $b_i$ 不一定不小于 $0$，所以这样可能导致 $x_{i}$ 小于 $0$ 然后，介绍一个魔法，通过矩阵变化改变 $b$ 的值来求合法解/最优解 #### 转轴 ~~之所以叫这个名称，是因为可行域是一个超立方体，这个魔法是经过一条轴，但是我不讲~~ 这个魔法是指： 选择两个数 $u,v$，考虑将 $x_v$ 变为非 $0$ 数而 $x_{u+n}$ 变为 $0$，也就是让 $x_v$ 来代替 $x_{u+n}$ 作为等于 $b_u$ 的位置。 如果要实现，就要运用矩阵变换。 我们可以先让第 $u$ 行全部除一个 $A_{u,v}$，然后把其它所有行 $i$ 的位置 $v$ 变为 $0$（即使 $A_{i,v}=0$），可以直接把第 $i$ 行减去 $A_{i,v}$ 倍第 $u$ 行 然后就可以了！ 因为（不算 $u+n,v$ 列的话）后面 $m$ 列只有 $A_{i,i+n}=1$，前面 $n$ 列都是 $0$，所以矩阵变化不会影响其它的 $x$。 当然，最好把第 $u+n,v$ 行换一下，这样会方便很多 （还有，因为矩阵后 $m$ 列值的特殊性，一般不存这 $m$ 列） 顺便一提，这对答案的贡献为： $$ \dfrac{A_{0,v}A_{u,0}}{A_{u,v}} $$ ```cpp int id[N<<1]; double a[N][N]; void pivot(int u,int v) { double k=a[u][v]; swap(id[u+n],id[v]); a[u][v]=1.;//已经知道第v列要被消掉了，就可以直接变为u+n列了（尽管要用一遍，就临时存一下） for(int i=0;i<=n;i++) a[u][i]/=k; for(int i=0;i<=m;i++) if(i!=u&&fabs(a[i][v])>eps) { k=a[i][v]; a[i][v]=0;//同上 for(int j=0;j<=n;j++) a[i][j]-=k*a[u][i]; } } ``` 知道了基本的变（魔）化（法），就来看怎么求一个合法解 #### 初始化 引用一下 Candy? 的博客 > 算法导论上有一个辅助线性规划的做法 > > 但我发现好多人都用了**随机初始化**的黑科技 > > 在所有 $a_{i,0}<0$ 的约束中随机选一个作为 $u$，再随机选一个 $a_{u,v}<0$ 作为 $v$，然后 `pivot(u,v)` 后 $a_{i,0}$ 就变正了... 但这可能出现循环，然后不断损精度，然后 WA/TLE，总之很玄学就对了 但是有「优先前面/后面」的选择倾向再加以小随机效果似乎不错 ```cpp int init() { while (1) { int u=0,v=0; for(int i=1;i<=m;i++) if(a[i][0]<-eps&&(!u||rand()%3==0)) u=i; if(!u)//说明这个解是可行解 return 1; for(int i=1;i<=n;i++) if(a[u][i]<-eps&&(!v||rand()%3==0)) v=i; if(!v) {//说明无解 printf("Infeasible\n"); return 0; } pivot(u,v); } } ``` #### 求解 每次找到一个 $v$ 使得 $A_{0,v}>0$，再找到一个 $u$ 使得 $A_{u,v}>0$，然后 `pivot(u,v)` 即可。 ##### 退化 在旋转的过程中，可能会存在保持答案不变的情况，这种现象称为退化。 ##### Bland 规则 在选择替入变量和替出变量的时候，我们总是选择满足条件的下标最小值。 - 替入变量 $v$ ：目标条件中，系数为正数（如果求最小值就是负数）的第一个作为替入变量 - 替出变量 $u$ ：对所有的约束条件中，选择对 $u$ 约束最紧的第一个（也就是 $\dfrac{A_{u,0}}{A_{u,v}}$ 最小的那个） ```cpp int simplex() { while (1) { int u=0,v=0; double mi=inf; for(int i=1;i<=n;i++) if(a[0][i]>eps) { v=i; return; } if(!v) return 1; for(int i=1;i<=m;i++) if(a[i][v]>eps&&a[i][0]/a[i][v]<mi) {//bland法则 mi=a[i][0]/a[i][v]; u=i; } if(!u) { //没有约束 printf("Unbounded\n"); return 0; } pivot(u,v); } } ``` #### 其它说明 ##### 为什么这能对 每次矩阵变化后的子问题是等价于初始问题的（易证） 并且这类问题有一个特殊性：它的可行域是凸的 意味着不会有除最优解外的局部最优解出现 ##### 这能跑多大 ~~心有多大，梦想有多大，它就能跑多大~~ 对于这种玄学算法，鬼才知道它能跑多大。 不过一般来说， $nm\le10^6$ 是~~不~~可行的
