【算法笔记】 压位高精度整除法

算法介绍

主要是模拟竖式计算的过程。可以扩展为 BASE 进制数的高精度整除除问题。不考虑负数问题。

记：
BASE为x（一般情形下为10）
被除数为 (A = a_p x^p + a_{p-1} x^{p-1} + ...... + a_1 x + a_0)
除数为 (B = b_q x^q + ...... +b_1 x + b_0)

当 (A<B) 时，我们有 (A/B=0)

此处，我们主要考虑 (A geq B)

这种情形下，必有 (p>q space or space a_p > b_q)

让我们观察竖式计算的过程。

不难发现，竖式计算原理如下：
要求 C 使得 (A = C imes B + R) 且 (R < B)
根据乘法分配律，可分解为 (A = (c_t x^t B + ... + c_1 x B + c_0 B) + R)
换句话说，我们只需要针对每个 (x^d B) 求出对应的 (c_d) 就行了
实际上，(x^d B) 就相当于 x 进制下的 B 移动 d 位。
而答案 C 可以表示为x进制下的 c1c2...ct

接下来考虑具体实现。

首先，考虑每一次计算时，要移动几位的问题。
这里我采用了一种比较简单的方式：移动到最高位重叠。
对于可能出现的移动后被除数最高位比除数最高位小的问题，我的做法是将此处的 (c_i) 记为 0，再将被除数的最高位乘以 x 后并入次高位里。
由于 (a_i < b_j) ，我们有 (b_j imes x > a_i imes x + a_{i-1})
也就是说，我们算出的 (c_k*b_j) 的 (c_k) 必然小于 x
仍然符合 x 进制数的特点。

接下考虑如何估计 (c_i)
在笔算时，我们有时不能马上看出答案。我们模拟竖式计算时也同理。
我采用的是利用最高位相除来估算的方法。
例如 (900 div 199) 估算百位时，直接看成 c = 9/1 ，再慢慢减小直到 199*c < 900

上面的这个做法在不压位(BASE较小或为10)的情形下需要减的次数不高，可视为常数。
但在压位下，需要减多少次就显得十分玄学。
对于这两个数 A 和 B，假定 p = q：
(A = a_p x^p + ...... + a_1 x + a_0)
(B = b_q x^q + ...... +b_1 x + b_0)
那么需要减的次数就是 (a_p div b_p - A div B)（整除）。
这个数的范围我不会求，有可能与 BASE 的规模相近。
我下文代码库中使用的是这种写法，但是 ZZY大师 给出了一种优化的思路，就是在求 c 的时候不是一个一个地减，而是从预估出的 c 的最高位开始减。这样就可以将规模缩小为与压位的位数相近，而不是与 BASE 相近。

代码库

数据生成器

