Category: 算法研究


好吧好久没写,果然变懒了 = = 简略点好了

N x N 方阵,左上角入,左下角出,只能上下左右相邻移动,每格仅访问一次遍历所有格,求遍历方法数。

暴力 + 剪枝,没想出什么别的方法。最多就位运算优化一下:

   Test 1: TEST OK [0.000 secs, 2900 KB]
   Test 2: TEST OK [0.000 secs, 2900 KB]
   Test 3: TEST OK [0.000 secs, 2900 KB]
   Test 4: TEST OK [0.000 secs, 2900 KB]
   Test 5: TEST OK [0.000 secs, 2900 KB]
   Test 6: TEST OK [0.000 secs, 2900 KB]
   Test 7: TEST OK [0.162 secs, 2900 KB]

All tests OK.

N = 7 的时候在 0.162 秒貌似还不错,嗯。

代码嘛……算了不贴了。在 wordpress 上贴代码各种麻烦。研究一下再说

Advertisements

平面最近点对与最远点对

给定平面上 N 个点的坐标,求其中最近点对与最远点对。看似相似的问题,也都可以在 O(N logN) 时间内解出,却是两种完全不同的思路

http://dantvt.is-programmer.com/posts/7974.html
复制算了。MSN 居然能认出全部格式…不错

通过 USACO 4.2.1 Ditch 学习一下最大流算法 。可惜它给的测试数据几乎没有任何杀伤力,后面测试时我们采用 DD_engi 写的程序生成的加强版数据。

总体上来说,最大流算法分为两大类:增广路 (Augmenting Path) 和预流推进重标号 (Push Relabel)。也有算法同时借鉴了两者的长处,如 Improved SAP。本篇主要介绍增广路类算法,思想、复杂度及实际运行效率比较,并试图从中选择一种兼顾代码复杂度和运行效率的较好方案。以下我们将会看到,有时理论分析的时间复杂度并不能很好的反映一种算法的实际效率。

1. Ford – Fulkerson 方法

所有增广路算法的基础都是 Ford – Fulkerson 方法。称之为方法而不是算法是因为 Ford – Fulkerson 只提供了一类思想,在此之上的具体操作可有不同的实现方案。

给定一个有向网络 G(V,E) 以及源点 s 终点 t ,FF 方法描述如下:

Ford-Fulkerson 方法 (G,s,t)
1 将各边上流量 f 初始化为 0
2 while 存在一条增广路径 p
3     do 沿路径 p 增广流量 f
4 return f

假设有向网络 G 中边 (i,j) 的容量为 c(i,j) ,当前流量为 f(i,j) ,则此边的剩余流量即为 r(i,j) = c(i,j) – f(i,j) ,其反向边的剩余流量为 r(j,i) = f(i,j) 。有向网中所有剩余流量 r(i,j) > 0 的边构成残量网络 Gf ,增广路径 p 即是残量网络中从源点 s 到终点 t 的路径。

沿路径 p 增广流量f的操作基本都是相同的,各算法的区别就在于寻找增广路径 p 的方法不同。例如可以寻找从 s 到 t 的最短路径,或者流量最大的路径。

2. Edmonds – Karp 算法

Shortest Augmenting Path (SAP) 是每次寻找最短增广路的一类算法,Edmonds – Karp 算法以及后来著名的 Dinic 算法都属于此。SAP 类算法可统一描述如下:

Shortest Augmenting Path
1 x <– 0
2 while 在残量网络 Gx 中存在增广路 s ~> t
3     do 找一条最短的增广路径 P
4        delta <– min{rij:(i,j) 属于 P}
5        沿 P 增广 delta 大小的流量
6        更新残量网络 Gx
7 return x

在无权边的有向图中寻找最短路,最简单的方法就是广度优先搜索 (BFS),E-K 算法就直接来源于此。每次用一遍 BFS 寻找从源点 s 到终点 t 的最短路作为增广路径,然后增广流量 f 并修改残量网络,直到不存在新的增广路径。

E-K 算法的时间复杂度为 O(VE2),由于 BFS 要搜索全部小于最短距离的分支路径之后才能找到终点,因此可以想象频繁的 BFS 效率是比较低的。实践中此算法使用的机会较少。

3. Dinic 算法

BFS 寻找终点太慢,而 DFS 又不能保证找到最短路径。1970年 Dinic 提出一种思想,结合了 BFS 与 DFS 的优势,采用构造分层网络的方法可以较快找到最短增广路,此算法又称为阻塞流算法 (Blocking Flow Algorithm)。

首先定义分层网络 AN(f)。在残量网络中从源点 s 起始进行 BFS,这样每个顶点在 BFS 树中会得到一个距源点 s 的距离 d,如 d(s) = 0,直接从 s 出发可到达的点距离为 1,下一层距离为2 … 。称所有具有相同距离的顶点位于同一层,在分层网络中,只保留满足条件 d(i) + 1 = d(j) 的边,这样在分层网络中的任意路径就成为到达此顶点的最短路径。

Dinic 算法每次用一遍 BFS 构建分层网络 AN(f),然后在 AN(f) 中一遍 DFS 找到所有到终点 t 的路径增广;之后重新构造 AN(f),若终点 t 不在 AN(f) 中则算法结束。DFS 部分算法可描述如下:

 1 p <– s
 2 while s 的出度 > 0 do
 3     u <– p.top
 4     if u != t then
 5         if u 的出度 > 0 then
 6             设 (u,v) 为 AN(f) 中一条边
 7             p <– p, v
 8         else
 9             从 p 和 AN(f) 中删除点 u 以及和 u 连接的所有边
