渗漏问题算法

给定由随机分布的绝缘材料和金属材料构成的组合系统：金属材料占多大比例才能使组合系统成为电导体？给定一个表面有水的多孔景观（或下面有油），水将在什么条件下能够通过底部排出（或油渗透到表面）？科学家们已经定义了一个称为渗透（percolation）的抽象过程来模拟这种情况。

模型

我们使用N×N网格点来模型一个渗透系统。每个格点或是open格点或是blocked格点。一个full site是一个open格点，它可以通过一连串的邻近（左，右，上，下）open格点连通到顶行的一个open格点。如果在底行中有一个full site格点，则称系统是渗透的。（对于绝缘/金属材料的例子，open格点对应于金属材料，渗透系统有一条从顶行到底行的金属路径，且full sites格点导电。对于多孔物质示例，open格点对应于空格，水可能流过，从而渗透系统使水充满open格点，自顶向下流动。）

问题

在一个著名的科学问题中，研究人员对以下问题感兴趣：如果将格点以概率p独立地设置为open格点（因此以概率1-p被设置为blocked格点），系统渗透的概率是多少？当p=0时，系统不会渗出; 当p=1时，系统渗透。下图显示了20×20随机网格（左）和100×100随机网格（右）的格点空置概率p与渗滤概率。

当N足够大时，存在阈值p*，使得当p <p*，随机N´N网格几乎不会渗透，并且当p>p*时，随机N´N网格几乎总是渗透。尚未得出用于确定渗滤阈值p*的数学解。我们的任务是编写一个计算机程序来估计p*。

我们可以使用深度优先搜索、广度优先搜索和合并-查找（union-find）数据结构，编写程序通过蒙特卡罗模拟（Monte Carlo simulation）来估计渗透阈值的值。在问题的核心部分，我将以三种方法解决这道问题，它们分别是深度优先搜索、Quick-Find和Quick-Union。

编程环境

C#, .Net Framework 4.6.2, Microsoft Visual Studio 2015 in Windows 10

深度优先搜索

本题可以将问题转化为用于解决迷宫的出口路径的深度优先搜索问题，即把block住的部分看作迷宫的墙，把open的地方看成是通路，把灌水的地方看做是走过的地方，然后通过深度优先搜索寻找一条从最上方某点到最下方某点的联通的路径。如果存在这样一条路径，那么，这个地图是渗漏的，否则是不渗漏的。

Quick-Find和Quick-Union

这两种方法的思路相同，都是使用了并查集的思想来解决这个问题。从不渗漏的地图中，随机打通一点，然后进行Union-Find，如果，地图的最上面一层的open点，与最下面一层的open点是同一个连通集的话，那么，此时这个地图是渗漏的，否则是不渗漏的。

不论用什么方法，当地图不渗漏的时候，随机打通一点让其为open状态，然后再判断其是否渗漏，直到渗漏为止，然后记录相关信息，计算open/all的比例，并计算多次运行的平均值，标准差和置信区间。

事实上，在这道题中，在地图尺寸不是非常大的情况下，(大概在100*100之内)，使用深度优先搜索获得的运行速度比Union-Find快得多，并且此时Quick-Union比Quick-Find也要快一些。

关键代码解释

public class Class_Percolation

使用了深度优先搜索算法的类

private bool dfs(int[][] v_map, point p_start, point p_end)

使用深度优先搜索对指定的有一定open点的地图进行路径搜索，搜索从起始点到终止点是否存在一个路径，使得地图渗漏。如果最后路径栈空了，说明不存在这样的路径，那么就是不连通的。

public void openblocks_ramdom(int num)

随机挖空一定数量的砖块，以改变地图

public class Class_Percolation_union_find

使用了Quick-Find和Quick-Union算法的类

HashSet<point> pset_unblock = new HashSet<point>();

用来存放已经打通的点的集合类，HashSet

public void open(int x,int y, int algorithms)

打开一个点，并使用指定算法将该点与其周围的点联通

