天天看点

多项式相乘的高效分治算法一

        多项式相乘最容易想到的方法就是采用2个for循环,依次遍历相乘累加即可。但时间复杂度为O(n2)。那么是否可用其他的方法来降低遍历比较的次数呢?可用分治法来进行对应切分,然后相乘,再相加。但是这好像并不能减少遍历的次数,因为与之前的合并排序不同的是,合并排序中,如果左边的数若比右边的数小,则左边的左边就不需要与右边的所有数进行比较,因为数是sorted,所以默认比右边的小,无需再比较;而两个多项式相乘,不存在排序,而且每个数都必须与对面的每个数相乘,因为数都是不同的,所以必须每个数与另一个多项式的每一项相乘(注:一个加号连接两个项),那么根据以上的分析,无论如何都要进行N2次乘法,这样即使采用分治法,最后的时间复杂度也为O(n2)。这个问题让我不得其解,因为不管怎么运用何种算法都不可能减少相乘的次数。然后高斯采用了一个出其不意的思想方法来减少了相乘的次数,变4次大相乘为3次大相乘,而每个数又都相乘了。

         该方法的理论依据是:

  1. 在计算上,乘法比加法更耗时,所以可用加法来代替乘法;
  2. 该方法来源于反思考宏观的推理。

        因此,多项式相乘高效方法来源于高斯的推理。那么首先从多项式的4次分治相乘的编程开始。如何采用分治法进行分割相乘呢?基本的步骤如下:

      例如有多项式A=[2,3,5,4]和B=[3,1,3,5]  ,注:多项式只存系数,系数的位序代表了该项的幂次。位序==幂次是编程实现的关键。如何关键呢?

     Step1. 将A,B分别从中间分割成A0和A1,B0和B1;

     Step2.  A*B=A0*B0+A0*B1+A1*B0+A1*B1,对这四项进行递归分割相乘,然后叠加;

      首先,根据分治法,编写一个主体的递归程序如下:

vector<double> DivideMerge::polynomialsMultiply(vector<double>&numA, int lowA, int highA, vector<double>&numB, int lowB, int highB)
{
	vector<double>numXY,numXYZ;
	vector<double>numX, numY, numZ;                                                            

	if (lowA < highA && lowB<highB)
	{
		int midA = (highA + lowA) / 2;
		int midB = (highB + lowB) / 2;

		numX = polynomialsMultiply(numA, lowA, midA, numB, lowB, midB);

		vector<double> numa, numb;
		int dy;                 //交叉相乘递归
		numa = polynomialsMultiply(numA, lowA, midA, numB, midB + 1, highB);
		numb = polynomialsMultiply(numA, midA + 1, highA, numB, lowB, midB);
		numY = crossAdd(numa, midA, highB, numb, highA, midB,dy);  //A0*B1+A1*B0

		int dxy;
		numXY = crossAdd(numX, midA, midB, numY, 0, dy, dxy);     //A0B0+(A0B1+A1B0)
		numZ = polynomialsMultiply(numA, midA + 1, highA, numB, midB + 1, highB);              

		int dxyz;
		numXYZ = crossAdd(numXY, 0, dxy, numZ, highA, highB,dxyz);
	}
else
	{
	 numXYZ.push_back(numA[lowA] * numB[lowB]);                                                 //调用的这个函数返回的是一个void类型
	}
}
           

      这个递归中只考虑if (lowA < highA && lowB<highB),当递归到终止条件lowA==highA或lowB==highB时,若用此一个if语句,当时出现了两个bug,一是地址越界:即numa = polynomialsMultiply(numA, lowA, midA, numB, midB + 1, highB); 其中midB+1==1,highB==0,则递归不会出错,但是在终止条件返回时则出现地址越界,因为这时的lowB>hgihB>numB.size()-1,这个错误是由于采用一个统筹的if语句引起的,前面numA、numB都是low<high,而mid肯定是比high小,所以在分治递归上没有问题,但是当一个满足if条件而另一个不满足时,这时仍然采用之前的mid计算方法,则就会使得mid==high,而分治递归时,mid+1>high导致终止条件后面的地址越界。当然可用条件判断来根据实际的low和high控制是否调用交叉相乘递归。但是本文采用了另外一种逻辑更清晰的结构,即归并排序的while循环结构,整个递归编程如下:

