技术图文：如何改进算法的运行效率？

背景

前段时间，一位好友发给我如下的文件：

每个CSV文件中的数据由三个属性组成，第一个属性为ID，第二个属性为X坐标，第三个属性为Y坐标。由于是二维数据，可以绘制出每个文件的散点图，把这些散点图连接起来，构成如下的GIF图像。类似于视频中的帧，每一帧对应一个CSV文件。

要做的事情是什么呢？

在每一帧数据中，快速找到以每个数据点为圆心，r = 0.004为半径的范围内包含的其它数据点。

技术分析

刚拿到这个问题的第一个反应就是利用暴力的方法，循环计算每个点与其它点的距离，把小于r的记录下来就可以了。可是这样做会耗费大量的计算时间。

如何改进呢，首先想到的是把计算欧氏距离的根号去掉，只计算欧氏距离的平方，这样可以避免开根号，从而提升运算效率。于是写了“暴力求解算法”的代码，发现执行效率很低，执行 10次 的平均时间约为 2754毫秒。

static float Distance(PointF p1, PointF p2)
{float d = (p1.X - p2.X)*(p1.X - p2.X) + (p1.Y - p2.Y)*(p1.Y - p2.Y);return d;
}for (int i = 0; i < lst.Count; i++)
{PointF p1 = lst[i].Original;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].Original;if (Distance(p1, p2) <= r * r){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}
}

怎么办呢？

我们知道平面上的点属于R2空间，是基于e1=(1,0),e2=(0,1)这组基的，为了方便比较我们可以做一个坐标变换，让这些数据基于a1=(r,0),a2=(0,r)这组基。这样给定一个点，我们就能很快判断出是否超出单位圆的范围。由于做基变换的过渡矩阵是可逆的，新的坐标和初始坐标可以相互线性表示，这样就把原问题简化了。于是写了“利用坐标变换进行计算”的代码，发现效率的提升并不明显，仅仅缩短了一半左右的执行时间，执行 10次 的平均时间约为 1296毫秒。

for (int i = 0; i < lst.Count; i++)
{PointF p1 = lst[i].AfterConversion;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].AfterConversion;if (p1.X + 1 < p2.X || p1.X - 1 > p2.X || p1.Y + 1 < p2.Y || p1.Y - 1 > p2.Y)continue;if (Distance(p1, p2) <= 1){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}
}

怎么办呢？

我们发现如果对数据的某一个维度进行排序，这样可以减少计算距离的次数，于是我们先对数据中的x维度进行排序，然后再来计算距离。于是写了“利用坐标变换和排序进行计算”的代码，发现效率提升了90%以上，执行 10次 的平均时间约为 235毫秒。

for (int i = 0; i < lst.Count; i++)
{PointF p1 = lst[i].AfterConversion;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].AfterConversion;if (p1.X + 1 < p2.X)break;if (Distance(p1, p2) <= 1){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}
}

这里还有一个插曲，我曾经想利用KD Tree的方式来建立每一帧数据的空间索引结构，这是加速K 邻近算法的方法。由于每帧数据都是变化的，虽然查找的速度很快，但构建这个Tree的过程需要耗费大量的时间，于是放弃了这种方法。

代码实现

记录时间的代码如下：

public class TimeRecord
{[DllImport("Kernel32.dll")]private static extern bool QueryPerformanceCounter(out long lpFrequency);[DllImport("Kernel32.dll")]private static extern bool QueryPerformanceFrequency(out long lpFrequency);private long _startTime, _stopTime;private readonly long _freq;/// <summary>/// 创建一个TimeRecord类的新实例/// </summary>public TimeRecord(){_startTime = 0;_stopTime = 0;if (QueryPerformanceFrequency(out _freq) == false){throw new Win32Exception();}}/// <summary>/// 开始计时/// </summary>public void Start(){Thread.Sleep(0);QueryPerformanceCounter(out _startTime);}/// <summary>/// 停止计时/// </summary>public void Stop(){QueryPerformanceCounter(out _stopTime);}/// <summary>/// 获取从Start()->Stop()中间的精确时间间隔/// 单位：毫秒/// 计数次数/计数频率/// </summary>public double DurationMs{get{return (double)(_stopTime - _startTime) * 1000 / _freq;}}
}

构造数据存储结构的代码如下：

public class Data
{/// <summary>/// 编号/// </summary>public int Id;/// <summary>/// 原先点/// </summary>public PointF Original;/// <summary>/// 转换后的点/// </summary>public PointF AfterConversion;/// <summary>/// 范围内的点/// </summary>public List<PointF> LstWithin = new List<PointF>();
}