public void union_find(point p,int q_x, int q_y, int algorithms)

用来将指定两集合的点做Quick-Find或Quick-Union的函数

public bool isPercolation(int algorithms)

用来通过指定算法判断该地图是否渗漏

public class QuickFind

封装了Quick-Find算法的类

public class QuickUnion

封装了Quick-Union算法的类

private point find(point p)

做Find操作的函数

public void union(point p, point q)

做Union操作的函数

运行结果

从结果中我们可以看出，大约在渗透的open比例为0.56-0.58左右的时候，开始出现上下全连通的。

我最开始的时候，首先使用的是深度优先搜索算法实现的，发现速度还是非常快的，然后使用了Union-Find的两种算法来做，发现这两种算法的速度竟然远比DFS慢。后来经过多次测试，发现，在地图尺寸比较小的时候，DFS的效率比较高，省去了很多没用的Union和Find操作。不过，当地图尺寸较大时，DFS带来的搜索空间的快速增长，使得该算法的效率急剧下降，大约在100*100尺寸的时候，三种算法的效率已经差不多了，而且Quick-Union算法的时间开销也在上升，此时Quick-Find似乎是一个在大尺寸地图下的较优算法。

那么我们是不是可以这样，在地图尺寸较小的时候，使用DFS算法，将问题转换成寻找一条连通的路径，在地图尺寸中等大小的时候，改用Quick-Union算法，通过判断顶部和底部是否在同一个连通集来求解，当地图尺寸较大时，使用Quick-Find算法，以提高调用次数更多的Find方法的效率。

实现代码

