zoukankan      html  css  js  c++  java
  • 从零开始学算法:高精度计算

    注:转载请注明:http://www.cnblogs.com/ECJTUACM-873284962/ 

    前言:由于计算机运算是有模运算,数据范围的表示有一定限制,如整型int(C++中int 与long相同)表达范围是(-2^31~2^31-1),unsigned long(无符号整数)是(0~2^32-1),都约为几十亿.如果采用实数型,则能保存最大的double只能提供15~16位的有效数字,即只能精确表达数百万亿的数.因此,在计算位数超过十几位的数时,不能采用现有类型,只能自己编程计算.
    高精度计算通用方法:高精度计算时一般用一个数组来存储一个数,数组的一个元素对应于数的一位(当然,在以后的学习中为了加快计算速度,也可用数组的一个元素表示数的多位数字,暂时不讲),表示时,由于数计算时可能要进位,因此为了方便,将数由低位到高位依次存在数组下标对应由低到高位置上,另外,我们申请数组大小时,一般考虑了最大的情况,在很多情况下,表示有富余,即高位有很多0,可能造成无效的运算和判断,因此,我们一般将数组的第0个下标对应位置来存储该数的位数.如数:3485(三千四百八十五),表达在数组a[10]上情况是:
    下标  0    1    2    3     4    5    6    7    8    9  
    内容  4    5    8    4     3    0    0    0    0    0
    说明:位数   个位  十位  百位 千位
    具体在计算加减乘除时方法就是小学时采用的列竖式方法.
    注:高精度计算时一般用正数,对于负数,通过处理符号位的修正.
    一.高精度数的存储
    1.如对数采用的字符串输入

     1 #include <iostream>
     2 #include <cstring>
     3 using namespace std;
     4 const int N=100;//最多100位
     5 int main()
     6 {
     7 int a[N+1],i;
     8 string s1;
     9 cin>>s1;//数s1
    10 memset(a,0,sizeof(a)); //数组清0
    11 a[0]=s1.length(); //位数
    12 for(i=1;i<=a[0];i++) a[i]=s1[a[0]-i]-'0';//将字符转为数字并倒序存储.
    13 return 0;
    14 }

    2.直接读入

     1 #include <iostream>
     2 using namespace std;
     3 const int N=100;//最多100位
     4 int main()
     5 {
     6 int a[N+1],i,s,key;
     7 cin>>key;//数key
     8 memset(a,0,sizeof(a)); //数组清0
     9 i=0;//第0位
    10 while(key)  //当key大于0
    11 {
    12   a[++i]=key%10;//取第i位的数
    13   key=key/10;
    14 }
    15 a[0]=i; //共i位数
    16 return 0;
    17 }

    3.直接初始化(用a[]存储)
    初始化为0: memset(a,0,sizeof(a));
    初始化为1: memset(a,0,sizeof(a));a[0]=1;a[1]=1;

    以下程序都只写函数,不写完整程序,所有高精度数存储都满足上述约定。
    二.高精度数比较

    1 int compare(int a[],int b[])   //比较a和b的大小关系,若a>b则为1,a<b则为-1,a=b则为0
    2 {int i;
    3 if (a[0]>b[0]) return 1;//a的位数大于b则a比b大
    4 if (a[0]<b[0]) return -1;//a的位数小于b则a比b小
    5 for(i=a[0];i>0;i--)  //从高位到低位比较
    6      {if (a[i]>b[i]) return 1;
    7       if (a[i]<b[i]) return -1;}
    8 return 0;//各位都相等则两数相等。
    9 }

    三、高精度加法

     1 int plus(int a[],int b[]) //计算a=a+b
     2 {int i,k;
     3 k=a[0]>b[0]?a[0]:b[0]; //k是a和b中位数最大的一个的位数
     4 for(i=1;i<=k;i++)
     5     {a[i+1]+=(a[i]+b[i])/10;  //若有进位,则先进位
     6     a[i]=(a[i]+b[i])%10;}  //计算当前位数字,注意:这条语句与上一条不能交换。
     7 if(a[k+1]>0) a[0]=k+1;  //修正新的a的位数(a+b最多只能的一个进位)
     8                else a[0]=k;
     9 return 0;
    10 }

    四、高精度减法

     1 int gminus(int a[],int b[]);//计算a=a-b,返加符号位0:正数 1:负数
     2 { int flag,i
     3   flag=compare(a,b); //调用比较函数判断大小
     4 if (falg==0)//相等
     5   {memset(a,0,sizeof(a));return 0;}  //若a=b,则a=0,也可在return前加一句a[0]=1,表示是 1位数0
     6 if(flag==1) //大于  
     7   {  for(i=1;i<=a[0];i++)
     8       {  if(a[i]<b[i]){ a[i+1]--;a[i]+=10;} //若不够减则向上借一位
     9         a[i]=a[i]-b[i];}
    10      while(a[a[0]]==0) a[0]--; //修正a的位数
    11     return 0;}
    12 if (flag==-1)//小于  则用a=b-a,返回-1
    13     { for(i=1;i<=b[0];i++)       {  if(b[i]<a[i]){ b[i+1]--;b[i]+=10;} //若不够减则向上借一位
    14         a[i]=b[i]-a[i];}
    15       a[0]=b[0];
    16      while(a[a[0]]==0) a[0]--; //修正a的位数
    17     return -1;}
    18 }

    五、高精度乘法1(高精度乘单精度数,单精度数是指通常的整型数)

    1 int multi1(int a[],long  key) //a=a*key,key是单精度数  
    2 {int i,k;
    3 if (key==0){memset(a,0,sizeof(a));a[0]=1;return 0;} //单独处理key=0
    4 for(i=1;i<=a[0];i++)a[i]=a[i]*key;//先每位乘起来
    5 for(i=1;i<=a[0];i++){a[i+1]+=a[i]/10;a[i]%=10;} //进位
    6 //注意上一语句退出时i=a[0]+1
    7 while(a[i]>0) {a[i+1]=a[i]/10;a[i]=a[i]%10;i++;a[0]++];}  //继续处理超过原a[0]位数的进位,修正a的位数
    8 return 0;
    9 }

    六、 高精度除以低精度;
    算法:按照从高位到低位的顺序,逐位相除。在除到第j位时,该位在接受了来自第j+1位的余数后与除数相除,如果最高位为零,则商的长度减一。源程序如下:

     1 #include  <stdio.h>
     2 #define   N  500
     3 main()
     4 {
     5   int  a[N] = {0}, c[N] = {0};
     6   int  i, k, d, b;
     7   char  a1[N];  
     8   printf("Input 除数:");
     9   scanf("%d", &b);
    10   printf("Input 被除数:");
    11   scanf("%s", a1);
    12   k = strlen(a1);
    13   for(i = 0; i < k; i++)  a[i] = a1[k - i - 1] - '0';
    14   d = 0;
    15   for(i = k - 1; i >= 0 ; i--)
    16   {
    17      d = d * 10 + a[i];
    18      c[i] = d / b;
    19      d = d % b;      
    20   }   
    21   while(c[k - 1] == 0 && k > 1)  k--;  
    22   printf("商=");
    23   for(i = k - 1; i >= 0; i--)  printf("%d", c[i]);
    24   printf("
    余数=%d", d);   
    25 }     

    七、高精度乘以高精度(要求用尽可能少的存储单元);
    算法:用数组保存两个高精度数,然后逐位相乘,注意考虑进位和总位数。源程序如下:

     1 #include  <stdio.h>
     2 main()
     3 {
     4   int  a[240] = {0}, b[240] = {0}, c[480] = {0};
     5   int  i, j, ka, kb, k;
     6   char  a1[240], b1[240];
     7   gets(a1);   
     8   ka = strlen(a1);
     9   gets(b1);   
    10   kb = strlen(b1);
    11   k = ka + kb;
    12   for(i = 0; i < ka; i++)  a[i] = a1[ka-i-1] - '0';
    13   for(i = 0; i < kb; i++)  b[i] = b1[kb-i-1] - '0';
    14   for(i = 0; i < ka; i++)
    15     for(j = 0; j < kb; j++)
    16     {
    17       c[i + j] = c[i + j] + a[i] * b[j];
    18       c[i + j +1] = c[i + j +1] + c[i + j]/10;
    19       c[i + j] = c[i + j] % 10;
    20     }
    21   if(!c[k])  k--;
    22   for(i = k-1; i >= 0; i--)  printf("%d", c[i]);        
    23 }     

    八、高精度除以高精度(要求用尽可能少的存储单元);
    算法:用计算机模拟手算除法,把除法试商转化为连减。

      1 #include  <stdio.h>
      2 #define   N  500
      3 int  bj(int a[], int b[], int k1, int k2)   /*比较大小函数*/
      4 {
      5    int i, t, flag;       /*flag作标志位*/
      6    if(k1 < k2)  
      7      flag = 0;           /*被除数小于除数返回0*/
      8    else if(k1 > k2)  
      9           flag = 1;      /*被除数大于除数返回1*/
     10         else
     11           {              /*被除数和除数位数相等则逐位进行比较*/
     12             i = k1;
     13             t = 0;
     14             while(t == 0 && i > 0)
     15             {
     16               if(a[i] > b[i]) {t = 1; flag = 1;}
     17               else if(a[i] == b[i])  i--;
     18               else  {t = 1; flag = 0;}        
     19             }
     20             if(i == 0 && t == 0)  flag = 2;     /*被除数等于除数返回2*/
     21           }
     22   return flag;           
     23 }
     24 int  jf(int a[], int b[], int k1, int k2)       /*减法运算*/
     25 {
     26   int  i, k, d[N];
     27   for(i = 0; i < k2; i++)  d[i] = b[i];        /*把除数赋给数组d*/
     28   for(i = k2; i < N; i++)  d[i] = 0;          /*d数组无数据的高位置0*/
     29   k = k1 - k2 - 1;                            /*计算减法起始位置*/
     30   if(k < 0)  k = 0;
     31   if(k > 0)
     32   {
     33     for(i = k2 - 1; i >= 0; i--)  d[i + k] = d[i];  /*移动减数位数与被减数对齐*/
     34     for(i = 0; i < k; i++)  d[i] = 0;            /*移动后的其余位置0*/
     35   }  
     36   for(i = 0; i < k1; i++)
     37   {
     38     if(a[i] >= d[i])  a[i] -= d[i];
     39     else
     40     {
     41       a[i + 1] = a[i + 1] - 1;
     42       a[i] = 10 + a[i] - d[i];  
     43     }     
     44   }   
     45   return k;
     46 }
     47 main()
     48 {
     49   int  a[N] = {0}, b[N] = {0}, c[N] = {0}, d[N] = {0};
     50   int  i, ka, kb, m, t, t1, t2, k, x, kd, kk;
     51   char  a1[N], b1[N];  
     52   printf("Input 被除数:");
     53   scanf("%s", a1);
     54   ka = strlen(a1);
     55   for(i = 0; i < ka; i++)  a[i] = a1[ka - i -1] - '0';
     56   printf("Input 除数:");
     57   scanf("%s", b1);
     58   kb = strlen(b1);
     59   for(i = 0; i < kb; i++)  b[i] = b1[kb - i -1] - '0';
     60   kd = ka;    /*保存被除数位数  */
     61   t2 = bj(a, b, ka, kb);
     62   m = 0;
     63   do
     64   {
     65     while(a[ka - 1] == 0)  ka--;
     66     t = bj(a, b, ka, kb);   
     67     if(t >= 1)
     68     {
     69       k = jf(a, b, ka, kb);
     70       c[k]++;      
     71       if(k > m)  m = k;
     72       t1 = 0;
     73       for(i = k; i <= m; i++)
     74       {
     75         x = c[i] + t1;
     76         c[i] = x % 10;
     77         t1 = x / 10;     
     78       }
     79       if(t1 > 0)  {m++; c[m] = t1;  }     
     80     }   
     81   }while(t == 1);
     82   if(t2 == 0)  
     83   {
     84     printf("商=0");  
     85     printf("
    余数=");
     86     for(i = kd - 1; i >= 0; i--)  printf("%d", a[i]);
     87     exit(1);  
     88   }
     89   if(t2 == 2)
     90   {
     91     printf("商 = 1");  
     92     printf("
    余数 = 0");
     93     exit(1);  
     94   }
     95   kk = kd;
     96   while(!c[kd - 1])  kd--;
     97   printf("商 = ");
     98   for(i = kd - 1; i >= 0; i--)  printf("%d", c[i]);
     99   while(!a[kk])  kk--;
    100   printf("
    余数 = ");
    101   if(kk < 0)  
    102   {
    103     printf("0");  
    104     exit(1);
    105   }
    106   for(i = kk; i >= 0; i--)  printf("%d", a[i]);
    107 } 

    下面给出一些案例:

    问题1. N!,要求精确到P位(0〈P〈1000〉。
    算法:结果用数组a保存,开始时a[0]=1,依次乘以数组中各位,注意进位和数组长度的变化。源程序如下:

     1 #include   <stdio.h>
     2 #define    M   1000
     3 main()
     4 {
     5   int a[M], i, n, j, flag = 1;
     6   printf("n=");
     7   scanf("%d",&n);
     8   printf("n!=");
     9   a[0] = 1;
    10   for(i = 1; i < M; i++) a[i] = 0;
    11    for(j = 2; j <= n; j++)
    12    {
    13      for(i = 0; i < flag; i++) a[i] *= j;
    14      for(i = 0; i < flag; i++)
    15        if(a[i] >= 10)
    16        {
    17          a[i+1] += a[i]/10;
    18          a[i] = a[i] % 10;
    19          if(i == flag-1)  flag++;
    20        }
    21     }
    22   for(j = flag - 1; j >= 0; j--)
    23     printf("%d", a[j]);
    24 } 

    问题2. 麦森数
    【问题描述】形如2P-1的素数称为麦森数,这时P一定也是个素数。但反过来不一定,即如果P是个素数,2P-1不一定也是素数。到1998年底,人们已找到了37个麦森数。最大的一个是P=3021377,它有909526位。麦森数有许多重要应用,它与完全数密切相关。
    任务:从文件中输入P(1000<P<3100000),计算2P-1的位数和最后500位数字(用十进制高精度数表示)
    【输入格式】
    文件中只包含一个整数P(1000<P<3100000)
    【输出格式】
    第一行:十进制高精度数2P-1的位数。
    第2-11行:十进制高精度数2P-1的最后500位数字。(每行输出50位,共输出10行,不足500位时高位补0)
    不必验证2P-1与P是否为素数。
    【输入样例】
    1279
    【输出样例】
    386
    00000000000000000000000000000000000000000000000000
    00000000000000000000000000000000000000000000000000
    00000000000000104079321946643990819252403273640855
    38615262247266704805319112350403608059673360298012
    23944173232418484242161395428100779138356624832346
    49081399066056773207629241295093892203457731833496
    61583550472959420547689811211693677147548478866962
    50138443826029173234888531116082853841658502825560
    46662248318909188018470682222031405210266984354887
    32958028878050869736186900714720710555703168729087
    算法:2的幂可以转化成左移运算,为了提高运算速度,可每次左移10位,即每次乘210。对于个位单独考虑,每次左移一位。源程序如下:

     1 #include <stdio.h>
     2 #include <math.h>
     3 #define  MAX  100000
     4 main()
     5 {
     6    int p;
     7    int i, j;
     8    scanf("%d", &p);
     9    printf("%d
    ", (int)(p * log10(2.0)) + 1);
    10    long  store[110] = {0};
    11    store[0] = 1;
    12    int left = p % 10;
    13    p /= 10;
    14     for(i = 1; i <= p; i++)
    15     {
    16       for(j = 0; j <= 100; j++)
    17         store[j] <<= 10;
    18       for(j = 0; j <= 100; j++)
    19       {
    20         if(store[j] >= MAX)
    21         {
    22           store[j + 1] += store[j] / MAX;
    23           store[j] %= MAX;
    24         }
    25       }
    26     }
    27     for(i = 1; i <= left; i++)
    28     {
    29       for(j = 0; j <= 100; j++)
    30         store[j] <<= 1;
    31       for(j = 0; j <= 100; j++)
    32       {
    33         if(store[j] >= MAX)
    34         {
    35           store[j + 1] += store[j] / MAX;
    36           store[j] %= MAX;
    37         }
    38       }
    39     }
    40     store[0] -= 1;
    41     for(i = 1; i < 100; i++)
    42     {
    43       if(store[i - 1] < 0)
    44       {
    45          store[i] -= 1;
    46          store[i - 1] += MAX;
    47       }
    48       else
    49         break;
    50     }
    51     for(i = 99; i >= 0; i--)
    52     {
    53       printf("%05d", store[i]);
    54       if((100 - i) % 10 == 0)
    55           printf("
    ");
    56     }
    57 } 

    问题3. 有一个正整数N(N可能达到120位),它是由若干个不大于65535的正整数相乘而得到的。请把这个数分解成素数因子(质因子)的乘积。
    输入:输入文件只有一行为N的值。
    输出:(1)素数因子由小到大分行输出;
    (2)每一行输出一个素数因子和该素数因子的个数,用一个空格分开;
    (3)如果正整数N的分解中有一个以上的大于65535的素数,请按照(1)、(2)的要求输出分解中的小于65535的素数后,在下一行输出
    “DATA  ERROR!”。
    算法:先将2到65535之间的所有素数保存在数组中,用这个数去除数组中的每一个数,得到一个质因数就打印出来。源程序如下:

     1 #include  <stdio.h>
     2 #include  <math.h>
     3 int length, temp[120];
     4 int sushu(int a[])
     5 {
     6   int i, j, k = 0, m;
     7   for(i = 2; i <= 65537; i++)
     8   {
     9     m = sqrt(i);
    10     for(j = 2; j <= m; j++)
    11       if(i % j == 0)  break;
    12     if(j > m)
    13     {
    14       a[k] = i;
    15       k++;   
    16     }         
    17   }
    18  return k;     
    19 }
    20 int divide(int a[], int k)
    21 {
    22   int i, d = 0;
    23   for(i = length - 1; i >= 0; i--)
    24   {
    25      d = d * 10 + a[i];
    26      temp[i] = d / k;
    27      d = d % k;      
    28   }   
    29   if(!d)
    30     {
    31       while(temp[length - 1] == 0 && length > 1)  length--;
    32       for(i = 0; i < length; i++)
    33       {
    34         a[i] = temp[i];
    35         temp[i] = 0;      
    36       }
    37       for(i = length; i < 120; i++) a[i] = 0;   
    38     }
    39     else  
    40       for(i = 0; i < length; i++)  temp[i] = 0;
    41   return d;     
    42 }
    43 main()
    44 {
    45   int i, k, s, d;       /*s计数器; d余数*/
    46   int a[6600], b[120] = {0}, c[120] = {0};
    47   char b1[120];
    48   gets(b1);
    49  length = strlen(b1);
    50   for(i = 0; i < length; i++)  b[i] = b1[length - i - 1] - '0';
    51   k = sushu(a);
    52   for(i = 0; i < k; i++)
    53   {
    54     s = 0;
    55     d = divide(b, a[i]);
    56     while(!d)
    57     {
    58       s++;
    59       d = divide(b, a[i]);         
    60     }
    61     if(i == k - 1)  
    62  
    63     {  
    64       printf("Data Error!");
    65       break;
    66     } 
  • 相关阅读:
    python3.7 打包(.exe)神器——pyinstaller 安装及用法
    python3.7下运行pyspider报错的问题及解决方案
    python3一键排版证件照(一寸照、二寸照),附源代码
    傻瓜式下载“喜马拉雅”音频文件
    windows 7 32位环境下安装Redis、安装桌面管理工具redis-desktop-manager
    python3爬虫之验证码的识别——selenium自动识别验证码并点击提交,附源代码
    python3爬虫之验证码的识别——第三方平台超级鹰
    python3爬虫之验证码的识别——图形验证码
    python3爬虫之图形验证码的识别——环境安装
    scrapy爬虫笔记(入门级案例)
  • 原文地址:https://www.cnblogs.com/ECJTUACM-873284962/p/6509429.html
Copyright © 2011-2022 走看看