10     else
11         沿 p 增广
12         令 p.top 为从 s 沿 p 可到达的最后顶点
13 end while

 实际代码中不必真的用一个图来存储分层网络,只需保存每个顶点的距离标号并在 DFS 时判断 dist[i] + 1 = dist[j] 即可。Dinic 的时间复杂度为 O(V2E)。由于较少的代码量和不错的运行效率,Dinic 在实践中比较常用。具体代码可参考 DD_engi 网络流算法评测包中的标程,这几天 dinic 算法的实现一共看过和比较过将近 10 个版本了,DD 写的那个在效率上是数一数二的,逻辑上也比较清晰。

 4. Improved SAP 算法

 本次介绍的重头戏。通常的 SAP 类算法在寻找增广路时总要先进行 BFS,BFS 的最坏情况下复杂度为 O(E),这样使得普通 SAP 类算法最坏情况下时间复杂度达到了 O(VE2)。为了避免这种情况,Ahuja 和 Orlin 在1987年提出了Improved SAP 算法,它充分利用了距离标号的作用,每次发现顶点无出弧时不是像 Dinic 算法那样到最后进行 BFS,而是就地对顶点距离重标号,这样相当于在遍历的同时顺便构建了新的分层网络,每轮寻找之间不必再插入全图的 BFS 操作,极大提高了运行效率。国内一般把这个算法称为 SAP…显然这是不准确的,毕竟从字面意思上来看 E-K 和 Dinic 都属于 SAP,我还是习惯称为 ISAP 或改进的 SAP 算法。

 与 Dinic 算法不同,ISAP 中的距离标号是每个顶点到达终点 t 的距离。同样也不需显式构造分层网络,只要保存每个顶点的距离标号即可。程序开始时用一个反向 BFS 初始化所有顶点的距离标号,之后从源点开始,进行如下三种操作:(1)当前顶点 i 为终点时增广 (2) 当前顶点有满足 dist[i] = dist[j] + 1 的出弧时前进 (3) 当前顶点无满足条件的出弧时重标号并回退一步。整个循环当源点 s 的距离标号 dist[s] >= n 时结束。对 i 点的重标号操作可概括为 dist[i] = 1 + min{dist[j] : (i,j)属于残量网络Gf}。具体算法描述如下:

algorithm Improved-Shortest-Augmenting-Path
 1 f <– 0
 2 从终点 t 开始进行一遍反向 BFS 求得所有顶点的起始距离标号 d(i)
 3 i <– s
 4 while d(s) < n do
 5     if i = t then    // 找到增广路
 6         Augment
 7         i <– s      // 从源点 s 开始下次寻找
 8     if Gf 包含从 i 出发的一条允许弧 (i,j)
 9         then Advance(i)
10         else Retreat(i)    // 没有从 i 出发的允许弧则回退
11 return f

procedure Advance(i)
1(i,j) 为从 i 出发的一条允许弧
2 pi(j) <– i    // 保存一条反向路径,为回退时准备
3 i <– j        // 前进一步,使 j 成为当前结点

procedure Retreat(i)
1 d(i) <– 1 + min{d(j):(i,j)属于残量网络Gf}    // 称为重标号的操作
2 if i != s
3     then i <– pi(i)    // 回退一步

procedure Augment
1 pi 中记录为当前找到的增广路 P
2 delta <– min{rij:(i,j)属于P}
3 沿路径 P 增广 delta 的流量
4 更新残量网络 Gf

 算法中的允许弧是指在残量网络中满足 dist[i] = dist[j] + 1 的弧。Retreat 过程中若从 i 出发没有弧属于残量网络 Gf 则把顶点距离重标号为 n 。

 虽然 ISAP 算法时间复杂度与 Dinic 相同都是 O(V2E),但在实际表现中要好得多。要提的一点是关于 ISAP 的一个所谓 GAP 优化。由于从 s 到 t 的一条最短路径的顶点距离标号单调递减,且相邻顶点标号差严格等于1,因此可以预见如果在当前网络中距离标号为 k (0 <= k < n) 的顶点数为 0,那么可以知道一定不存在一条从 s 到 t 的增广路径,此时可直接跳出主循环。在我的实测中,这个优化是绝对不能少的,一方面可以提高速度,另外可增强 ISAP 算法时间上的稳定性,不然某些情况下 ISAP 会出奇的费时,而且大大慢于 Dinic 算法。相对的,初始的一遍 BFS 却是可有可无,因为 ISAP 可在循环中自动建立起分层网络。实测加不加 BFS 运行时间差只有 5% 左右,代码量可节省 15~20 行。

 5. 最大容量路径算法 (Maximum Capacity Path Algorithm)

 1972年还是那个 E-K 组合提出的另一种最大流算法。每次寻找增广路径时不找最短路径,而找容量最大的。可以预见,此方法与 SAP 类算法相比可更快逼近最大流,从而降低增广操作的次数。实际算法也很简单,只用把前面 E-K 算法的 BFS 部分替换为一个类 Dijkstra 算法即可。USACO 4.2 节的说明详细介绍了此算法,这里就不详述了。

 时间复杂度方面。BFS 是 O(E),简单 Dijkstra 是 O(V2),因此效果可想而知。但提到 Dijkstra 就不能不提那个 Heap 优化,虽然 USACO 的算法例子中没有用 Heap ,我自己还是实现了一个加 Heap 的版本,毕竟 STL 的优先队列太好用了不加白不加啊。效果也是非常明显的,但比起 Dinic 或 ISAP 仍然存在海量差距,这里就不再详细介绍了。

 6. Capacity Scaling Algorithm

 不知道怎么翻比较好,索性就这么放着吧。叫什么的都有,容量缩放算法、容量变尺度算法等,反正就那个意思。类似于二分查找的思想,寻找增广路时不必非要局限于寻找最大容量,而是找到一个可接受的较大值即可,一方面有效降低寻找增广路时的复杂度,另一方面增广操作次数也不会增加太多。时间复杂度 O(E2logU) 实际效率嘛大约稍好于最前面 BFS 的 E-K 算法,稀疏图时表现较优,但仍然不敌 Dinic 与 ISAP。

 7. 算法效率实测!

 重头戏之二,虽然引用比较多,哎~

 首先引用此篇强文 《Maximum Flow: Augmenting Path Algorithms Comparison》

 对以上算法在稀疏图、中等稠密图及稠密图上分别进行了对比测试。直接看结果吧:

稀疏图:

ISAP 轻松拿下第一的位置,图中最左边的 SAP 应该指的是 E-K 算法,这里没有比较 Dinic 算法是个小遗憾吧,他把 Dinic 与 SAP 归为一类了。最大流量路径的简单 Dijkstra 实现实在是太失败了 – -,好在 Heap 优化后还比较能接受……可以看到 Scaling 算法也有不错的表现。

 中等稠密图:

 

 ISAP 依然领先。最大流量算法依然不太好过……几个 Scaling 类算法仍然可接受。

 稠密图:

 

 ISAP……你无敌了!这次可以看出 BFS 的 Scaling 比 DFS 实现好得多,而且几乎与 Improved Scaling 不相上下,代码量却差不少。看来除 ISAP 外 BFS Scaling 也是个不错的选择。

 最后补个我自己实测的图,比较算法有很多是 DD 网络流算法评测包里的标程,评测系统用的 Cena,评测数据为 DD ditch 数据生成程序生成的加强版数据:

 

 我一共写了 7 个版本的 ISAP ,Dinic 算法也写了几个递归版的但效率都不高,只放上来一个算了。从上图来看似乎 ISAP 全面超越了号称最大流最快速算法的 HLPP,但在另外一台机器上测试结果与此却不大相同,有时 ISAP 与 HLPP 基本持平,有时又稍稍慢一些。在这种差距非常小的情况下似乎编译器的效果也比较明显。那个 HLPP 是用 PASCAL 写的,我不知在 Win32 平台下编译代码效率如何,至少我的几个 ISAP 用 VC2008 + SP1 编译比用 g++ 要强不少,也可能是参数设置的问题。