Program.cs

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace 渗漏问题
{
    class Program
    {
        /// <summary>
        /// 渗漏结果
        /// </summary>
        public struct myresult
        {
            public int num_open;
            public int num;
            public bool isPercolation;
        }

        /// <summary>
        /// 置信区间
        /// </summary>
        public struct myconfidence
        {
            public double low;
            public double high;
        }

        static void Main(string[] args)
        {
            Program p = new Program();
            //p.test();
            //p.test(16,5);
            Console.WriteLine("DFS >>>");
            p.control_test(200, 0);
            Console.WriteLine("UnionFind >>>");
            p.control_test(200, 1);
            Console.WriteLine("UnionQuick >>>");
            p.control_test(200, 2);
            //myresult r= p.test_UnionFind();
            //Console.WriteLine("{0} {1} {2} {3}", r.isPercolation, r.num_open, r.num, (double)r.num_open / (double)r.num);
            System.Console.ReadKey();
        }

        public Program()
        {
            
        }

        private void control_test(int num, int algorithms)
        {
            List<myresult> r = new List<myresult>();

            System.Console.WriteLine("开始测试...");
            //int num = 200;

            DateTime beforDT = System.DateTime.Now;//开始时间

            for (int i = 0; i < num; i++)
            {
                //System.Console.WriteLine("正在测试第 " + i + " 项...");
                myresult tmp;
                if (algorithms == 0)
                    tmp = test();
                else
                    tmp = test_UnionFind(algorithms - 1);
                System.Threading.Thread.Sleep(2);//加一点延时会使得结果比较稳定
                r.Add(tmp);
            }

            DateTime afterDT = System.DateTime.Now;
            TimeSpan ts = afterDT.Subtract(beforDT);

            System.Console.WriteLine("测试完成，开始统计...");
            myresult[] re = r.ToArray();
            List<double> list_rate = new List<double>();
            for (int i = 0; i < num; i++)
            {
                double tmp = (double)re[i].num_open / (double)re[i].num;
                list_rate.Add(tmp);
            }
            
            //均值
            double rate_means = statistic_means(list_rate);
            System.Console.WriteLine("平均渗透open比例：" + rate_means);

            //标准差
            double rate_stvdev = statistic_stddev(rate_means, list_rate);
            System.Console.WriteLine("标准差：" + rate_stvdev);

            //置信区间
            myconfidence con = statistic_confidence(rate_means, rate_stvdev, num);
            System.Console.WriteLine("95%置信区间：  [ " + con.low + " , " + con.high + " ]");

            Console.WriteLine("DateTime总共花费{0}ms.", ts.TotalMilliseconds - 2 * num);

            System.Console.WriteLine("===============\n");

            //System.Console.ReadKey();
        }
        
        /// <summary>
        /// 随机尺寸测试
        /// </summary>
        /// <returns></returns>
        public myresult test()
        {
            Class_Percolation p = new Class_Percolation();

            int[][] map;
            //System.Console.WriteLine("正在生成地图...");
            //Random ram = new Random(DateTime.Now.Minute * 10000 + DateTime.Now.Second * 100 + DateTime.Now.Millisecond);
            Random ram = new Random();
            int size_x = 5 + ram.Next() % 45;
            int size_y = 5 + ram.Next() % 45;
            //System.Console.WriteLine(size_x + " " + size_y);
            int num_open = 1;
            map = make_map(size_y, size_x, num_open);

            p.map = map;
            //p.print_map();
            //System.Console.WriteLine("挖空一个...");
            //p.openblocks_ramdom(1);
            //p.print_map();

            //System.Console.WriteLine("正在测试是否渗漏...");
            bool result_test = p.test();
            while (result_test != true)
            {
                p.openblocks_ramdom(1);
                num_open++;
                result_test = p.test();
            }
            //p.print_map();
            // System.Console.WriteLine("实验结果：" + result);
            //System.Console.Read();
            myresult r = new myresult();
            r.num_open = num_open;
            r.num = size_x * size_y;
            r.isPercolation = true;
            return r;
        }

        /// <summary>
        /// 自定义尺寸的测试
        /// </summary>
        /// <param name="x"></param>
        /// <param name="y"></param>
        /// <returns></returns>
        public myresult test(int x,int y)
        {
            Class_Percolation p = new Class_Percolation();

            int[][] map;
            //System.Console.WriteLine("正在生成地图...");
            
            int size_x = x;
            int size_y = y;
            //System.Console.WriteLine(size_x + " " + size_y);
            int num_open = 1;
            map = make_map(size_y, size_x, num_open);

            p.map = map;
            //p.print_map();
            //System.Console.WriteLine("挖空一个...");
            //p.openblocks_ramdom(1);
            p.print_map();

            //System.Console.WriteLine("正在测试是否渗漏...");
            bool result_test = p.test();
            while (result_test != true)
            {
                //System.Console.WriteLine("openblocks_ramdom()前...");
                p.openblocks_ramdom(1);
                num_open++;
                //System.Console.WriteLine("test()前...");
                result_test = p.test();
            }
            p.print_map();
            // System.Console.WriteLine("实验结果：" + result);
            //System.Console.Read();
            myresult r = new myresult();
            r.num_open = num_open;
            r.num = size_x * size_y;
            r.isPercolation = true;
            return r;
        }
        
        /// <summary>
        /// 随机尺寸测试
        /// </summary>
        /// <returns></returns>
        public myresult test_UnionFind(int algorithms)
        {
            Class_Percolation_union_find p = new Class_Percolation_union_find();

            point[][] map;
            //System.Console.WriteLine("正在生成地图...");
            //Random ram = new Random(DateTime.Now.Minute * 10000 + DateTime.Now.Second * 100 + DateTime.Now.Millisecond);
            Random ram = new Random();
            int size_x = 5 + ram.Next() % 45;
            int size_y = 5 + ram.Next() % 45;
            //System.Console.WriteLine(size_x + " " + size_y);
            //int num_open = 1;
            //map = make_map(size_y, size_x, num_open);
            map = new point[size_x][];
            for(int i=0;i<map.Length;i++)
            {
                map[i] = new point[size_y];
            }

            p.map = map;
            //p.print_map();
            //System.Console.WriteLine("挖空一个...");
            //p.openblocks_ramdom(1);
            //p.print_map();

            //System.Console.WriteLine("正在测试是否渗漏...");

            if(algorithms== p.ALGORITHMS_QUICK_FIND)
            {
                while (!p.isPercolation(p.ALGORITHMS_QUICK_FIND))
                {
                    p.open(p.ALGORITHMS_QUICK_FIND);
                }
            }

            if (algorithms == p.ALGORITHMS_QUICK_UNION)
            {
                while (!p.isPercolation(p.ALGORITHMS_QUICK_UNION))
                {
                    p.open(p.ALGORITHMS_QUICK_UNION);
                }
            }

                



            //p.print_map();
            // System.Console.WriteLine("实验结果：" + result);
            //System.Console.Read();
            myresult r = new myresult();
            r.num_open = p.GetOpenedCount();
            r.num = size_x * size_y;
            r.isPercolation = true;
            return r;
        }
        /// <summary>
        /// 生成随机地图
        /// </summary>
        /// <param name="x"></param>
        /// <param name="y"></param>
        /// <param name="count_open"></param>
        /// <returns></returns>
        public int[][] make_map(int x,int y,int count_open)
        {
            int size_x = x, size_y = y;
            int[][] map = new int[size_y][];
            //Random r = new Random(DateTime.Now.Hour * 10000 + DateTime.Now.Minute * 100 + DateTime.Now.Second);
            Random r = new Random();
            for (int i = 0; i < map.Length; i++)
            {
                map[i] = new int[size_x];
                for (int j = 0; j < map[i].Length; j++)
                {
                    map[i][j] = -1;
                }
            }
            int x_open = 0, y_open = 0;
            for (int i=0;i<count_open;i++)
            {
                x_open = r.Next() % size_x;
                y_open = r.Next() % size_y;
                map[y_open][x_open] = 0;
            }
            return map;
        }

        /// <summary>
        /// 统计均值
        /// </summary>
        /// <param name="list_rate"></param>
        /// <returns></returns>
        public double statistic_means(List<double> list_rate)
        {
            double rate_means = 0.0;
            double[] rate = list_rate.ToArray();
            for (int i = 0; i < rate.Length; i++)
            {
                rate_means += rate[i];
            }
            rate_means = rate_means / rate.Length;

            return rate_means;
        }

        /// <summary>
        /// 标准差
        /// </summary>
        /// <param name="means"></param>
        /// <param name="list_rate"></param>
        /// <returns></returns>
        public double statistic_stddev(double means,List<double> list_rate)
        {
            double rate_stddev = 0.0;//标准差
            double[] rate = list_rate.ToArray();

            for(int i=0;i<rate.Length;i++)
            {
                //rate_stddev += Math.Pow(rate[i] - means, 2.0);
                rate_stddev += (rate[i] - means) * (rate[i] - means);
            }
            rate_stddev /= (rate.Length - 1);
            return rate_stddev;
        }

        /// <summary>
        /// 计算置信区间
        /// </summary>
        /// <param name="means">均值</param>
        /// <param name="stddev">方差</param>
        /// <returns></returns>
        public myconfidence statistic_confidence(double means,double stddev,int T)
        {
            double tmp = 0.0;
            tmp = 1.96 * stddev / Math.Pow(T, 0.5);

            myconfidence con = new myconfidence();
            con.low = means - tmp;
            con.high = means + tmp;
            return con;
        }
    }
}

