这几天在折腾算阶乘的方法，想到了因式分解法算，但因式分解又不是真的因式分解，而是找到一个素数就看这个阶乘里面总共会出现多少次，把N！里面所有的素数和出现次数都算出来以后再计算就会很快了！
首先给初始的素数{2,3,5,7}，这里面还有一个技巧，就是素数5出现的次数，就是最终有多少个0的个数，因为2出现的次数一定是比5出现的次数多的。
而求N！中素数q出现的次数其实很简单，就一个while循环就好
count=0;
m=N;
while(m>=q)
{
m/=f;
count+=m;
}
经过这几行代码后count里面就是N！中素数q出现的次数。
通过计算10000！中列举了最前面的几个素数和各自出现的次数如下：
1:2,9995
2:3,4996
3:5,2499
4:7,1665
5:11,998
6:13,832
7:17,624
8:19,554
9:23,452
因为2出现9995次，而5出现2499次，那么10000！最后面一定是2499个0，所以5不用参加计算，而2参加计算的是（9995-2499）次方。
为了计算大数专门写了一个类bigint
其中最主要的成员就一个：
private:
vector<int> numbers;
这个成员保存一堆的int，每个int让它小于10的9次方，这样留一点加法不越界。乘法在64位不越界。#define MAXINT (1000000000)
三个构造函数：
bigint()
{
}
bigint(const int value)
{
numbers.push_back(value);
}
bigint(const bigint& b)
{
for(int i=0;i<b.len();i++)
{
numbers.push_back(b.getnum(i));
}
}
to_string函数，除最高位以外不满9位的前面补0：
string to_string()
{
int i=len();
string r,s0="000000000",s1;
i--;
r=std::to_string(numbers[i]);
while(i>0)
{
i--;
s1 = std::to_string(numbers[i]);
if(s1.length()<9)
s1=s0.substr(0,9-s1.length())+s1;
r+=s1;
}
return r;
}
重载*乘法：
bigint operator*(const bigint& b)
{
bigint r;
unsigned long long rare=0;
int i,j,k,ib,ie,m,n;
m=len()-1;
n=b.len()-1;
for(k=0;k<=m+n;k++)
{
unsigned long long ull=rare,ull2;
rare=0;
ib=(k>n)?k-n:0;
ie=(m>k)?k:m;
for(i=ib;i<=ie;i++)
{
j=k-i;
ull2=getnum(i);
ull2*=b.getnum(j);
ull+=ull2;
rare+=ull/MAXINT;
ull=ull%MAXINT;
}
r.setnum(k,(int)ull);
}
while(rare>0)
{
r.setnum(k,rare%MAXINT);
rare=rare/MAXINT;
k++;
}
return r;
}
被乘数[0,m]，乘数[0,n]，乘后变成[0,m+n]，有进位还会增加。第一个循环是以k从0到m+n，执行i是被乘数的变量，它的循环开始和结束有些关键，看下代码要好好想想。
乘方：
bigint power(int n)
{
bigint p(*this);
static int masks[]={0x00,
0x01,
0x02,0x02,
0x04,0x04,0x04,0x04,
0x08,0x08,0x08,0x08,0x08,0x08,0x08,0x08};
int mask=0,caln=n;
if(n==0)
return 1;
if(caln&0x7fff0000)
{
mask=16;
caln>>=16;
}
caln&=0xffff;
if(caln&0xff00)
{
mask+=8;
caln/=0x100;
}
caln&=0xff;
if(caln&0xf0)
{
mask+=4;
caln/=0x10;
}
caln&=0xf;
mask=(masks[caln]<<mask);
mask>>=1;
while(mask)
{
p=p*p;
if(n&mask)
p=p*(*this);
mask>>=1;
}
return p;
}
代码看上去很多，其实关键的代码就后面这个while循环，在此之前都是在找罪高位的1，
当找到这个高位以后，继续往后面看，有后面一位首先自己平方一下，如果下面一位是1还需要再乘上一个（*this）。
读写长度和数字代码如下：
int len() const
{
return numbers.size();
}
int getnum(int i) const
{
if(i<0||i>=numbers.size())
return 0;
return numbers[i];
}
bool setnum(int i,int value)
{
if(i<0||i>numbers.size())
return false;
if(i==numbers.size())
numbers.push_back(value);
else
numbers[i]=value;
return true;
}
主函数：
因式分解：
vector<int> factors={2,3,5,7};
vector<int> counts;
N！中包含这几个素数的次数找出来并把它push_back到counts里面
for(i=0;i<factors.size();i++)
{
c=0;
f=factors[i];
m=N;
while(m>=f)
{
m/=f;
c+=m;
}
counts.push_back(c);
cout<<counts.size()<<":"<<f<<","<<c<<endl;
}
要找后面的素数就从11，开始找奇数，因为11以后的偶数肯定不是素数。一边找一边在素数vector里面插入，先用一个bool bf=true;表示是素数，当下面一个循环中发现不是素数再把bf置false，当查询的素数平方大于要查的数的时候表面它本身就是素数也可以停下来了。
查找到一个素数后马上计算这个素数有多少次，并把这两个值插在最后面。
for(i=11;i<=N;i+=2)
{
bool bf=true;
for(j=0;j<factors.size();j++)
{
f=factors[j];
if(f*f>i) break;
if(i%f==0)
{
bf=false;
break;
}
}
if(bf)
{
c=0;
m=N;
f=i;
while(m>=f)
{
m/=f;
c+=m;
}
factors.push_back(i);
counts.push_back(c);
cout<<counts.size()<<":"<<f<<","<<c<<endl;
}
}
完整代码如下：
#include <string>
#include <iostream>
#include <vector>
#include <ctime>
using namespace std;
#define MAXINT (1000000000)
class bigint{
public:
bigint()
{
}
bigint(const int value)
{
numbers.push_back(value);
}
bigint(const bigint& b)
{
for(int i=0;i<b.len();i++)
{
numbers.push_back(b.getnum(i));
}
}
bigint operator*(const bigint& b)
{
bigint r;
unsigned long long rare=0;
int i,j,k,ib,ie,m,n;
m=len()-1;
n=b.len()-1;
for(k=0;k<=m+n;k++)
{
unsigned long long ull=rare,ull2;
rare=0;
ib=(k>n)?k-n:0;
ie=(m>k)?k:m;
for(i=ib;i<=ie;i++)
{
j=k-i;
ull2=getnum(i);
ull2*=b.getnum(j);
ull+=ull2;
rare+=ull/MAXINT;
ull=ull%MAXINT;
}
r.setnum(k,(int)ull);
}
while(rare>0)
{
r.setnum(k,rare%MAXINT);
rare=rare/MAXINT;
k++;
}
return r;
}
bigint operator*(const int& b)
{
long long ll;
int rare=0;
int i;
bigint t(0);
for(i=0;i<numbers.size();i++)
{
ll=numbers[i];
ll*=b;
ll+=rare;
t.setnum(i,(int)(ll%MAXINT));
rare=int(ll/MAXINT);
}
while(rare>0)
{
t.setnum(i,rare%MAXINT);
rare/=MAXINT;
}
return t;
}
bigint power(int n)
{
bigint p(*this);
static int masks[]={0x00,
0x01,
0x02,0x02,
0x04,0x04,0x04,0x04,
0x08,0x08,0x08,0x08,0x08,0x08,0x08,0x08};
int mask=0,caln=n;
if(n==0)
return 1;
if(caln&0x7fff0000)
{
mask=16;
caln>>=16;
}
caln&=0xffff;
if(caln&0xff00)
{
mask+=8;
caln/=0x100;
}
caln&=0xff;
if(caln&0xf0)
{
mask+=4;
caln/=0x10;
}
caln&=0xf;
mask=(masks[caln]<<mask);
mask>>=1;
while(mask)
{
p=p*p;
if(n&mask)
p=p*(*this);
mask>>=1;
}
return p;
}
int len() const
{
return numbers.size();
}
int getnum(int i) const
{
if(i<0||i>=numbers.size())
return 0;
return numbers[i];
}
bool setnum(int i,int value)
{
if(i<0||i>numbers.size())
return false;
if(i==numbers.size())
numbers.push_back(value);
else
numbers[i]=value;
return true;
}
string to_string()
{
int i=len();
string r,s0="000000000",s1;
i--;
r=std::to_string(numbers[i]);
while(i>0)
{
i--;
s1 = std::to_string(numbers[i]);
if(s1.length()<9)
s1=s0.substr(0,9-s1.length())+s1;
r+=s1;
}
return r;
}
private:
vector<int> numbers;
};
int main()
{
int N=10000;
cout<<"input N:";
cin>>N;
clock_t begin_time = clock();
vector<int> factors={2,3,5,7};
vector<int> counts;
int i,j,c,f,m;
for(i=0;i<factors.size();i++)
{
c=0;
f=factors[i];
m=N;
while(m>=f)
{
m/=f;
c+=m;
}
counts.push_back(c);
cout<<counts.size()<<":"<<f<<","<<c<<endl;
}
for(i=11;i<=N;i+=2)
{
bool bf=true;
for(j=0;j<factors.size();j++)
{
f=factors[j];
if(f*f>i) break;
if(i%f==0)
{
bf=false;
break;
}
}
if(bf)
{
c=0;
m=N;
f=i;
while(m>=f)
{
m/=f;
c+=m;
}
factors.push_back(i);
counts.push_back(c);
cout<<counts.size()<<":"<<f<<","<<c<<endl;
}
}
cout<<"depart:"<<(clock()-begin_time)<<endl;
m=counts[2];
counts[0]=counts[0]-m;
counts[2]=0;
i=counts.size();
bigint fact(1);
int curcount=0;
while(i-->0)
{
if(counts[i]==0)continue;
if(curcount!=counts[i])
{
bigint fs(factors[i]);
curcount=counts[i];
while(i-->0)
{
if(curcount!=counts[i])
break;
fs=fs*factors[i];
}
i++;
if(curcount!=1)
fs=fs.power(curcount);
fact=fact*fs;
}
}
cout<<"calc:"<<(clock()-begin_time)<<endl;
if(m>=10)
cout<<fact.to_string()<<"e+"<<m<<endl;
else
{
cout<<fact.to_string();
for(i=0;i<m;i++)
cout<<"0";
cout<<endl;
}
cout<<"disp:"<<(clock()-begin_time)<<endl;
getchar();
return 0;
}
扫一扫
425
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
