针对 P1001 A+B Problem 的大炮使用手册
第一次在洛谷交题解。
大炮维护日志：修改了大炮的运行逻辑，特别优化了其在对蚊子实施打击时的使用体验。同时，同步更新了《大炮使用手册》，润色了一些句子，并增添了对大炮新功能的说明。
思路：
观察数据范围，发现答案竟然在 \([-2^{31},2^{31}-1]\) 范围内，为了防止溢出，我们考虑使用高精度算法，然而朴素的高精度算法时间复杂度为 \(O(n)\)，也许需要一定程度的优化。
我们知道，FFT（快速傅里叶变换）可以高效地计算多项式卷积，因此被用来优化乘法的计算。借用这种算法，我们尝试对加法进行优化。具体思路借用了这道题。
FFT 具有线性性质。这也是我们应用其优化运算的基础。具体地，对于两个序列 \(a\) 和 \(b\)，有 \(\operatorname{FFT}(a\times b)=\operatorname{FFT}(a)\times \operatorname{FFT}(b)\)，同理，有 \(\operatorname{FFT}(a+ b)=\operatorname{FFT}(a)+\operatorname{FFT}(b)\)，考虑利用这个性质，对 \(a\) 和 \(b\) 进行 FFT，然后将结果相加，最后再对相加后的结果逆 FFT 即得答案。
常数优化：
高精度运算的一种简单且有效的优化方式就是压位，然而 FFT 在计算过程中涉及浮点运算，这会带来严重的精度问题，尤其是在使用压位高精的时候。所以我们只压两位。对于本题，\(\dfrac{1}{2}\) 的常数已经足够了。
使用体验优化：
为了方便大炮的使用，我对其进行了封装，同时内附一份标准的 FFT 优化高精度乘法运算模板。
#include<cmath>
#include<iostream>
#include<cstring>
//********一些辅助函数与FFT实现
const double PI=4*atan(1);
template<typename T>
void Swap(T &a,T &b){
T c=a;
a=b;
b=c;
return;
}
template<typename T>
T Max(const T &a,const T &b){
return a<b?b:a;
}
typedef long long ll;
struct comp{
double real,imag;
comp operator+(const comp &x)const{
return {real+x.real,imag+x.imag};
}
comp operator-(const comp &x)const{
return {real-x.real,imag-x.imag};
}
comp operator*(const comp &x)const{
return {real*x.real-imag*x.imag,real*x.imag+x.real*imag};
}
comp operator/(const unsigned &x)const{
return {real/(double)x,imag/(double)x};
}
};
void FFT(comp *f,unsigned n,int rev){
for(unsigned i=1,j=n>>1,k;i<n-1;i++){//位逆序置换
if(i<j)
Swap(f[i],f[j]);
k=n>>1;
while(j>=k){
j-=k;
k>>=1;
}
j+=k;
}
for(unsigned l=2;l<=n;l<<=1){//蝶形运算
double arg=2*PI*rev/l;
comp wn={cos(arg),sin(arg)};
for(unsigned i=0;i<n;i+=l){
comp w={1,0};
for(unsigned j=0;j<(l>>1);j++){
comp f1=f[i+j];
comp f2=f[i+j+(l>>1)];
f[i+j]=f1+w*f2;
f[i+j+(l>>1)]=f1-w*f2;
w=w*wn;
}
}
}
if(!~rev)
for(unsigned i=0;i<n;i++)
f[i]=f[i]/n;
}
//********高精度整数类
#define BASE 100//因为FFT的精度问题严重，我们只压2位
template<const unsigned Size>
class bigint{
private:
unsigned len;
int num[Size];
void init(){
memset(num,0,sizeof(num));
len=1;
}
bool abs_greater_equal(const bigint &a)const{
if(len!=a.len)
return len>a.len;
for(int i=len;i;i--)
if(num[i]!=a.num[i])
return num[i]>a.num[i];
return 1;
}
public:
bigint(){
init();
}
void get_num(std::string s){
init();
int f=0;
unsigned slen=s.length();
if(s[0]=='-')
num[0]=f=1;
len=0;
unsigned temp=0,w=1;
for(int i=slen-1;i>=f;i--){
temp+=(s[i]^48)*w;
w=(w<<1)+(w<<3);
if(w==BASE||i==f){
num[++len]=(int)temp;
temp=0;
w=1;
}
}
if(temp||len==0)
num[++len]=temp;
}
bool operator<(const bigint &a)const{
if(num[0]&&!a.num[0])
return 1;
if(!num[0]&&a.num[0])
return 0;
if(num[0]){
if(len!=a.len)
return len>a.len;
for(int i=len;i;i--)
if(num[i]!=a.num[i])
return num[i]>a.num[i];
}
else{
if(len!=a.len)
return len<a.len;
for(int i=len;i;i--)
if(num[i]!=a.num[i])
return num[i]<a.num[i];
}
return 0;
}
bigint operator+(const bigint &a)const{
bigint res;
if(len==1&&num[1]==0){
res=a;
return res;
}
if(a.len==1&&a.num[1]==0){
res=*this;
return res;
}
if(num[0]==a.num[0]){
res.num[0]=num[0];
unsigned len_sum=1;
while(len_sum<len+a.len)
len_sum<<=1;
comp *fa=new comp[len_sum]();
comp *fb=new comp[len_sum]();
for(unsigned i=0;i<len;i++)
fa[i]={(double)num[i+1],0};
for(unsigned i=0;i<a.len;i++)
fb[i]={(double)a.num[i+1],0};
FFT(fa,len_sum,1);
FFT(fb,len_sum,1);
for(unsigned i=0;i<len_sum;i++)
fa[i]=fa[i]+fb[i];
FFT(fa,len_sum,-1);
res.len=Max(len,a.len);
ll temp=0;
for(unsigned i=0;i<res.len;i++){
ll val=(ll)round(fa[i].real)+temp;
res.num[i+1]=(int)(val%BASE);
temp=val/BASE;
}
if(temp)
res.num[++res.len]=temp;
while(res.len>1&&res.num[res.len]==0)
res.len--;
delete[] fa;
delete[] fb;
}
else{
if(abs_greater_equal(a)){
res.num[0]=num[0];
unsigned len_sum=1;
while(len_sum<len+a.len)
len_sum<<=1;
comp *fa=new comp[len_sum]();
comp *fb=new comp[len_sum]();
for(unsigned i=0;i<len;i++)
fa[i]={(double)num[i+1],0};
for(unsigned i=0;i<a.len;i++)
fb[i]={(double)a.num[i+1],0};
FFT(fa,len_sum,1);
FFT(fb,len_sum,1);
for(unsigned i=0;i<len_sum;i++)
fa[i]=fa[i]-fb[i];
FFT(fa,len_sum,-1);
res.len=Max(len,a.len);
ll temp=0;
for(unsigned i=0;i<res.len;i++){
ll val=(ll)round(fa[i].real)+temp;
if(val<0){
val+=BASE;
temp=-1;
}
else
temp=0;//借位
res.num[i+1]=(int)(val%BASE);
}
if(temp)
res.num[++res.len]=temp;
while(res.len>1&&res.num[res.len]==0)
res.len--;
delete[] fa;
delete[] fb;
}
else{
res.num[0]=a.num[0];
unsigned len_sum=1;
while(len_sum<len+a.len)
len_sum<<=1;
comp *fa=new comp[len_sum]();
comp *fb=new comp[len_sum]();
for(unsigned i=0;i<len;i++)
fa[i]={(double)num[i+1],0};
for(unsigned i=0;i<a.len;i++)
fb[i]={(double)a.num[i+1],0};
FFT(fa,len_sum,1);
FFT(fb,len_sum,1);
for(unsigned i=0;i<len_sum;i++)
fa[i]=fb[i]-fa[i];
FFT(fa,len_sum,-1);
res.len=Max(len,a.len);
ll temp=0;
for(unsigned i=0;i<res.len;i++){
ll val=(ll)round(fa[i].real)+temp;
if(val<0){
val+=BASE;
temp=-1;
}
else
temp=0;
res.num[i+1]=(int)(val%BASE);
}
if(temp)
res.num[++res.len]=temp;
while(res.len>1&&res.num[res.len]==0)
res.len--;
delete[] fa;
delete[] fb;
}
if(res.len==1&&res.num[1]==0)
res.num[0]=0;
}
return res;
}
bigint operator*(const bigint &a)const{
bigint res;
if((len==1&&num[1]==0)||(a.len==1&&a.num[1]==0))
return res;
res.num[0]=num[0]^a.num[0];
unsigned len_sum=1;
while(len_sum<len+a.len)
len_sum<<=1;
comp *fa=new comp[len_sum]();
comp *fb=new comp[len_sum]();
for(unsigned i=0;i<len;i++)
fa[i]={(double)num[i+1],0};
for(unsigned i=0;i<a.len;i++)
fb[i]={(double)a.num[i+1],0};
FFT(fa,len_sum,1);
FFT(fb,len_sum,1);
for(unsigned i=0;i<len_sum;i++)
fa[i]=fa[i]*fb[i];
FFT(fa,len_sum,-1);
res.len=len+a.len;
ll temp=0;
for(unsigned i=0;i<res.len;i++){
ll val=(ll)(fa[i].real+0.5)+temp;
res.num[i+1]=(int)(val%BASE);
temp=val/BASE;
}
if(temp)
res.num[++res.len]=temp;
while(res.len>1&&res.num[res.len]==0)
res.len--;
delete[] fa;
delete[] fb;
return res;
}
void read(){
init();
std::string s;
char ch=getchar();
while(ch<'0'||ch>'9'){
if(ch=='-')
s.push_back('-');
ch=getchar();
}
while(ch>='0'&&ch<='9'){
s.push_back(ch);
ch=getchar();
}
get_num(s);
}
void print(){
if(num[0])
putchar('-');
bool leading_zero=1;
for(int i=len;i;i--){
if(leading_zero)
printf("%d",num[i]);
else
printf("%02d",num[i]);
leading_zero=0;
}
putchar('\n');
return;
}
};
//********程序主体
const int N=1<<10,M=100;
int n,m;
bigint<114514> a,b,c;
int main(){
a.read();
b.read();
c=a+b;
c.print();
return 0;
}
写在结尾：
我们成功把 \(O(1)\) 的朴素算法优化到了 \(O(n)\)，接着进一步将 \(O(n)\) 的高精度算法优化到了 FFT 的 \(O(n \log n)\)，这体现了我们对 OI 知识的灵活迁移运用。希望大家能（可千万别）将这种思想运用到学习生活中。
浙公网安备 33010602011771号