Class_Percolation.cs

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace 渗漏问题
{
    public class Class_Percolation
    {
        /// <summary>
        /// 坐标点
        /// </summary>
        struct point
        {
            public int x;
            public int y;
        }
        /// <summary>
        /// 方格为空
        /// </summary>
        public const int STATUS_OPEN = 0;
        /// <summary>
        /// 方格阻塞
        /// </summary>
        public const int STATUS_BLOCKS = -1;
        /// <summary>
        /// 方格灌水
        /// </summary>
        public const int STATUS_FULL = 1;

        private int[][] _map;
        /// <summary>
        /// 地图
        /// </summary>
        public int[][] map
        {
            get { return _map; }
            set { _map = value; }
        }

        public Class_Percolation()
        {
            _map = null;
        }
        public Class_Percolation(int[][] v_map)
        {
            _map = v_map;
        }
        /// <summary>
        /// 测试是否渗漏
        /// </summary>
        /// <returns>渗漏则返回true</returns>
        public bool test()
        {
            Exception e = new Exception();
            if (_map == null) throw e;
            int height = _map.Length;
            int width = _map[0].Length;
            int[][] tmp_map = add_head_bottom(_map);

            point p_start, p_end;
            p_start.x = 0;
            p_start.y = 0;
            p_end.y = tmp_map.Length - 1;
            p_end.x = tmp_map[0].Length - 1;

            return dfs(tmp_map, p_start, p_end);
        }

        /// <summary>
        /// 给地图添加全连通的顶部和底部，以方便计算
        /// </summary>
        /// <param name="v_map">原始地图</param>
        /// <returns>添加顶部底部后的地图</returns>
        private int[][] add_head_bottom(int[][] v_map)
        {
            int height = v_map.Length;
            int width = v_map[0].Length;
            int[][] tmp_map = new int[height + 2][];

            tmp_map[0] = new int[width];
            tmp_map[tmp_map.Length - 1] = new int[width];

            //顶部底部全部连通
            for (int i = 0; i < tmp_map[0].Length; i++)
            {
                tmp_map[0][i] = STATUS_OPEN;
                tmp_map[tmp_map.Length - 1][i] = STATUS_OPEN;
            }
            //中间依次复制
            for (int i = 0; i < v_map.Length; i++)
            {
                tmp_map[i + 1] = new int[width];
                v_map[i].CopyTo(tmp_map[i + 1], 0);
            }
            return tmp_map;
        }

        /// <summary>
        /// 深度优先搜索，判断是否联通
        /// </summary>
        /// <param name="v_map"></param>
        /// <param name="r_start"></param>
        /// <param name="p_end"></param>
        /// <returns>返回是否联通</returns>
        private bool dfs(int[][] v_map, point p_start, point p_end)
        {
            Stack<point> point_path = new Stack<point>();

            int height = v_map.Length;
            int width = v_map[0].Length;
            //初始化
            point_path.Push(p_start);
            point p_current = p_start;
            v_map[p_current.y][p_current.x] = STATUS_FULL;//灌水
            //开始行走
            while ((p_current.x != p_end.x || p_current.y != p_end.y) && point_path.Count > 0 )
            {
                //System.Console.WriteLine(p_current.x+"\t"+p_current.y);
                //向下走
                if (p_current.y + 1 < height && v_map[p_current.y + 1][p_current.x] == STATUS_OPEN)
                {
                    p_current.y = p_current.y + 1;
                    v_map[p_current.y][p_current.x] = STATUS_FULL;//灌水
                    point_path.Push(p_current);
                    continue;
                }
                //向左走
                if (p_current.x > 0 && v_map[p_current.y][p_current.x - 1] == STATUS_OPEN)
                {
                    p_current.x = p_current.x - 1;
                    v_map[p_current.y][p_current.x] = STATUS_FULL;//灌水
                    point_path.Push(p_current);
                    continue;
                }
                //向右走
                if (p_current.x + 1 < width && v_map[p_current.y][p_current.x + 1] == STATUS_OPEN)
                {
                    p_current.x = p_current.x + 1;
                    v_map[p_current.y][p_current.x] = STATUS_FULL;//灌水
                    point_path.Push(p_current);
                    continue;
                }
                //向上走
                if (p_current.y > 0 && v_map[p_current.y - 1][p_current.x] == STATUS_OPEN)
                {
                    p_current.y = p_current.y - 1;
                    v_map[p_current.y][p_current.x] = STATUS_FULL;//灌水
                    point_path.Push(p_current);
                    continue;
                }
                //否则，所有的路都走不通了，回滚吧
                p_current = point_path.Pop();
            }

            //如果栈不为空，则存在一条渗漏路径
            if(point_path.Count > 0)
            {
                return true;
            }

            return false;
        }
        
        /// <summary>
        /// 随机挖空一定数量的砖块，以改变地图
        /// </summary>
        /// <param name="num">挖空数量</param>
        public void openblocks_ramdom(int num)
        {
            //先统计一下有多少block
            //System.Console.WriteLine("openblocks_ramdom()中开始...");
            int num_blocks = 0;
            for(int i=0;i<_map.Length;i++)
            {
                for(int j=0;j<_map[i].Length;j++)
                {
                    if (_map[i][j] == STATUS_BLOCKS)
                        num_blocks++;
                }
            }

            //System.Console.WriteLine("openblocks_ramdom()中 1");
            //避免数量超了导致死循环
            if (num > num_blocks)
            {
                num = num_blocks;
            }

            //System.Console.WriteLine("openblocks_ramdom()中 2");
            //如果是全部挖空，则直接置0
            if (num == num_blocks)
            {
                for (int i = 0; i < _map.Length; i++)
                {
                    for (int j = 0; j < _map[i].Length; j++)
                    {
                        _map[i][j] = STATUS_OPEN;
                    }
                }
            }

            //System.Console.WriteLine("openblocks_ramdom()中 3");
            //print_map(_map);
            //System.Console.WriteLine("===3===  "+ num);
            //其他情况
            //Random r = new Random(DateTime.Now.Hour * 10000 + DateTime.Now.Minute * 100 + DateTime.Now.Second);
            Random r = new Random();
            int x_open = 0, y_open = 0;
            for (int j = 0/*j记录当前次的挖空数量*/; j < num; )
            {
                y_open = r.Next() % _map.Length;
                x_open = r.Next() % _map[0].Length;
                if(_map[y_open][x_open]==STATUS_BLOCKS)
                {
                    _map[y_open][x_open] = STATUS_OPEN;
                    j++;
                }
                
            }
            //System.Console.WriteLine("openblocks_ramdom()结束...");
        }
        /// <summary>
        /// 打印地图
        /// </summary>
        public void print_map()
        {
            if (_map == null) return;
            for(int i = 0; i < _map.Length; i++)
            {
                string tmp = "";
                for(int j = 0; j < _map[i].Length; j++)
                {
                    tmp += _map[i][j] + "\t";
                }
                System.Console.WriteLine(tmp);
            }
        }
        /// <summary>
        /// 打印地图
        /// </summary>
        public void print_map(int[][] v_map)
        {
            if (v_map == null) return;
            for (int i = 0; i < v_map.Length; i++)
            {
                string tmp = "";
                for (int j = 0; j < v_map[i].Length; j++)
                {
                    tmp += v_map[i][j] + "\t";
                }
                System.Console.WriteLine(tmp);
            }
        }
    }
}