不过这些都是小事,关键是见证了 ISAP 的实际效率。从上面可以看出不加 GAP 优化的 ISAP 有几个测试点干脆无法通过,而不加 BFS 却无甚大影响,递归与非递归相差在 7% 左右的样子。综合以上表现,推荐采用 ISAP 不加 BFS,非递归 + GAP 优化的写法,Ditch 这道题一共也才 80 行左右代码。要提的一点是 GAP 优化用递归来表现的话不如 while 循环来得直接。另外,看起来 HLPP 表现确实很优秀,有机会也好好研究一下吧,预流推进重标号算法也是一大类呢……

最后附上一个 ISAP + GAP + BFS 的非递归版本代码,网络采用邻接表 + 反向弧指针:

  1. #include<cstdio>
  2. #include<memory>
  3. using namespace std;
  4.  
  5. const int maxnode = 1024;
  6. const int infinity = 2100000000;
  7.  
  8. struct edge{
  9.     int ver;    // vertex
  10.     int cap;    // capacity
  11.     int flow;   // current flow in this arc
  12.     edge *next; // next arc
  13.     edge *rev;  // reverse arc
  14.     edge(){}
  15.     edge(int Vertex, int Capacity, edge *Next)
  16.         :ver(Vertex), cap(Capacity), flow(0), next(Next), rev((edge*)NULL){}
  17.     void* operator new(size_t, void *p){
  18.         return p;
  19.     }
  20. }*Net[maxnode];
  21. int dist[maxnode]= {0}, numbs[maxnode] = {0}, src, des, n;
  22.  
  23. void rev_BFS(){
  24.     int Q[maxnode], head = 0, tail = 0;
  25.     for(int i=1; i<=n; ++i){
  26.         dist[i] = maxnode;
  27.         numbs[i] = 0;
  28.     }
  29.  
  30.     Q[tail++] = des;
  31.     dist[des] = 0;
  32.     numbs[0] = 1;
  33.     while(head != tail){
  34.         int v = Q[head++];
  35.         for(edge *e = Net[v]; e; e = e->next){
  36.             if(e->rev->cap == 0 || dist[e->ver] < maxnode)continue;
  37.             dist[e->ver] = dist[v] + 1;
  38.             ++numbs[dist[e->ver]];
  39.             Q[tail++] = e->ver;
  40.         }
  41.     }
  42. }
  43.  
  44. int maxflow(){
  45.     int u, totalflow = 0;
  46.     edge *CurEdge[maxnode], *revpath[maxnode];
  47.     for(int i=1; i<=n; ++i)CurEdge[i] = Net[i];
  48.     u = src;
  49.     while(dist[src] < n){
  50.         if(u == des){    // find an augmenting path
  51.             int augflow = infinity;
  52.             for(int i = src; i != des; i = CurEdge[i]->ver)
  53.                 augflow = min(augflow, CurEdge[i]->cap);
  54.             for(int i = src; i != des; i = CurEdge[i]->ver){
  55.                 CurEdge[i]->cap -= augflow;
  56.                 CurEdge[i]->rev->cap += augflow;
  57.                 CurEdge[i]->flow += augflow;
  58.                 CurEdge[i]->rev->flow -= augflow;
  59.             }
  60.             totalflow += augflow;
  61.             u = src;
  62.         }
  63.  
  64.         edge *e;
  65.         for(e = CurEdge[u]; e; e = e->next)
  66.             if(e->cap > 0 && dist[u] == dist[e->ver] + 1)break;
  67.         if(e){    // find an admissible arc, then Advance
  68.             CurEdge[u] = e;
  69.             revpath[e->ver] = e->rev;
  70.             u = e->ver;
  71.         } else {    // no admissible arc, then relabel this vertex
  72.             if(0 == (–numbs[dist[u]]))break;    // GAP cut, Important!
  73.             CurEdge[u] = Net[u];
  74.             int mindist = n;
  75.             for(edge *te = Net[u]; te; te = te->next)
  76.                 if(te->cap > 0)mindist = min(mindist, dist[te->ver]);
  77.             dist[u] = mindist + 1;
  78.             ++numbs[dist[u]];
  79.             if(u != src)
  80.                 u = revpath[u]->ver;    // Backtrack
  81.         }
  82.     }
  83.     return totalflow;
  84. }
  85.  
  86. int main(){
  87.     int m, u, v, w;
  88.     freopen("ditch.in", "r", stdin);
  89.     freopen("ditch.out", "w", stdout);
  90.     while(scanf("%d%d", &m, &n) != EOF){    // POJ 1273 need this while loop
  91.         edge *buffer = new edge[2*m];
  92.         edge *data = buffer;
  93.         src = 1; des = n;
  94.         while(m–){
  95.             scanf("%d%d%d", &u, &v, &w);
  96.             Net[u] = new((void*) data++) edge(v, w, Net[u]);
  97.             Net[v] = new((void*) data++) edge(u, 0, Net[v]);
  98.             Net[u]->rev = Net[v];
  99.             Net[v]->rev = Net[u];
  100.         }
  101.         rev_BFS();
  102.         printf("%d\n", maxflow());
  103.         delete [] buffer;
  104.     }
  105.     return 0;
  106. }

最小环 与 Floyd

Floyd 的 改进写法可以解决最小环问题,时间复杂度依然是 O(n^3),储存结构也是邻接矩阵
 
int mincircle = infinity;
Dist = Graph;