暴力求解算法的代码如下：

static List<Data> LoadData(string[] fileNames)
{List<Data> result = new List<Data>();foreach (string file in fileNames){string[] strs = File.ReadAllLines(file);foreach (string str in strs){string[] temp = str.Split(new char[] {','});Data data = new Data{Id = int.Parse(temp[0]),Original = new PointF(float.Parse(temp[1]), float.Parse(temp[2]))};result.Add(data);}}return result;
}static float Distance(PointF p1, PointF p2)
{float d = (p1.X - p2.X)*(p1.X - p2.X) + (p1.Y - p2.Y)*(p1.Y - p2.Y);return d;
}static double Test01(string[] files)
{TimeRecord tr = new TimeRecord();tr.Start();float r = 0.004f;List<Data> lst = LoadData(files);for (int i = 0; i < lst.Count; i++){PointF p1 = lst[i].Original;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].Original;if (Distance(p1, p2) <= r * r){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}}tr.Stop();return tr.DurationMs;
}

利用坐标变换进行计算的代码如下：

static List<Data> LoadData(string[] fileNames,float d)
{List<Data> result = new List<Data>();foreach (string file in fileNames){string[] strs = File.ReadAllLines(file);foreach (string str in strs){string[] temp = str.Split(new char[] { ',' });Data data = new Data();data.Id = int.Parse(temp[0]);float x = float.Parse(temp[1]);float y = float.Parse(temp[2]);data.Original = new PointF(x, y);data.AfterConversion = new PointF(x*d, y*d);result.Add(data);}}return result;
}static double Test02(string[] files)
{TimeRecord tr = new TimeRecord();tr.Start();List<Data> lst = LoadData(files, 250);for (int i = 0; i < lst.Count; i++){PointF p1 = lst[i].AfterConversion;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].AfterConversion;if (p1.X + 1 < p2.X || p1.X - 1 > p2.X || p1.Y + 1 < p2.Y || p1.Y - 1 > p2.Y)continue;if (Distance(p1, p2) <= 1){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}}tr.Stop();return tr.DurationMs;
}

利用坐标变换和排序进行计算的代码如下：

static double Test03(string[] files)
{TimeRecord tr = new TimeRecord();tr.Start();List<Data> lst = LoadData(files, 250);lst = lst.OrderBy(a => a.AfterConversion.X).ToList();for (int i = 0; i < lst.Count; i++){PointF p1 = lst[i].AfterConversion;for (int j = i + 1; j < lst.Count; j++){PointF p2 = lst[j].AfterConversion;if (p1.X + 1 < p2.X)break;if (Distance(p1, p2) <= 1){lst[i].LstWithin.Add(p2);lst[j].LstWithin.Add(p1);}}}tr.Stop();return tr.DurationMs;
}

测试代码如下：

delegate double Test(string[] files);static double Run(Test test)
{string[] files = new string[]{"40000weizhi.csv","41000weizhi.csv","41500weizhi.csv","42000weizhi.csv"};double[] r = new double[10];for (int i = 0; i < 10; i++){r[i] = test(files);}return r.Average();
}static void Main(string[] args)
{double r1 = Run(Test01);double r2 = Run(Test02);double r3 = Run(Test03);Console.WriteLine("方法1耗费时间:{0}", r1);Console.WriteLine("方法2耗费时间:{0}", r2);Console.WriteLine("方法3耗费时间:{0}", r3);
}

测试结果如下：

总结

这篇图文记录了我尝试提升算法执行效率的过程，不知道能否达到这位朋友的要求，毕竟比起他们现在的算法，缩短了90%以上的时间。今天就这样吧！See You！大家如果有什么好的方法，欢迎给我留言啊，我们一起把这个问题解决掉。

往期活动

LSGO软件技术团队会定期开展提升编程技能的刻意练习活动，希望大家能够参与进来一起刻意练习，一起学习进步！

• Python基础刻意练习活动即将开启，你参加吗？
• Task01：变量、运算符与数据类型
• Task02：条件与循环
• Task03：列表与元组
• Task04：字符串与序列
• Task05：函数与Lambda表达式
• Task06：字典与集合
• Task07：文件与文件系统
• Task08：异常处理
• Task09：else 与 with 语句（1day）
• Task10：类与对象
• Task11：魔法方法
• Task12：模块