vector<double> DivideMerge::polynomialsMultiply(vector<double>&numA, int lowA, int highA, vector<double>&numB, int lowB, int highB)
{
	vector<double>numXY,numXYZ;
	vector<double>numX, numY, numZ;              

	if (lowA < highA && lowB<highB)
	{
		int midA = (highA + lowA) / 2;
		int midB = (highB + lowB) / 2;

		numX = polynomialsMultiply(numA, lowA, midA, numB, lowB, midB);

		vector<double> numa, numb;
		int dy;
		numa = polynomialsMultiply(numA, lowA, midA, numB, midB + 1, highB);
		numb = polynomialsMultiply(numA, midA + 1, highA, numB, lowB, midB);
		numY = crossAdd(numa, midA, highB, numb, highA, midB,dy);

		int dxy;
		numXY = crossAdd(numX, midA, midB, numY, 0, dy, dxy);
		numZ = polynomialsMultiply(numA, midA + 1, highA, numB, midB + 1, highB);               

		int dxyz;
		numXYZ = crossAdd(numXY, 0, dxy, numZ, highA, highB,dxyz);
	
	}
	else if (lowA<highA)
	{
		int midA = (highA + lowA) / 2;
		numX = polynomialsMultiply(numA, lowA, midA, numB, lowB, highB);                          //只有一个就不用交叉
		numZ = polynomialsMultiply(numA, midA + 1, highA, numB, lowB, highB);

		int dxz;
		numXYZ = crossAdd(numX, midA, highB, numZ, highA, highB,dxz);
	}
	else if (lowB < highB)
	{
		int midB = (lowB+highB) / 2;
		numX = polynomialsMultiply(numA, lowA, highA, numB, lowB, midB);                          //只有一个就不用交叉
		numZ = polynomialsMultiply(numA, lowA, highA, numB, midB+1, highB);

		int dxz;
		numXYZ = crossAdd(numX, highA, midB, numZ, highA, highB,dxz);
	}
	else
	{
		numXYZ.push_back(numA[lowA] * numB[lowB]);                                                 //调用的这个函数返回的是一个void类型
	}

	return numXYZ;
}
           

        其次,4次相乘依次结束后,就需要对同幂次的项进行叠加,如何对所有递归回来的项进行正确的叠加是实现该算法的核心和关键,前前后后花了2个星期的思考才彻底解决了这个问题:(numX,numY和numZ 分别表示A0*B0,A0B1+B0A1和A1*B1递归相乘返回的结果)

        (思考1)最先想到的是对numX的最后一项与numY的第一项相加,这样仅仅试用对3×3的多项式,但当多项式较多时,numX的最后一项的幂次就与numY的幂次就不对应了。

        (思考2)根据4×4的例子,采用n/2的划分来叠加,即元素项的总个数n/2对应numC = numX+numY,n/2+n/2对应于numC+numZ,测试用例得出,该方法也不能普遍试用所有的情况,当项数较多时,结果叠加不正确,也不适用于1×3的情况。

       那么,到底应该如果叠加呢,当多项式的项数增多时,且两个多项式的项数不一致时,得出的numX,numY和numZ的情况仅仅通过n/2的方法进行叠加,你会发现会出现叠加错误,根本原因是这里根据相乘结果的位置来叠加,而numX 和numY位置上的元素的幂次不能保证对应相等,比如numX最后一位为3次幂,而numY的第1项可能使2次幂,这样就没有对应幂次相加。到底该如何通用的叠加呢?这个问题一直思考了1个星期。

         最终回归到问题的本身,依据什么叠加,依据是同幂次才能相加,那么只要保证相加的数是同幂次,则就能保证叠加的正确性。那么如何保证相加的数是同幂次呢?就要根据幂次在程序中的定义来编程。还记得本文开头加黑的位序==幂次是编程实现的关键,多项式的存储是按照幂次的高低依次存储在vector中的。因此只有根据总体的low和high的值就能准确计算出整个递归过程中的结果项和中间项的幂次,即通过vector的下标就能遍历该项,这样再次印证了“数据善置,算法自成”,对数据能掌控的存和取,就能对数据进行心中想要的计算。

      幂次的数据叠加步骤如下:

      Step1. 对于给定的两个vector数据项numa 和numb 的最高次幂,就能推到出numa中前面项的幂次,因为幂次都是按顺序存储在numa中的,若最高次幂为n,则前一项为n-1,而第一项则是n-numa.size()+1,这样编程就能很容易实现:若numa的首项的幂次小于相对应的numb的首项。则当i的循环值还没有到达numa和numb的同幂次项时,则先让numa[i]进vector,这里有一个同幂次基准下标就是numb的首项幂次,所以可以通过if判断下标i(i定义为numa的下标)若小于低次幂的元素个数(低次幂的个数=numa的首项幂次-numa的首项幂次),就直接numa[i]进vector;而两项相加则注意i的值小于numa.size(),避免下标越界。