for(int k=0;k<nVertex;++k){
    //新增部分:
    for(int i=0;i<k;++i)
        for(int j=0;j<i;++j)
            mincircle = min(mincircle,Dist[i][j]+Graph[j][k]+Graph[k][i]);
    //通常的 floyd 部分:
    for(int i=0;i<nVertex;++i)
        for(int j=0;j<i;++j){
            int temp = Dist[i][k] + Disk[k][j];
            if(temp < Dist[i][j])
                Dist[i][j] = Dist[j][i] = temp;
        }
}

 
上面是对无向图的情况。
Floyd 算法保证了最外层循环到 k 时所有顶点间已求得以 0…k-1 为中间点的最短路径。一个环至少有3个顶点,设某环编号最大的顶点为 L ,在环中直接与之相连的两个顶点编号分别为 M 和 N (M,N < L),则最大编号为 L 的最小环长度即为 Graph(M,L) + Graph(N,L) + Dist(M,N) ,其中 Dist(M,N) 表示以 0…L-1 号顶点为中间点时的最短路径,刚好符合 Floyd 算法最外层循环到 k=L 时的情况,则此时对 M 和 N 循环所有编号小于 L 的顶点组合即可找到最大编号为 L 的最小环。再经过最外层 k 的循环,即可找到整个图的最小环。
 
若是有向图,只需稍作改动。注意考虑有向图中2顶点即可组成环的情况。
多进程对共享资源互斥访问及进程同步的经典问题
 
设有一文件F,多个并发读进程和写进程都要访问,要求:
(1)读写互斥
(2)写写互斥
(3)允许多个读进程同时访问
采用记录型信号量机制解决
 
较常见的写法:

semaphore fmutex=1, rdcntmutex=1;
//fmutex –> access to file; rdcntmutex –> access to readcount
int readcount = 0;
void reader(){
    while(1){
        wait(rdcntmutex);
        if(0 == readcount)wait(fmutex);
        readcount = readcount + 1;
        signal(rdcntmutex);
        //Do read operation …
        wait(rdcntmutex);
        readcount = readcount – 1;
        if(0 == readcount)signal(fmutex);
        signal(rdcntmutex);
    }
}
void writer(){
    while(1){
        wait(fmutex);
        //Do write operation …
        signal(fmutex);
    }
}

读进程只要看到有其他读进程正在访问文件,就可以继续作读访问;写进程必须等待所有读进程都不访问时才能写文件,即使写进程可能比一些读进程更早提出申请。所以以上解法实际是 读者优先 的解法。如果在读访问非常频繁的场合,有可能造成写进程一直无法访问文件的局面….
 
为了解决以上问题,需要提高写进程的优先级。这里另增加一个排队信号量:queue。读写进程访问文件前都要在此信号量上排队,通过区别对待读写进程便可达到提高写进程优先级的目的。另外再增加一个 writecount 以记录提出写访问申请和正在写的进程总数:

semaphore fmutex=1, rdcntmutex=1, wtcntmutex=1, queue=1;
//fmutex –> access to file; rdcntmutex –> access to readcount
//wtcntmutex –> access to writecount
int readcount = 0, writecount = 0;
void reader(){
    while(1){
        wait(queue);
        wait(rdcntmutex);
        if(0 == readcount)wait(fmutex);
        readcount = readcount + 1;
        signal(rdcntmutex);
        signal(queue);
        //Do read operation …
        wait(rdcntmutex);
        readcount = readcount – 1;
        if(0 == readcount)signal(fmutex);
        signal(rdcntmutex);
    }
}
void writer(){
    while(1){
        wait(wtcntmutex);
        if(0 == writecount)wait(queue);
        writecount = writecount + 1;
        signal(wtcntmutex);
        wait(fmutex);
        //Do write operation …
        signal(fmutex);
        wait(wtcntmutex);
        writecount = writecount – 1;
        if(0 == writecount)signal(queue);
        signal(wtcntmutex);
    }
}

每个读进程最开始都要申请一下 queue 信号量,之后在真正做读操作前即让出(使得写进程可以随时申请到 queue)。而只有第一个写进程需要申请 queue,之后就一直占着不放了,直到所有写进程都完成后才让出。等于只要有写进程提出申请就禁止读进程排队,变相提高了写进程的优先级。
 
通过类似思想即可实现读写进程的公平竞争:

semaphore fmutex=1, rdcntmutex=1, queue=1;
//fmutex –> access to file; rdcntmutex –> access to readcount
int readcount = 0;
void reader(){
    while(1){
        wait(queue);
        wait(rdcntmutex);
        if(0 == readcount)wait(fmutex);
        readcount = readcount + 1;
        signal(rdcntmutex);
        signal(queue);
        //Do read operation …
        wait(rdcntmutex);
        readcount = readcount – 1;
        if(0 == readcount)signal(fmutex);
        signal(rdcntmutex);
    }
}
void writer(){
    while(1){
        wait(queue);
        wait(fmutex);
        signal(queue);
        //Do write operation …
        signal(fmutex);
    }
}

读进程没变,写进程变成在每次写操作前都要等待 queue 信号量。
 
课本上一般只会写第一种解法吧。看了后两种方法即可发现,在第一个解法中,fmutex 信号量实际是双重身份,首先实现对文件的互斥访问,其次起到了和后面排队信号量 queue 相同的作用,只不过在那种排序下只能是读者优先。如果直接看过后两种解法,应该会有更清楚的理解吧。
有朋友回复说 在稠密图上,简单的dijkstra 并不比 heap + dijkstra 来得慢。仔细想了一下,理论上 简单dijkstra 时间复杂度为 O(V^2),heap + dijkstra 为 O((V+E)lgV),在完全图时前者仍为 O(V^2),后者变为 O( V^2 * lgV ),反而比 简单dijkstra 来得慢了。实际效果如何呢?还是来实验一下吧:
 
先把我的结果贴上来:
( Floyd 用来作对照;SPFA 顺带也加上了,这里是优化版的效率最高;HpDijkstra 表示二叉堆优化的 dijkstra;NvDijkstra 为简单dijkstra;NvDijkstra2 类似,只是用一个 list 来维护还未加入最终集合的顶点列表,避免循环所有顶点 )
一共五组数据,分别是 100、200、500、1000 和 1500 个顶点的完全图。上次实验电脑是老笔记本了,CPU PM 1.4G;这回换个台机 CPU P4 515 2.93G,测试数据可以变得稍大一些了。
——————————
Number of vertices: 100
 
Floyd consumed:       0.016 secs
SPFA(opt) consumed:   0.015 secs
HpDijkstra consumed:  0.031 secs
NvDijkstra consumed:  0.047 secs
NvDijkstra2 consumed: 0.094 secs
—————————–
Number of vertices: 200
 
Floyd consumed:       0.125 secs
SPFA(opt) consumed:   0.125 secs
HpDijkstra consumed:  0.188 secs
NvDijkstra consumed:  0.375 secs
NvDijkstra2 consumed: 0.796 secs
—————————–
Number of vertices: 500
 
