|
发表于 2007-6-29 16:24:30
|
显示全部楼层
头文件Karatsuba.h:
---------------------------------------------------------------------------------------------
#define MAX_INDEX 1000
#define NULL_VALUE -1e15
-----------------------------------------------------------------------------------------------
cpp文件Karatsuba.cpp:
------------------------------------------------------------------------------------------------
#include <stdio.h>
#include "math.h"
#include "Karatsuba.h"
#include <afx.h>
//得到最高指数
int Get_index(double f[MAX_INDEX])
{
int i;
for(i=0;i<MAX_INDEX;i++)
{
if(fabs(f-NULL_VALUE)<1e-6)
return i-1;
}
return -1;
}
//得到指数的幂指数
int Get_n(double f[MAX_INDEX],double g[MAX_INDEX])
{
double n1,n2;
n1=(double)Get_index(f);
n2=(double)Get_index(g);
double maxn;
if(n1>n2)
maxn=n1;
else
maxn=n2;
int ret;
if(maxn==0)
return 0;
ret=(int)(log(maxn)/log(2))+1;
return ret;
}
//判断是否有值
bool IsNull(double a)
{
if(fabs(a-NULL_VALUE)<1e-6)
return true;
else
return false;
}
//加法
bool Add(double f[MAX_INDEX],double g[MAX_INDEX],double result[MAX_INDEX])
{
for(int i=0;i<MAX_INDEX;i++)
{
if(IsNull(f)&&IsNull(g))
{
return true;
}
else if(IsNull(f))
{
result=g;
}
else if(IsNull(g))
{
result=f;
}
else
{
result=g+f;
}
}
return true;
}
//减法
bool Sub(double f[MAX_INDEX],double g[MAX_INDEX],double result[MAX_INDEX])
{
for(int i=0;i<MAX_INDEX;i++)
{
if(IsNull(f)&&IsNull(g))
{
return true;
}
else if(IsNull(f))
{
result=0.0-g;
}
else if(IsNull(g))
{
result=f;
}
else
{
result=f-g;
}
}
return true;
}
//Karatsuba算法
bool Karatsuba(double f[MAX_INDEX],double g[MAX_INDEX], double result[MAX_INDEX])
{
int i;
int n;
n=Get_n(f,g);
if(n==0)
{
result[0]=f[0]*g[0];
return true;
}
double f1[MAX_INDEX],f0[MAX_INDEX],g1[MAX_INDEX],g0[MAX_INDEX];
double f0g0[MAX_INDEX],f1g1[MAX_INDEX],f0af1[MAX_INDEX],g0ag1[MAX_INDEX],
f0af1g0ag1[MAX_INDEX];
double res1[MAX_INDEX],res2[MAX_INDEX],res3[MAX_INDEX],res4[MAX_INDEX],res5[MAX_INDEX];
for(i=0;i<MAX_INDEX;i++)
{
f1=NULL_VALUE;
f0=NULL_VALUE;
g1=NULL_VALUE;
g0=NULL_VALUE;
f0g0=NULL_VALUE;
f1g1=NULL_VALUE;
f0af1=NULL_VALUE;
g0ag1=NULL_VALUE;
f0af1g0ag1=NULL_VALUE;
res1=NULL_VALUE;
res2=NULL_VALUE;
res3=NULL_VALUE;
res4=NULL_VALUE;
res5=NULL_VALUE;
}
int sub_n;
sub_n=1<<(n-1);
n=1<<n;
for(i=0;i<sub_n;i++)
{
f0=f;
g0=g;
}
i=0;
while(!IsNull(f[i+sub_n]))
{
f1=f[i+sub_n];
i++;
}
if(IsNull(f1[0]))
f1[0]=0.0;
i=0;
while(!IsNull(g[i+sub_n]))
{
g1=g[i+sub_n];
i++;
}
if(IsNull(g1[0]))
g1[0]=0.0;
Karatsuba(f0,g0,f0g0);
Karatsuba(f1,g1,f1g1);
Add(f0,f1,f0af1);
Add(g0,g1,g0ag1);
Karatsuba(f0af1,g0ag1,f0af1g0ag1);
i=0;
while(!IsNull(f1g1))
{
res1[i+n]=f1g1;
i++;
}
for(i=0;i<n;i++)
{
res1=0.0;
}
Sub(f0af1g0ag1,f0g0,res2);
Sub(res2,f1g1,res3);
i=0;
while(!IsNull(res3))
{
res4[i+sub_n]=res3;
i++;
}
for(i=0;i<sub_n;i++)
{
res4=0.0;
}
Add(res1,res4,res5);
Add(res5,f0g0,result);
//去掉最高位的0
i=0;
while(!IsNull(result))
{
i++;
}
i--;
while(result==0.0)
{
result=NULL_VALUE;
i--;
}
return true;
}
main()
{
double f[MAX_INDEX];
double g[MAX_INDEX];
double result[MAX_INDEX];
int i;
//数组清空
for(i=0;i<MAX_INDEX;i++)
{
f=NULL_VALUE;
g=NULL_VALUE;
result=NULL_VALUE;
}
//初始化
///f(x)=1.0x0+2.0x1
f[0]=1.0;
f[1]=2.0;
///g(x)=1.0x0+1.0x1+2.0x2
g[0]=1.0;
g[1]=1.0;
g[2]=3.0;
//乘法
Karatsuba(f,g,result);
//显示结果
CString msg;
CString string_res;
string_res="(";
i=0;
while(!IsNull(f[i+1]))
{
msg.Format("%.2fx%d+",f,i);
string_res+=msg;
i++;
}
msg.Format("%.2fx%d",f,i);
string_res+=msg+")*(";
i=0;
while(!IsNull(g[i+1]))
{
msg.Format("%.2fx%d+",g,i);
string_res+=msg;
i++;
}
msg.Format("%.2fx%d",g,i);
string_res+=msg+")=";
i=0;
while(!IsNull(result[i+1]))
{
msg.Format("%.2fx%d+",result,i);
string_res+=msg;
i++;
}
msg.Format("%.2fx%d\n",result,i);
string_res+=msg;
printf(string_res);
return 0;
}------------------------------------------------------------------------------------------------------------------------------------
输出结果:
(1.00x0+2.00x1)×(1.00x0+1.00x1+3.00x2)=1.00x0+3.00x1+5.00x2+6.00x3
字符串部分偷了点懒,用了MFC的CString类 :)
ps: 这个算法好好玩:)第一次见到这样算多项式的。
贴源代码的文字时,系统好像自动去掉了一些中括号,楼主还是下载附件看吧。 |
本帖子中包含更多资源
您需要 登录 才可以下载或查看,没有账号?注册
×
|