Class_Percolation_union_find.cs

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace 渗漏问题
{
    /// <summary>
    /// 坐标点
    /// </summary>
    public struct point
    {
        public int x;
        public int y;
        public point(int _x,int _y)
        {
            x = _x;
            y = _y;
        }
        public static bool operator ==(point rec1, point rec2)
        {
            //注意我们调用Equals来判断是否相等。而不是在自己的函数中判、断。这是因为如果在自己的函数中判断。比如有rec2=null的情况。如果是这种情况。我们要判断if(rec2==null) {…}。其中rec2==null也是调用一个等号运算符，这里面有一个递归的过程，造成了死循环。  
            return rec1.x == rec2.x && rec1.y == rec2.y;
        }
        public static bool operator !=(point rec1, point rec2)
        {
            return !(rec1.x == rec2.x && rec1.y == rec2.y);
        }
    }

    public class Class_Percolation_union_find
    {
        public int ALGORITHMS_QUICK_FIND { get; private set; } = 0;
        public int ALGORITHMS_QUICK_UNION { get; private set; } = 1;

        private point[][] _map;
        /// <summary>
        /// 地图
        /// </summary>
        public point[][] map
        {
            get { return _map; }
            set { _map = _map = add_head_bottom(value); }
        }

        /// <summary>
        /// 方格为空
        /// </summary>
        //public const int STATUS_OPEN = 0;
        /// <summary>
        /// 方格阻塞
        /// </summary>
        //public const int STATUS_BLOCKS = -1;

        public point STATUS_BLOCKS = new point(-1, -1);
        
        public Class_Percolation_union_find()
        {
            _map = null;
        }

        public Class_Percolation_union_find(point[][] v_map)
        {
            //_map = v_map;
            _map = add_head_bottom(v_map);
        }

        /// <summary>
        /// 给地图添加全连通的顶部和底部，以方便计算
        /// </summary>
        /// <param name="v_map">原始地图</param>
        /// <returns>添加顶部底部后的地图</returns>
        private point[][] add_head_bottom(point[][] v_map)
        {
            int height = v_map.Length;
            int width = v_map[0].Length;
            point[][] tmp_map = new point[height + 2][];

            tmp_map[0] = new point[width];
            tmp_map[tmp_map.Length - 1] = new point[width];
            for(int i=0;i<tmp_map[tmp_map.Length - 1].Length; i++)
            {
                tmp_map[tmp_map.Length - 1][i] = new point(tmp_map.Length - 1, 0);
            }

            //顶部底部全部连通
            for (int i = 0; i < tmp_map[0].Length; i++)
            {
                tmp_map[0][i] = tmp_map[0][0];
                tmp_map[tmp_map.Length - 1][i] = tmp_map[tmp_map.Length - 1][0];

                point p0 = new point(0, i);
                point p1 = new point(tmp_map.Length - 1, i);
                pset_unblock.Add(p0);
                pset_unblock.Add(p1);
            }

            //中间依次复制
            for (int i = 0; i < v_map.Length; i++)
            {
                tmp_map[i + 1] = new point[v_map[0].Length];
                for (int j = 0; j < v_map[0].Length; j++)
                {
                    tmp_map[i + 1][j] = new point(i + 1, j);
                    //tmp_map[i + 1][j] = STATUS_BLOCKS;
                    //v_map[i].CopyTo(tmp_map[i + 1], 0);
                }
            }
            
            return tmp_map;
        }

        /// <summary>
        /// 存放已经打通的点
        /// </summary>
        HashSet<point> pset_unblock = new HashSet<point>();

        

        /// <summary>
        /// 打开一个点，并使用特定算法将该点与其周围的点联通
        /// </summary>
        /// <param name="algorithms"></param>
        public void open(int algorithms)
        {
            Random ran = new Random();
            point p;
            int a, b;
            do
            {
                a = ran.Next() % _map.Length;
                b = ran.Next() % _map[0].Length;
                p = new point(a, b);
            } while (pset_unblock.Contains(p) && pset_unblock.Count < (_map.Length * _map[0].Length));
            
            open(a, b, algorithms);
        }

        /// <summary>
        /// 打开一个点，并使用特定算法将该点与其周围的点联通
        /// </summary>
        /// <param name="x"></param>
        /// <param name="y"></param>
        /// <param name="algorithms"></param>
        public void open(int x,int y, int algorithms)
        {
            point p = new point(x, y);
            pset_unblock.Add(p);

            //point tmp_p = new point(-1, -1);
            
            if (p.x - 1 >= 0 && pset_unblock.Contains(_map[x - 1][y]))
            {
                //tmp_p = _map[x - 1][y];
                union_find(p, x - 1, y, algorithms);
            }
            if (p.x + 1< _map.Length && pset_unblock.Contains(_map[p.x + 1][y]))
            {
                //tmp_p = _map[p.x + 1][y];
                union_find(p, p.x + 1,y, algorithms);
            }
            if (p.y - 1 >= 0 && pset_unblock.Contains(_map[x][y - 1]))
            {
                //tmp_p = _map[x][y - 1];
                union_find(p, x,y-1, algorithms);
            }
            if (p.y + 1 < _map[0].Length && pset_unblock.Contains(_map[x][p.y + 1]))
            {
                //tmp_p = _map[x][p.y + 1];
                union_find(p, x, p.y + 1, algorithms);
            }
            
        }
        public void union_find(point p,int q_x, int q_y, int algorithms)
        {

            if (algorithms == ALGORITHMS_QUICK_FIND)
            {
                QuickFind qf = new QuickFind(_map);
                point q = _map[q_x][q_y];
                qf.union(p, q);
            }
            else if (algorithms == ALGORITHMS_QUICK_UNION)
            {
                QuickUnion qu = new QuickUnion(_map);
                point q = new point(q_x, q_y);
                qu.union(p, q);
            }

        }
        public bool isPercolation(int algorithms)
        {
            if (algorithms == ALGORITHMS_QUICK_FIND)
            {
                if (_map[0][0] == _map[_map.Length - 1][0])
                {
                    return true;
                }
                return false;
            }
            else if (algorithms == ALGORITHMS_QUICK_UNION)
            {
                QuickUnion qu = new QuickUnion(_map);
                return qu.connected(_map[0][0], _map[_map.Length - 1][0]);
            }
            return false;
        }
        
        /// <summary>
        /// 获取当前打通的砖块的数量
        /// </summary>
        /// <returns></returns>
        public int GetOpenedCount()
        {
            //int count = 0;
            //for(int i=1;i<_map.Length-1;i++)
            //{
            //    for(int j=0;j<_map[i].Length;j++)
            //    {
            //        if (_map[i][j] != STATUS_BLOCKS)
            //            count++;
            //    }
            //}
            //return count;
            return pset_unblock.Count - _map[0].Length * 2;
        }
        /// <summary>
        /// 打印地图
        /// </summary>
        public void print_map()
        {
            if (_map == null) return;
            for (int i = 0; i < _map.Length; i++)
            {
                string tmp = "";
                for (int j = 0; j < _map[i].Length; j++)
                {
                    tmp += _map[i][j] + "\t";
                }
                System.Console.WriteLine(tmp);
            }
        }
        /// <summary>
        /// 打印地图
        /// </summary>
        public void print_map(int[][] v_map)
        {
            if (v_map == null) return;
            for (int i = 0; i < v_map.Length; i++)
            {
                string tmp = "";
                for (int j = 0; j < v_map[i].Length; j++)
                {
                    tmp += v_map[i][j] + "\t";
                }
                System.Console.WriteLine(tmp);
            }
        }
    }

    public class QuickFind
    {
        private point[][] _map;

        public QuickFind(point[][] map)
        {
            _map = map;
        }

        public point find(point p)
        {
            return _map[p.x][p.y];
        }

        public void union(point p, point q)
        {
            point pID = find(p);
            point qID = find(q);

            if (pID == qID)
            {
                return;
            }

            for(int i = 0; i < _map.Length; i++)
            {
                for(int j = 0; j < _map[0].Length; j++)
                {
                    if(_map[i][j] == pID)
                    {
                        _map[i][j] = qID;
                    }
                }
            }
        }

    }

    public class QuickUnion
    {
        private point[][] _map;
        public QuickUnion(point[][] map)
        {
            _map = map;
        }

        private point find(point p)
        {
            while (p != _map[p.x][p.y])
                p = _map[p.x][p.y];
            return p;
        }

        public void union(point p, point q)
        {
            point pRoot = find(p);
            point qRoot = find(q);

            if (pRoot == qRoot)
                return;
            _map[pRoot.x][pRoot.y] = qRoot;

        }

        public bool connected(point p,point q)
        {
            return find(p) == find(q);
        }
    }
}

版权声明 本博客的文章除特别说明外均为原创，本人版权所有。欢迎转载，转载请注明作者及来源链接，谢谢。 本文地址: https://blog.ailemon.net/2018/07/29/algorithm-to-percolation/ All articles are under Attribution-NonCommercial-ShareAlike 4.0