Floyd consumed:       2.047 secs
SPFA(opt) consumed:   1.703 secs
HpDijkstra consumed:  2.469 secs
NvDijkstra consumed:  5.844 secs
NvDijkstra2 consumed: 4.406 secs
——————————
Number of vertices: 1000
 
Floyd consumed:       16.344 secs
SPFA(opt) consumed:   13.813 secs
HpDijkstra consumed:  18.250 secs
NvDijkstra consumed:  46.609 secs
NvDijkstra2 consumed: 35.734 secs
——————————
Number of vertices: 1500
 
Floyd consumed:       54.422 secs
SPFA(opt) consumed:   46.453 secs
HpDijkstra consumed:  59.156 secs
NvDijkstra consumed:  156.609 secs
NvDijkstra2 consumed: 116.296 secs
 
结果与预想的并不相同,heap + dijkstra 完胜 简单dijkstra。为什么,下面是我的分析,仅做参考吧
先把代码贴一下。程序与上次实验类似,只贴关键部分
先是公共的 relax函数,需要松弛时返回 true,否则返回 false

inline bool relax(int u, int v, vector<int> &d){
    int nlen=d[u]+vAdjMatrix[u][v];
    if(nlen<d[v]){
        d[v]=nlen;
        return true;
    }
    return false;
}

heap + dijkstra: (与上次实验的稍有不同,写法简化一些)

vector<bool> IsInset(N,false);
int setcount=0;
priority_queue< ver_weight,vector<ver_weight>,greater<ver_weight> > pq;
 
dist[iSource]=0;
pq.push(ver_weight(iSource,0));
 
while(setcount<N){  //O(V)
    ver_weight uw=pq.top();
    pq.pop();   // < O(lgV) worst: 2*lgV
    if(IsInset[uw.vertex])continue;
    IsInset[uw.vertex]=true;
    ++setcount;
    for(int v=0;v<N;++v)if(relax(uw.vertex,v,dist))  // Sum{V} -> O(E)
        pq.push(ver_weight(v,dist[v]));  // O(1) or < O(lgV)
}

简单dijkstra:

vector<bool> IsInset(N,false);
int setcount=0;
 
dist[iSource]=0;
 
while(setcount<N){  //O(V)
    int min=infinity,minvertex=0;
    for(int i=0;i<N;++i)if(!IsInset[i] && dist[i]<min){  //O(V)
        min=dist[i];
        minvertex=i;
    }
    IsInset[minvertex]=true;
    ++setcount;
    for(int v=0;v<N;++v)    //Sum{V} -> O(E)
        relax(minvertex,v,dist);
}

简单dijkstra + list 与上面类似就不贴了,只是多维护一个 list,时间复杂度还是一样的。
N 小的时候,加 list 的简单dijkstra 反而比较慢,可见维护 list 多了不少开销,不过 N 大了以后效果就体现出来了。
 
加 heap 效率比较高的原因,首先我想到的是 N 的大小还不够大。毕竟在分析时间复杂度的时候都是 N 趋近于无穷时的极限值。但在上面的数据中,随着 N 的增大,并没有看到加 heap 的效率有下降的趋势,它与不加 heap 两者差距反而越拉越大。因此这个理由大概是站不住脚的。
其次,我的 简单dijkstra 算法代码有问题以至于效率明显下降?暂时没看出来,一般 dijkstra 应该都是这么写的吧。
最后,还是老老实实从时间复杂度上分析。为此特地把《算法导论》又请出来,好好翻了下 dijkstra 算法时间复杂度分析的那一节,比较花时间的几步操作都在上面代码中加了注释。
先来看 简单dijkstra,这个没什么说的,while 循环对每个顶点都要进入一次,后面寻找最近顶点又是对 N 个顶点循环,平均也要 N/2 次,再后面的 relax 操作对所有顶点加起来其实就是对所有边的操作。所以总的时间复杂度是 V*V/2 + E = O(V^2 + E)。这个复杂度是实打实的,几乎对所有图都需要花费这么多时间,相当稳定。
然后再来看 heap + dijkstra。最外层 while 循环一样,仍是 N 次。下面是 heap 的 pop 操作,复杂度 lgN 也没有争议。最大问题就出在最后:如果 relax 有松弛操作,则将该点 push 入 heap。虽然 push 复杂度是 lgN 也没什么问题,但这个 push 并不是所有边都需要的操作,只有前面 relax 返回 true 的时候才需要。在计算时间复杂度的时候是按最坏情况,即所有边都需要 push 来算的。我认为这个放大得有些过了,好好想一下,要构造一个这样所有边都需要 push 的图的最坏情况都不是一件容易的事情,更别说普通的图了。那么在这 N*N 条边中,到底有多少真正做了这个操作呢?我也做了统计,程序中用一个变量计算了一下,结果如下:

N        push    push/E
50        171    6.84 %
300     1588    1.76 %
800     5010    0.78 %
1500  10335    0.46 %


可见需要 push 的边数远比 E=N*N 来得小。而且随着 N 的增大,push 操作数反而与 E 的比例变小了。如果假定随着边数 N 增大,需要 push 的次数只按 lgN/N 的程度增大,那么估算 heap + dijkstra 的复杂度即为:V*lgV + E*(lgV/V)*lgV + E = O(V*lgV*lgV + E),就比 简单dijkstra 来得小了。
补充: 标准的 dijkstra 写法中,最后 relax 返回 true 后应该是一个 decrease-key 操作,这里却用的 push,这样改动基于两点原因:(1)要完成 decrease-key 就需要先在 heap 中定位这个顶点,也即根据 value 找 key。而一个节点的位置在 heap 中是不固定的,为了效率还得先维护一个 value-key 的反向查找表;要么就在 heap 中保存指针,实际的 ver_weight 组保存在位置固定的 vector 中,不管哪个写起来都稍显复杂;其实还是有点懒,因为 STL 中的 heap 没有封装这个操作,呵呵。(2)由 push 代替 decrease-key 对效率影响不大。二者的时间复杂度都是 lgN,区别在于前者需要在 heap 中多维护一些无效数据,对时间上影响类似于 2lgN 与 lgN 的差别,空间上就是多费点内存啦,现在内存也不值钱 = = ….要是真觉得不爽就自己也写个 heap 吧…
 
总之,不管怎么说,即便是稠密图,加 heap 的效果还是很明显的,所以一般条件允许时还是尽量用 heap + dijkstra 代替 简单dijkstra 吧。而且使用 STL 的优先队列,加 heap 的代码反而比简单dijkstra 来得还少,连循环都省了,还是挺好写的 呵呵
 