#include <cstdio>
#include <cstdlib>
#include <ctime>
int main(){
    char str[]="0.in";
    srand(time(NULL));
    for(int i=1;i<=9;i++){
        str[0]=i+48;
        freopen(str,"w",stdout);
        putchar(i+48);
        for(int j=2;j<=1000;j++) putchar(rand()%10+48);
        putchar('
');
        putchar(10-i+48);
        for(int j=2;j<=233;j++) putchar(rand()%10+48);
        putchar('
');
        fclose(stdout);
    }
    return 0;
}

数据生成器 python 部分

root = "D:\Data_Testing\data\Bint Div"
for i in range(1,10) :
    f = open(root+"\%d.in"%(i),"r");
    a = f.readline(); b = f.readline();
    a1 = int(a); b1 = int(b);
    f.close();
    f = open(root+"\%d.out"%(i),"w+");
    f.write(str(a1//b1));
    f.close;

数据验证器 python

a = int(input()); 
b = int(input());
print(a//b);

压位高精度整除

#include <cstdio>
#include <cstring>
typedef long long ll;
#define REG register
#define rep(i,a,b) for(REG int i=a;i<=b;i++)
#define Rep(i,a,b) for(REG int i=a;i>=b;i--)
const int N=2000,BASE=100000;
inline char getc(){
    static char buf[1<<14],*p1=buf,*p2=buf;
    return p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<14,stdin),p1==p2)?EOF:*p1++;
}
struct Bint{
    ll a[N+2];
    Bint(){ memset(a,0,sizeof(a)); }
    int getSiz() const{
        Rep(i,N,0) if(a[i]) return i;
        return -1;
    }
    void maintain(){
        rep(i,0,N){
            if(a[i]<0) a[i+1]-=(-a[i]/BASE+1),a[i]+=(-a[i]/BASE+1)*BASE;
            else if(a[i]>=BASE) a[i+1]+=a[i]/BASE,a[i]%=BASE;
        }
    }
    void test()const{
        rep(i,0,15) printf("%lld ",a[i]);
        putchar('
');
    }
    bool operator<=(const Bint&B)const{
        Rep(i,N,0) if(a[i]!=B.a[i]) return a[i]<B.a[i];
        return 1;
    }
    bool operator<(const Bint&B)const{
        Rep(i,N,0) if(a[i]!=B.a[i]) return a[i]<B.a[i];
        return 0;
    }
    Bint operator<<(const int&x)const{
        Bint C;
        rep(i,0,N-x) C.a[i+x]=a[i];
        return C;
    }
    Bint operator*(const int&B)const{
        Bint C;
        rep(i,0,N) C.a[i]=a[i]*B;
        C.maintain();
        return C;
    }
    Bint operator/(const Bint&B)const{
        Bint C=(*this),temp,ans,QQ; REG int t1,t2=B.getSiz();
        while(B<=C){
            t1=C.getSiz(); temp=B<<(t1-t2);
            REG ll div=C.a[t1]/temp.a[t1];
            //printf("p:%d
",div);
            while(1){
                QQ=temp*div;
                int st=QQ.getSiz();
                if(st==t1+1) QQ.a[t1]+=QQ.a[st]*BASE,QQ.a[st]=0;
                if(!(C<=QQ)) break;
                div--;
            }
            //printf("%d
",div);
            rep(i,0,N) C.a[i]-=temp.a[i]*div;
            C.maintain();
            C.a[t1-1]+=C.a[t1]*BASE; C.a[t1]=0;
            //C.test(); putchar('
');
            ans.a[t1-t2]=div;
        }
        return ans;
    }
    ll getll()const{
        ll x=0;
        Rep(i,N,0) x=x*BASE+a[i];
        return x;
    }
};
inline Bint scan(){
    REG int c=0;
    Bint T; static char ch=0,s[N];
    memset(s,0,sizeof(s));
    while(ch<48) ch=getc();
    while(ch>=48) s[++c]=ch-48,ch=getc();
    #define BAS 5
    REG int x,pos=c,cnt=0;
    while(pos-BAS+1>0){
        x=0;
        rep(i,pos-BAS+1,pos) x=x*10+s[i];
        T.a[cnt++]=x; pos-=BAS;
    }
    if(pos>0){
        x=0;
        rep(i,1,pos) x=x*10+s[i];
        T.a[cnt++]=x;
    }
    return T;
}
inline void print(const Bint&x){
    REG bool f=0;
    Rep(i,N,0){
        //别忘了 !f 因为咱压了位
        if(!f&&x.a[i]){
            f=1; printf("%lld",x.a[i]);
        }else if(f) printf("%.5lld",x.a[i]);
    }
    if(!f) putchar('0');
}
int main(){
    //freopen("1.in","r",stdin);
    #define pc putchar('
')
    Bint A=scan(),B=scan(),C;
    //A.test(); pc;
    //B.test(); pc;
    C=A/B;
    print(C); pc;
    //C.test(); pc;
    //原先t1 t2写反了 出现了奇怪的问题
    //print(A); pc;
    //printf("%lld
",A.getll()/B.getll());
    return 0;
}
//通过观察发现 答案基本是对的，就是丢了几个0
//然后查出是输出的问题

END