最后的同幂次叠加算法如下:(注:该算法不针对特定的某个项,是一个普适的只要知道元素项和最高次幂,因此不管是累加还是交叉叠加都可以)

std::vector<double> DivideMerge::crossAdd(vector<double> &numa, int midA, int highB, vector<double>&numb, int highA, int midB,int &d)  //这里的int类型的参数都是最高次幂
{
	vector<double> numc;
	int i = 0;
	int j = 0;
	int na = midA + highB;                             //numa的最高幂次 2
	int nb = highA + midB;                             //numb的最高幂次 3
	int diffA = na - numa.size() + 1;
	int diffB = nb - numb.size() + 1;
	int flag = 0;                                     //相等就为0

	d = na > nb ? na : nb;                            //原则:返回最高次幂

	if (diffA > diffB)
	{
		flag = 1;
	}
	else if (diffA < diffB)
	{
		flag = 2;
	}

	while (i < numa.size() || j < numb.size())
	{
		switch (flag)
		{
			case 0:{
					   if (i < numa.size() && j < numb.size()) 
					   {
							numc.push_back(numa[i] + numb[j]);
							i++;
							j++;
					   }
					   else if (i<numa.size())
					   {
						   numc.push_back(numa[i]);
						   i++;
					   }
					   else
					   {
						   numc.push_back(numb[j]);
						   j++;
					   }
					  
			}continue;

			case 1:{
					   if (j < diffA - diffB)
					   {
						   numc.push_back(numb[j]);
						   j++;
					   }
					   else if (j < numb.size())
					   {
						   numc.push_back(numa[i] + numb[j]);
						   i++;
						   j++;
					   }
					   else
					   {
						   numc.push_back(numa[i]);
						   i++;
					   }

			}continue;

			default:{

						if (i < diffB - diffA)
						{
							numc.push_back(numa[i]);
							i++;
						}
						else if (i<numa.size())
						{
							numc.push_back(numa[i] + numb[j]);
							i++;
							j++;
						}
						else
						{
							numc.push_back(numb[j]);
							j++;
						}
			}continue;
		}//switch
	}//while

	return numc;
}
           

       总结:本文基于分治法的多项式相乘算法的时间复杂度仍为O(n2),而降低时间复杂度的方法就是采用高斯方法,即减少相乘的次数,这个将在下一篇文中写到。

      在叠加算法中,第一次采用了switch... case ... continue的搭配,这是由于外层有循环条件,因此具体用break 还是continue 要根据实际的功能需求而定,break跳出整个循环,而continue则仅仅是终止接下来语句的执行,并不一定是循环,思维一定要开阔。在测试中,仅仅针对于项数2×3,4×4,1×3来测试 或8×8等,这样的测试用例没有全考虑涵盖,所以最后根据1×2、1×3、1×4... 2×2、2×3、2×4 等等,根据遍历的方法既能都考虑到,又避免重复的全测试方法。

继续阅读