总结下:
  稠密图 稀疏图 有负权边
单源问题 Dijkstra+heap SPFA (或Dijkstra+heap,根据稀疏程度) SPFA
APSP(无向图) SPFA(优化)/Floyd SPFA(优化) SPFA(优化)
APSP(有向图) Floyd SPFA (或Dijkstra+heap,根据稀疏程度) SPFA
虽然时间复杂度都清楚,不过实际运行起来如何心里还是没底,实践才是检验真理的标准啊
 
前面在USACO 3.2.6 Sweet Butter 那道题中已经看到,稀疏图中对单源问题来说 SPFA 的效率略高于 Heap+Dijkstra ;对于无向图上的 APSP (All Pairs Shortest Paths)问题,SPFA算法加上优化后效率更是远高于 Heap+Dijkstra。今次再来看一看稠密图上的实验结果。
 
(注:以下提到的Dijkstra,如无特别说明,均指加入二叉堆(Binary-Heap)优化的Dijkstra算法,程序中由STL的优先队列(piority_queue)实现)
 
程序随机生成一个 N 顶点的无向完全图,之后分别用三种算法求所有顶点对之间的最短路径并统计各算法耗费的时间。程序代码附后
 
(1)不加优化
试验一:
Number of vertexes: 100
Floyd consumed:     0.070 secs
SPFA consumed:      0.240 secs
Dijkstra consumed:  0.040 secs
 
试验二:
Number of vertexes: 300
Floyd consumed:     0.911 secs
SPFA consumed:      2.183 secs
Dijkstra consumed:  0.891 secs
 
试验三:
Number of vertexes: 600
Floyd consumed:     6.719 secs
SPFA consumed:      19.709 secs
Dijkstra consumed:  6.589 secs
 
可以看到完全图果然是Floyd算法的老本行,写起来复杂得多的Heap+Dijkstra比起区区三行的Floyd并不具有太多优势;虽然没有写进程序,不过可以预见不加Heap的Dijkstra是肯定要败给Floyd的。
比起前两个,SPFA就不那么好过了,基本耗费了2-3倍的时间才完成计算,可见在稠密图的单源最短路径问题上SPFA比起Dijkstra确实有很大劣势。
 
(2)SPFA针对无向图的APSP问题加优化,优化方法见前面文章所述
试验四:
Number of vertexes: 100
Floyd consumed:     0.070 secs
SPFA consumed:      0.070 secs
Dijkstra consumed:  0.080 secs
 
试验五:
Number of vertexes: 300
Floyd consumed:     0.981 secs
SPFA consumed:      0.641 secs
Dijkstra consumed:  0.902 secs
 
试验六:
Number of vertexes: 600
Floyd consumed:     6.820 secs
SPFA consumed:      5.077 secs
Dijkstra consumed:  6.590 secs
 
SPFA加上优化后效率有了大幅提高,不但超过了Floyd,比起Dijkstra还略高20%左右。我不打算继续针对完全图的情况做任何优化,因为这里要比较的应该是对一般图都适用的普通算法。
 
总结一下:
(1)对于稀疏图,当然是SPFA的天下,不论是单源问题还是APSP问题,SPFA的效率都是最高的,写起来也比Dijkstra简单。对于无向图的APSP问题还可以加入优化使效率提高2倍以上。
(2)对于稠密图,就得分情况讨论了。单源问题明显还是Dijkstra的势力范围,效率比SPFA要高2-3倍。APSP问题,如果对时间要求不是那么严苛的话简简单单的Floyd即可满足要求,又快又不容易写错;否则就得使用Dijkstra或其他更高级的算法了。如果是无向图,则可以把Dijkstra扔掉了,加上优化的SPFA绝对是必然的选择。
 
 

稠密图

稀疏图

有负权边

单源问题

Dijkstra+heap

SPFA (或Dijkstra+heap,根据稀疏程度)

SPFA

APSP(无向图)

SPFA(优化)/Floyd

SPFA(优化)

SPFA(优化)

APSP(有向图)

Floyd

SPFA (或Dijkstra+heap,根据稀疏程度)

SPFA
 
OK,今天总算彻底弄清这仨的恩怨情仇,至少一段时间内不用再为选择哪个而头大了……
 
最后附上实验程序的代码:
#include<iostream>
#include<fstream>
#include<vector>
#include<cstdlib>
#include<ctime>
#include<queue>
using namespace std;
 
const int infinity=1000000000;
int N;
vector< vector<int> > vAdjMatrix;
vector< vector<int> > dijkstra;
vector< vector<int> > spfa;
vector< vector<int> > floyd;
 
inline bool relax(int u, int v, vector<int> &d){
    int nlen=d[u]+vAdjMatrix[u][v];
    if(nlen<d[v]){
        d[v]=nlen;
        return true;
    }
    return false;
}
 
void DoFloyd(){
    for(int i=0;i<N;++i)for(int j=0;j<N;++j)
        floyd[i][j]=vAdjMatrix[i][j];
 
    int newlen;
    for(int k=0;k<N;++k)
        for(int u=0;u<N;++u)
            for(int v=0;v<N;++v)
                if((newlen=floyd[u][k]+floyd[k][v])<floyd[u][v])
                    floyd[u][v]=floyd[v][u]=newlen;
}
 
void DoSPFA(){
    for(int iSource=0;iSource<N;++iSource){
        queue<int> Q;
        vector<bool> IsInQ(N,false);
 
        //optimizing begin:
        if(iSource>0){
            int iShort=0;
            for(int k=1;k<iSource;++k)if(spfa[k][iSource]<spfa[iShort][iSource])
                iShort=k;
            for(int iDes=0;iDes<N;++iDes)
                spfa[iSource][iDes]=spfa[iShort][iSource]+spfa[iShort][iDes];
        }
        //optimizing end.
 
        Q.push(iSource);
        IsInQ[iSource]=true;
        spfa[iSource][iSource]=0;
        while(!Q.empty()){
            int u=Q.front();
            Q.pop();
            IsInQ[u]=false;
            for(int v=0;v<N;++v)
                if(relax(u,v,spfa[iSource]) && !IsInQ[v]){
                    Q.push(v);
                    IsInQ[v]=true;
                }
        }
    }
}
 
