转自：http://www.cnblogs.com/kkgreen/archive/2011/06/12/2078668.html
参考 《算法设计与分析》 第四章 分治法 Anany Levitin著 翻译版 清华大学出版社
在上一篇文章中，介绍了分治策略的思想，主定理，以及几个用分治策略的经典案例。这一篇文章将继续探讨分治算法的其他应用，包括大整数乘法和Strassen矩阵乘法，最近点对问题和凸包问题这4个算法，一般来说常规的数据结构教程上不包括这些内容。
————————————————————————————————————————————————–
4.5 大整数乘法和Strassen矩阵乘法
两个算法都试图以增加少量的加法为代价，减少乘法运算的次数。
1）大整数乘法
书上讲述的是2个位数都为n位的数相乘，我这里扩展一下，写一个更为一般的方法，任意位数的2个数相乘：
对于任意位数的2个数相乘a * b，写成：
a = a1 * 10^(n1/2) + a0 —–n1为a的位数
b = b1 * 10^(n2/2) + b0 —–n2为b的位数
分治策略就是基于以上变换，将a，b写成前一半数字和后一半数字相加的形式，例如若a = 5423678，那么a1 = 542，a0 = 3678（注意若不是偶数截取较小一半）
这样a和b相乘就可以写为：a * b = { a1 * 10^(n1/2) + a0 } * { b1 * 10^(n2/2) + b0 }
展开后整理得： a * b = a1*b1 * 10^[ (n1+n2)/2 ] + a1*b0 * 10^(n1/2) + a0*b1 * 10^(n2/2) + a0*b0 四项
这样就很容易递归的来求a * b，如果你嫌分解后的数还太大，就可以继续分解。（你可以自己规定在何时结束递归）
实现方法：我们定义一个支持方法Mul(String s1,String s2)，用于在结束递归时（在本例中，我定义有一个数是1位时结束递归，直接用普通乘法）计算两个字符串的乘积（为了表示大数，用字符串来接受参数）。有了这个支持方法，分治递归实现两个大数乘法的实现如下：
public static long Mutiply(String a,String b) // 用字符串读入2个大整数
{
long result = 0 ;
if (a.length() == 1 || b.length() == 1 ) // 递归结束的条件
result = Mul(a,b);
else // 如果2个字符串的长度都 >= 2
{
String a1 = a.substring( 0 , a.length() / 2 ); // 截取前一半的字符串(较短的一半)
String a0 = a.substring(a1.length(), a.length()); // 截取后一半的字符串
// System.out.println(a1);
// System.out.println(a0);
String b1 = b.substring( 0 , b.length() / 2 );
String b0 = b.substring(b1.length(), b.length());
// 分治的思想将整数写成这样： a = a1 * 10^(n1/2) + a0, b = b1 * 10^(n2/2)，相乘展开得到以下四项
// 其中n1，n2为2个整数a，b的位数
result = ( long ) (Mutiply(a1,b1) * Math.pow( 10 , a0.length() + b0.length())
+ Mutiply(a1,b0) * Math.pow( 10 , a0.length()) + Mutiply(a0,b1) * Math.pow( 10 , b0.length())
+ Mutiply(a0,b0));
}
return result;
}
可以看到如果我们2个字符串都大于1位，那么我们就继续分解a和b。
完整的代码和测试的例子：
BigIntegerMutiply
95211154和903923的乘积：
86063551957142 86063551957142
true
注意，写这个算法仅仅是看看分治是怎样实现的，我的测试用例也算不上是大整数，真正的大整数在普通的计算机上会溢出，计算不出来。
一些研究人员的实验显示，在大于600位的整数相乘时，分治策略的性能才超过普通的算法
，显然，这么大的数字在常规的计算机上无法运算。
扩展一点，大整数加法怎么做？———-大整数乘法是为了提高效率，减少乘法运算次数而设计的，这里所说的大整数加法不是效率的问题，而是计算机无法表示的问题。例如对于100位甚至1000位的数字，计算机的表示范围显然无法达到，那么应该用数组来存放大整数的每一位，对于进位，直接向前一位加一即可。
2）Strassen矩阵乘法
理论上的分析就不多写了，涉及很多数学上的东西，计算麻烦，直接上书上的图：
分治策略将2个二阶矩阵采用下列方式来计算：
其中：
数一下，这样来计算2个二阶矩阵的乘法用了7次乘法，18次加法。而蛮力法用了8次乘法和4次加法。当然，这还不能体现出它的优越性，它的优越性表现在当矩阵的阶趋于无穷大时的渐进效率。
由上面的式子，根据矩阵的相关知识，对于分块矩阵，也可以写成上述形式：
用分治策略来计算矩阵乘法的方法如下：（支持方法较多，详见完整代码）
public static int [][] StrassenMul( int [][] a, int [][] b){ // a，b均是2的乘方的方阵
int [][] result = new int [a.length][a.length];
if (a.length == 2 ) // 如果a,b均是2阶的，递归结束条件
result = StrassMul(a,b);
else // 否则（即a，b都是4,8,16...阶的）
{
// a的四个子矩阵
int [][] A00 = copyArrays(a, 1 );
int [][] A01 = copyArrays(a, 2 );
int [][] A10 = copyArrays(a, 3 );
int [][] A11 = copyArrays(a, 4 );
// b的四个子矩阵
int [][] B00 = copyArrays(b, 1 );
int [][] B01 = copyArrays(b, 2 );
int [][] B10 = copyArrays(b, 3 );
int [][] B11 = copyArrays(b, 4 );
// 递归调用
int [][] m1 = StrassenMul(addArrays(A00,A11),addArrays(B00,B11));
int [][] m2 = StrassenMul(addArrays(A10,A11),B00);
int [][] m3 = StrassenMul(A00,subArrays(B01,B11));
int [][] m4 = StrassenMul(A11,subArrays(B10,B00));
int [][] m5 = StrassenMul(addArrays(A00,A01),B11);
int [][] m6 = StrassenMul(subArrays(A10,A00),addArrays(B00,B01));
int [][] m7 = StrassenMul(subArrays(A01,A11),addArrays(B10,B11));
// 得到result的四个子矩阵
int [][] C00 = addArrays(m7,subArrays(addArrays(m1,m4),m5)); // m1+m4-m5+m7
int [][] C01 = addArrays(m3,m5); // m3+m5
int [][] C10 = addArrays(m2,m4); // m2+m4
int [][] C11 = addArrays(m6,subArrays(addArrays(m1,m3),m2)); // m1+m3-m2+m6
// 也可以按照下列方法来求C
// C00 = addArrays(StrassenMul(A00,B00),StrassenMul(A01,B10));
// C01 = addArrays(StrassenMul(A00,B01),StrassenMul(A01,B11));
// C10 = addArrays(StrassenMul(A10,B00),StrassenMul(A11,B10));
// C11 = addArrays(StrassenMul(A10,B01),StrassenMul(A11,B11));
// 将四个子矩阵合并成result
Merge(result,C00, 1 );
Merge(result,C01, 2 );
Merge(result,C10, 3 );
Merge(result,C11, 4 );
}
return result;
}
完整代码：
Strassen
对给出的2个矩阵，乘积输出如下：
输出矩阵：
5 4 7 3
4 5 1 9
8 1 3 7
5 8 7 7
当然，矩阵的阶较大时才能体现出它的优势。
提示：其实二阶的时候结束递归直接蛮力就行了，而且你可以在四阶的时候就结束递归也行，取决于你的实现
————————————————————————————————————————————————–
4.6 用分治法解最近点对问题和凸包问题
问题描述在蛮力法里都介绍过，不再介绍了。
1）最近点对问题
在蛮力法里写了一个2重循环的蛮力算法，时间复杂度是n^2。
a，预排序！
b，在预排序基础上，将点递归的分为左右一半。
实现方式：按照x坐标升序对所求点集进行升序排序。然后递归的将点集分为左右一半，在点集中仅剩2或3个点时用蛮力法结束递归。
注意分治的策略：
划分子问题：每次将问题规模递归的划分为2个一半
合并子问题的解：分别求的2个一半中距离最近的2个点，假设在左一半是点A，点B，右边是点C，点D。再取其距离更小的2个作为暂时的最小值。注意此时AB或CD还不是最终要求的2个点，最小值可能还来自于一个点在左子集中一个点在右子集中。（位于左右子集分割边界上的2个点，一个在左子集中，一个在由子集中）
主方法如下：该方法接受一个点数组，返回其中距离最近的2个点。
public static Point[] getNearestPoints(Point[] Points) {
// 从一个点数组里面找到最近的两个点，并返回这两个点
Point[] result = new Point[ 2 ];
if (Points.length == 3 || Points.length == 2 ) // 递归结束的条件
result = getNear(Points);
else // 多于3个点，分治，分别找出两个子集合的最近点对，然后合并结果
{
Point[] left = Arrays.copyOfRange(Points, 0 , Points.length / 2 ); // 最后一个下标不包括
Point[] right = Arrays.copyOfRange(Points, Points.length / 2 ,
Points.length);
// 得到2个子集里面分别最短距离的2个点
Point[] result1 = getNearestPoints(left);
Point[] result2 = getNearestPoints(right);
double d1 = dPoints(result1[ 0 ], result1[ 1 ]);
double d2 = dPoints(result2[ 0 ], result2[ 1 ]);
// 忘了将result赋值
if (d1 <= d2)
result = result1;
else
result = result2;
// 合并结果：找到全局距离最短的两个点
double dmin = Math.min(d1, d2);
int x1 = left.length - 1 ; // 两个x的分界点
int x2 = x1 + 1 ;
// 在Points.length/2是一个整数时是错误的
// int x1 = Points[Points.length/2 - 1].x; // 两个x的分界点
// int x2 = Points[Points.length/2].x;
for ( int i = x1; i >= 0 ; i -- ) {
// if(x2 - Points[i].x > dmin) // 直接导致调试很久都不知道错在哪！！！！！！
if (Points[x2].x - Points[i].x > dmin)
break ;
else
// for(int j = Points.length/2;j < Points.length;j++)
for ( int j = x2; j < Points.length; j ++ ) {
// System.out.println(Points[j].y);
// if(Points[j].x - x1 > dmin)
if (Points[j].x - Points[x1].x > dmin)
break ;
else {
double temp = dPoints(Points[i], Points[j]);
// System.out.println(temp);
if (temp < dmin) {
dmin = temp;
result[ 0 ] = Points[i];
result[ 1 ] = Points[j];
}
}
}
}
}
return result;
}
注意要得到左右边界的值，然后从左右边界处往两边扩展，看看有没有更小的值。一开始再这个地方犯了错误，见注释。
其他支持方法见完整代码：
NearestPoint
输出距离最近的两个点是：
(7,8) (7,9)
这个算法递归的关系式为：
根据主定理，时间复杂度为n * log（n），注意这仅仅是对于划分和合并子问题的解的递推式。还有一个预排序，由于预排序时间复杂也是n * log（n），故它不影响整体的复杂度。
提示：这个例子告诉我们，有时候合并子问题的解并不仅仅是简单的把子问题的解合并起来，它需要一些进一步的处理。最终的解可能并不在子问题的解中，但它来自于对子问题的解进一步的处理结果。
2）凸包问题
描述比较复杂，简单说下，先预排序，预排序后最左和最右的点肯定是凸包中的点。然后可以递归的从内向外扩展凸包，在当前直线的2侧寻找最高点，最高点肯定在凸包中，这里涉及到一些数学知识：
a，首先定义射线p1到p2的左侧：若p1 p2 p构成的顺序是逆时针，称p在射线的左侧
b，三角形p1 p2 p3的面积等于下列行列式的一半：
仅当p3在射线p1p2左侧时这个值才为正。
由此我们很容易求p1，p2左侧的最高点（离直线最远的点，这个点即凸包向外扩展得到的新顶点），得到一个最高点后，就得到了2条新边，继续向外扩展，如图：
注意我只定义了左侧，其实左右都一样，你把射线方向换一下，处理右侧不就等于处理左侧了吗，免得搞得那么麻烦。详见代码：
package Section4;
/* 第4章 分治法 凸包问题的分治解法 */
import java.util.Arrays;
import java.util.Iterator;
import java.util.LinkedList;
import java.util.Stack;
public class ConvexHull {
/**
* @param args
*/
public static void main(String[] args) {
// TODO Auto-generated method stub
Point[] Points = new Point[ 15 ];
Points[ 0 ] = new Point( - 5 , 7 );
Points[ 1 ] = new Point( 3 , - 6 );
Points[ 2 ] = new Point( 5 , 4 );
Points[ 3 ] = new Point( - 5 , - 5 );
Points[ 4 ] = new Point( 1 , 7 );
Points[ 5 ] = new Point( 6 , 0 );
Points[ 6 ] = new Point( 0 , 0 );
Points[ 7 ] = new Point( - 5 , 0 );
Points[ 8 ] = new Point( 3 , - 2 );
Points[ 9 ] = new Point( 3 , 4 );
Points[ 10 ] = new Point( 1 , 6 );
Points[ 11 ] = new Point( 5 , 3 );
Points[ 12 ] = new Point( - 4 , - 5 );
Points[ 13 ] = new Point( - 3 , 6 );
Points[ 14 ] = new Point( 2 , 5 );
Arrays.sort(Points); // 预排序处理
LinkedList < Point > list = new LinkedList < Point > ();
for ( int i = 0 ;i < Points.length;i ++ )
list.add(Points[i]); // list存放全部的顶点
LinkedList < Point > result = getConvexHulls(list); // result用来存放最终的结果顶点
System.out.println( " 一共有 " + result.size() + " 个顶点, " + " 凸包的顶点是： " );
Iterator it = result.iterator();
while (it.hasNext())
{
Point next = (Point) it.next();
System.out.print( " ( " + next.x + " , " + next.y + " ) " + " " );
}
}
public static LinkedList < Point > getConvexHulls(LinkedList < Point > list){
// 将凸包顶点以result链表返回
LinkedList < Point > result = new LinkedList < Point > ();
Point temp1 = list.removeFirst();
Point temp2 = list.removeLast();
result.add(temp1);
result.add(temp2);
// 递归的处理temp1 ---> temp2左右两侧的点
dealWithLeft(temp1,temp2,result,list);
dealWithLeft(temp2,temp1,result,list); // 注意每次要将result带着，存放结果集
return result;
}
private static void dealWithLeft(Point p1,Point p2,LinkedList result,LinkedList list){
// 递归的处理p1，p2构成的射线左边的点
Iterator it = list.iterator();
// 找出左边最高的点Pmax
Point Pmax = null ;
int max = 0 ;
while (it.hasNext())
{
Point next = (Point) it.next();
int x1 = p1.x,y1 = p1.y;
int x2 = p2.x,y2 = p2.y;
int x3 = next.x,y3 = next.y;
// int max = 0; // 小小的一个错误啊！！！！！！！
int compute = x1 * y2 + x3 * y1 + x2 * y3 - x3 * y2 - x2 * y1 - x1 * y3;
if (compute > max)
{
max = compute;
Pmax = next;
}
}
// 又找到了一个顶点
if (Pmax != null )
{
result.add(Pmax);
list.remove(Pmax);
// 递归
dealWithLeft(p1,Pmax,result,list);
dealWithLeft(Pmax,p2,result,list);
}
}
private static boolean onLeft(Point target,Point p1,Point p2){
// 判断target是否在p1--->p2射线的左侧
int x1 = p1.x,y1 = p1.y;
int x2 = p2.x,y2 = p2.y;
int x3 = target.x,y3 = target.y;
int compute = x1 * y2 + x3 * y1 + x2 * y3 - x3 * y2 - x2 * y1 - x1 * y3;
if (compute > 0 )
return true ;
else
return false ;
}
}
以上测试的14个点是精心构造的包含了各种特殊情况的14个点（包括有3点共线，有点的横坐标相等，等等），运行结果如下，有兴趣的可以在坐标系里画一下。
一共有 6 个顶点, 凸包的顶点是：
(-5,7) (6,0) (1,7) (5,4) (-5,-5) (3,-6)
提示：这个代码写的不好，其实每次判断一侧后，就可以将另一侧的点排除掉，这里还没有这样去做。
————————————————————————————————————————————————–
总结：
分治的策略
怎么划分子问题
怎么用递归
规模小的时候结束递归
怎么合并结果