void DoDijkstra(){
    struct ver_weight{
        int vertex;
        int weight;
        ver_weight(int Vertex, int Weight):vertex(Vertex),weight(Weight){}
        bool operator>(const ver_weight &comp)const{
            return weight>comp.weight;
        }
    };
 
    for(int iSource=0;iSource<N;++iSource){
        vector<bool> IsInset(N,false);
        int setcount=0;
        priority_queue< ver_weight,vector<ver_weight>,greater<ver_weight> > pq;
 
        IsInset[iSource]=true;
        ++setcount;
        for(int v=0;v<N;++v)if(vAdjMatrix[iSource][v]<infinity){
            dijkstra[iSource][v]=vAdjMatrix[iSource][v];
            pq.push(ver_weight(v,dijkstra[iSource][v]));
        }
        while(setcount<N){
            ver_weight uw=pq.top();
            pq.pop();
            if(IsInset[uw.vertex])continue;
            IsInset[uw.vertex]=true;
            ++setcount;
            for(int v=0;v<N;++v)if(relax(uw.vertex,v,dijkstra[iSource]))
                pq.push(ver_weight(v,dijkstra[iSource][v]));
        }
    }
}
 
int main(){
    ofstream fout("out.txt");
    cout<<"Input the number of vertexes: ";
    cin>>N;
    vAdjMatrix.resize(N,vector<int>(N,infinity));
    floyd.resize(N,vector<int>(N,infinity));
    spfa.resize(N,vector<int>(N,infinity));
    dijkstra.resize(N,vector<int>(N,infinity));
    fout<<"Number of vertexes: "<<N<<‘\n’<<endl;
 
    //create the graph
    srand(unsigned(time(NULL)));
    for(int i=0;i<N;++i)for(int j=0;j<i;++j)
        vAdjMatrix[i][j]=vAdjMatrix[j][i]=rand();
    fout.precision(3);
    fout.setf(ios::fixed);
    clock_t start,finish;
 
    //floyd
    cout<<"Floyd start … \n";
    start=clock();
    DoFloyd();
    finish=clock();
    cout<<"Floyd end.\n\n";
    fout<<"Floyd consumed:     "<<(double)(finish-start)/CLOCKS_PER_SEC<<" secs\n";
 
    //SPFA
    cout<<"SPFA start … \n";
    start=clock();
    DoSPFA();
    finish=clock();
    cout<<"SPFA end.\n\n";
    fout<<"SPFA consumed:      "<<(double)(finish-start)/CLOCKS_PER_SEC<<" secs\n";
 
    for(int i=0;i<N;++i)for(int j=0;j<i;++j)
        if(floyd[i][j]!=spfa[i][j]){
            fout<<"SPFA worked wrong!\n";
            goto SPFA_end;   // try…catch are much better than goto…
        }
 
SPFA_end:
 
    //Dijkstra
    cout<<"Dijkstra start … \n";
    start=clock();
    DoDijkstra();
    finish=clock();
    cout<<"Dijkstra end.\n";
    fout<<"Dijkstra consumed:  "<<(double)(finish-start)/CLOCKS_PER_SEC<<" secs\n";
 
    for(int i=0;i<N;++i)for(int j=0;j<i;++j)
        if(floyd[i][j]!=dijkstra[i][j]){
            fout<<"Dijkstra worked wrong!\n";
            goto Dijk_end;
        }
 
Dijk_end:
 
    //system("pause");
    return 0;
}
 
——————-我–是–分–割–线——————
好吧,既然有朋友提出来写个 naive 的 dijkstra 看看,我就也来比较一下
这里。
 
—————-PS.————————————–
再补充,关于 SPFA 的两个优化 SLF 和 LLL:   懒得再写一篇了…
仍然拿随机生成的图做了实验,对单源问题,SLF 可使速度提高 15 ~ 20 %; SLF + LLL 可提高约 50 % (都是与纯SPFA相比)。
不过 LLL 写起来加的东西比较多一点,SLF 简单能加就随便加吧反正多一行而已
即便如此在稠密图上与 Heap + Dijkstra 相比仍有差距
 
对无向图 APSP 如果加上前述优化的话再加 SLF 或 LLL 反正略有提高,但不明显,最多好像也就 7~8 % 的样子,忘了。加不加无所谓。

单源最短路径算法小结

这里就不写具体算法了,只将他们的时间复杂度、适用范围、代码复杂程度简单做个比较
待搜索的图都指有向图(无向图类似)。储存方式均为邻接表
 
一、广度优先搜索(BFS)
时间复杂度:O(V+E),效率很高
适用范围:(很窄)仅适于无权边的图。即每条边长度都为1的情况
代码复杂程度:一般,需队列
 
二、Bellman-Ford
时间复杂度:O(VE),效率一般
适用范围:(很广)允许存在负权边,能够判断图中是否存在从源点可到达的负权环路
代码复杂程度:较易
 
三、有向无环图算法
时间复杂度:O(V+E),效率很高
适用范围:(很窄)仅适于有向无环图
代码复杂程度:一般,需拓扑排序
 
四、Dijkstra
适用范围:(一般)不允许存在负权边
这个算法复杂度取决于"取最小"(Extract-min)操作使用的算法
       Extract-min操作                时间复杂度               代码复杂程度
顺序检测所有点决定最小值           O(V^2)                        一般
使用Binary-Heap(优先队列)       O((V+E)lgV)                 较复杂
    使用Fibonacci-Heap               O(VlgV+E)                  较复杂
 
五、SPFA (Shortest Path Faster Algorithm)
时间复杂度:O(kE),k为一较小常量。效率很高
适用范围:(较广),允许存在负权边,但不允许负权环路
代码复杂程度:较易,需队列
这个算法可算是 Bellman-Ford 的优化版本,去除冗余的Relax操作,效率有很大提升
有些奇怪的是这个算法在《算法导论》和 Wikipedia 上竟然都没有介绍,而只在国内的OI相关网站论坛才发现的解释,大多还不是很全面。从适用范围和代码复杂程度来看绝对是一个值得推荐的算法(比熟知的 Dijkstra 都要好),效率上也不比使用 Heap 优化的 Dijkstra 低,一般来说甚至还要高(还有待更多题目实践验证)。即使寻找所有顶点对之间的最短路径这个算法也是值得考虑的(对每个顶点运行一次SPFA)
 
之后再遇到最短路径的题目我会尽量把Dijkstra和SPFA都用一下做个比较
 
附:这一篇做了较详细的比较

USACO section 2.1

先写几个常用缩写:
AC = Accept
WA = Wrong Answer
TLE = Time Limit Exceeded
MLE = Memory Limit Exceeded
 
嗯,虽说这几题都比较简单,不过1个多小时AC掉3道题感觉还是不错滴 ^_^
 
USER: Jinlin Wang [wjl1241]
TASK: frac1
LANG: C++
Compiling…
Compile: OK
Executing…
   Test 1: TEST OK [0.022 secs, 2840 KB]
   Test 2: TEST OK [0.000 secs, 2844 KB]
   Test 3: TEST OK [0.000 secs, 2840 KB]
   Test 4: TEST OK [0.000 secs, 2844 KB]
   Test 5: TEST OK [0.011 secs, 2844 KB]
   Test 6: TEST OK [0.000 secs, 2840 KB]
   Test 7: TEST OK [0.000 secs, 2840 KB]
   Test 8: TEST OK [0.011 secs, 2844 KB]
   Test 9: TEST OK [0.011 secs, 2840 KB]
   Test 10: TEST OK [0.011 secs, 2840 KB]
   Test 11: TEST OK [0.065 secs, 2844 KB]
All tests OK.
YOUR PROGRAM (‘frac1’) WORKED FIRST TIME!  That’s fantastic
— and a rare thing.  Please accept these special automated
congratulations.
这题比较顺利一次PASS。没感觉需要什么特殊算法,把数据结构组织好就差不多了
 
USER: Jinlin Wang [wjl1241]
TASK: sort3
LANG: C++
Compiling…
Compile: OK
Executing…
   Test 1: TEST OK [0.011 secs, 2848 KB]
   Test 2: TEST OK [0.000 secs, 2844 KB]
   Test 3: TEST OK [0.000 secs, 2848 KB]
   Test 4: TEST OK [0.000 secs, 2848 KB]
   Test 5: TEST OK [0.011 secs, 2844 KB]
   Test 6: TEST OK [0.000 secs, 2844 KB]
   Test 7: TEST OK [0.000 secs, 2848 KB]
   Test 8: TEST OK [0.011 secs, 2844 KB]
All tests OK.
Your program (‘sort3’) produced all correct answers!  This is your
submission #3 for this problem.  Congratulations!
一开始想错了。。。结果WA了两次
 
USER: Jinlin Wang [wjl1241]
TASK: holstein
LANG: C++
Compiling…
Compile: OK
Executing…
   Test 1: TEST OK [0.000 secs, 2848 KB]
   Test 2: TEST OK [0.000 secs, 2844 KB]
   Test 3: TEST OK [0.000 secs, 2848 KB]
   Test 4: TEST OK [0.000 secs, 2844 KB]
   Test 5: TEST OK [0.022 secs, 2844 KB]
   Test 6: TEST OK [0.000 secs, 2844 KB]
   Test 7: TEST OK [0.000 secs, 2848 KB]
   Test 8: TEST OK [0.011 secs, 2844 KB]
   Test 9: TEST OK [0.000 secs, 2844 KB]
   Test 10: TEST OK [0.011 secs, 2844 KB]
All tests OK.
Your program (‘holstein’) produced all correct answers!  This is your
submission #2 for this problem.  Congratulations!
前一段写五子棋程序写过一个组合类,直接拿来就用了。可方便了
本来应该一次就对的,第一遍最后结果少输出一个回车就给判错了。。。
 
================ 华丽 滴 分割 =========================
 
可惜这次还是稍稍有点无聊的内容,辜负某些同学的期望了,不好意思啊。不过本人感情方面确实没啥可汇报的,见到女生都不知该说啥,更别说平时基本都见不到几个女生。就由它去吧……
嗯,下次尝试写点不同的话题

暂时放弃了。。。

一个高效的VCF算法看来比我开始预计的要复杂得多……
 
最容易,也是最开始想到的是广度优先搜索。实现起来相当方便,一棵树和一个队列即可解决。不过缺点也是显而易见的,当解存在于7-8层深度以上时,程序只有把7层以下所有冲四的全排列都列举完才能发现解。而全排列这个东东随着搜索深度增加增长太快了……中盘时局面上同时出现10个左右可以冲四的点是很常见的情况,这时就要搜索10!次。。。马上就不行了。我用Renlib软件自带的一个VCF问题谱来试验这个算法,31个VCF问题中解决了29个,大部分都是1秒内即给出解法,少数几个用了2-3秒。未解决的两个是搜索20秒后仍然没结果。看了一下这两个问题,解法都在12步以上,并且每一步都有很多点可以选择。。。怪不得
 
这样看来,貌似深度优先可行。最开始没用这个是因为深度优先不能保证找到最短的VCF路径,若局面上同时存在两条VCF解,一个3步另一个7步,则深度优先有可能找到的是9步的解——其中除了必要的7步外还有2步废棋。。。这样最后还得花时间来清除废棋,写起来比较麻烦。不过从结果来看,当存在VCF时,深度优先仍然比广度优先快很多,特别是解的路径比较深的时候。之后就是缺点了,同样明显。若一个比较复杂的局面却无解,则深度优先与广度优先一样都要遍历整个搜索树才能判断,也即同样存在上面的全排列问题。。。仍然用那31个VCF问题试验,这次全部解决掉了,速度也很快。不过败在了我自己随便布的几个局面下,局面上冲四点很多,但又没有解,它就不行了
 
这样看来,要改变思考方法了。想一下我们人在判断VCF时首先是找组合,比如这3个冲四能形成活三,则可能存在VCF;之后才考虑走棋顺序,比如会不会让对方形成反冲四啊、会不会形成禁手啊什么的。等于将全排列变为了先组合再选择性的排列,一下把要考虑的情况缩减到了几十或几百几千分之一。不过这时又有了新的问题,即有些冲四是需要其他几个冲四做铺垫才能出现的,这样就还得把有可能形成新冲四点的冲四组合都走掉。。。比前两种算法要复杂得多,然后在试验这些组合时还要考虑对方反冲四啦、通过冲四解禁手啦等等情况。最后克服各种困难终于写出一个试验版本。。。而且只判断黑棋的VCF(白棋类似,不过还要考虑抓禁手。。。更复杂),程序长度一下增长了3倍多,不过速度确实有明显提高。不管是有解还是无解情况都能很快判断。但同样用那31个VCF问题试验时却发现有几个问题没找到解。。。后来经过仔细分析发现有一种较少见的情况有可能出现解的在这个算法中被排除掉了,不仔细写了,总之要解决这个BUG得改动很多地方。
 
至此为了这个程序我已经花了太多时间在上面了,想想看应该还有更多的事情要做,有更多的东西要学。。。不能再在这上面耽误太久了。考虑放弃……就先随它去吧,以后说不定哪一天又想重新拿起再继续吧。希望那时我会有更多新想法可用。。。