高精度运算及其应用.doc

Jsoi2010春季函授讲义（B2） 13/13高精度运算及其应用一、引言利用计算机进行数值运算，经常会遇到数值太大，超出Longint、int64等系统标准数据类型的有效范围，如计算mn，而m、n100；有时又会遇到对运算的精度要求特别高的情况，如计算圆周率，要求精确到小数点后100位，此时real、double等数据类型也无能为力。这些情况下，我们都要用“高精度运算”来解决。一般我们将小数点后几百位或者更多，当然也可能是几千亿几百亿的大数字统称为高精度数。高精度运算首先要解决存储问题。一般都是定义一个一维数组来存储一个高精度数，用每一个数组元素存储该数的每一位或某几位。高精度数的读入可以采用两种方法，一是采用字符串（String,AnsiString）方式一起读入，再逐位处理成数字存储在数组中；另一种方法是一位一位读入并存储到数组中。在实际使用时，请大家注意比较各自的优、缺点。高精度运算一般都是采用模拟的方法解决。输出时一定要注意格式和精度。二、高精度运算1、编程实现高精度加法问题描述 输入两个正整数（最多250位），输出它们的和。 比如输入：99999999999999999999999999999999999999999999999999999912345678999999999999999999999999 输出：add=1000000000000000000000012345678999999999999999999999998问题分析只要模拟“加法运算”的过程，从低位（对齐）开始逐位相加，最后再统一处理进位即可。参考程序Program ex1(input,output);const max=250;var s1,s2:string; a,b,c:array1.max of byte; l1,l2,l,i:integer;begin writeln(input two large integer:); readln(s1); readln(s2); 用字符串方式读入两个高精度数 l1:=length(s1); l2:=length(s2); for i:=1 to max do begin ai:=0;bi:=0;ci:=0;end; 注意一定要初始化 for i:=1 to l1 do ai:=ord(s1l1+1-i)-48; for i:=1 to l2 do bi:=ord(s2l2+1-i)-48; 以上是把两个高精度数逐位处理并转存到a、b两个数组中 if l1l2 then l:=l1 else l:=l2; for i:=1 to l do ci:=ai+bi; 对应位相加 for i:=1 to l do 从低位到高位，统一处理进位 if ci=10 then begin ci:=ci-10; ci+1:=ci+1+1; end; if cl+10 then l:=l+1; write(add=); 输出 for i:=l downto 1 do write(ci); readln;end.思考和练习1、 如果要一边加一边进位，程序怎么修改？你觉得好不好？2、 如果输入的数再大一点，比如1000位，还好用String类型读入吗？程序怎么修改？3、 请你编写一个高精度减法的程序，注意结果的正负。2、高精度乘法例2、编程求n!的值，n=10 then begin aj+1:=aj+1+aj div 10; 和高精度加法有何区别？ aj:=aj mod 10; end; while ah+10 do 最高位的进位处理，和高精度加法有何区别？ begin h:=h+1; ah+1:=ah div 10; ah:=ah mod 10; end; if ah+10 then h:=h+1; end; write(n,!=); 输出for i:=h downto 1 do write(ai); readln;end.运行示例输入：input n:500输出：500!=1220136825991110068701238785423046926253574342803192842192413588385845373153881997605496447502203281863013616477148203584163378722078177200480785205159329285477907571939330603772960859086270429174547882424912726344305670173270769461062802310452644218878789465754777149863494367781037644274033827365397471386477878495438489595537537990423241061271326984327745715546309977202781014561081188373709531016356324432987029563896628911658974769572087926928871281780070265174507768410719624390394322536422605234945850129918571501248706961568141625359056693423813008856249246891564126775654481886506593847951775360894005745238940335798476363944905313062323749066445048824665075946735862074637925184200459369692981022263971952597190945217823331756934581508552332820762820023402626907898342451712006207714640979456116127629145951237229913340169552363850942885592018727433795173014586357570828355780158735432768888680120399882384702151467605445407663535984174430480128938313896881639487469658817504506926365338175055478128640000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000思考和练习1、 如果要问你n!的结果有多少位，要不要做高精度运算，怎么做最好？请编程完成。2、 如果要问你n!的结果后面有多少个连续的0，要不要做高精度运算，怎么做最好？请编程完成。例3、输入两个高精度正整数m,n（m,n都在200位以内），输出它们的乘积。问题分析本题为一个高精度数乘以另一个高精度数。算法思想仍然是模拟乘法的过程，用一个数的每一位（从低位开始）逐位与另一个数的每一位相乘，同时处理进位。程序如下：参考程序Program ex3(input,output);const max=200;var s1,s2:string; a,b:array1.max of integer; c:array1.max*max of integer; l1,l2,w,jw,h,i,j,f:longint;begin assign(input,t3.in); 文件输入输出 assign(output,t3.out); reset(input); rewrite(output); readln(s1); 读入及预处理 readln(s2); l1:=length(s1); l2:=length(s2); for i:=1 to max do begin ai:=0;bi:=0;end; for i:=1 to max*max do ci:=0; for i:=1 to l1 do ai:=ord(s1l1+1-i)-48; for i:=1 to l2 do bi:=ord(s2l2+1-i)-48; jw:=0; 逐位乘 for i:=1 to l1 do for j:=1 to l2 do begin f:=ai*bj; jw:=f div 10; f:=f mod 10; w:=i+j-1; cw:=cw+f; cw+1:=cw+1+jw+cw div 10; cw:=cw mod 10; end; h:=l1+l2; while ch=0 do h:=h-1; 判断最高位 f:=0; 输出 for i:=h downto 1 do begin write(ci); f:=f+1; if f mod 30=0 then writeln; 30位换一行 end; close(input); close(output)end.运行示例输入：234234535436546546098948569895349583947593532459999999999999999999999999999999999999999999999999999999999900000000354348594543988989898989898988999898989898898776655432123456789765432456789输出：830006784256044531573047212363730755471182948895939218737538376740503405551785126396123228045386811816943500803653644720410101010101011000101010101101223344567876543210234567543211000000003、高精度除法例4、计算n/m的值，设n，m是integer范围内的正整数，要求精确到小数点后500位。问题分析依然采用模拟的方法。由数学知识可知： 除法运算中被除数、除数和商、余数的关系为： 新的被除数 = 10 * 余数；商 = trunc(被除数除数)；余数 = 被除数mod 除数。程序如下：参考程序Program ex4(input,output);const e=500;var a,d,x:array0.e of integer; n,m,t:integer;begin write(input n,m:); readln(n,m); write(n,/,m,=); a0:=n; d0:=n div m; x0:= n mod m; write(d0,.); 初始化处理，得到整数部分的商 for t:=1 to e do begin if xt-1=0 then exit; at:=xt-1 * 10; 将余数扩大10倍 dt:=at div m; write(d t); 计算商 xt:=at mod m; 生成新余数 end; readln;end.运行示例输入:input n，m: 355 113 输出:355/113=3.14159292035398230088495575221238938053097345132743362831858407079646017699115044247787610619469026548672566371681415929203539823008849557522123893805309734513274336283185840707964601769911504424778761061946902654867256637168141592920353982300884955752212389380530973451327433628318584070796460176991150442477876106194690265486725663716814159292035398230088495575221238938053097345132743362831858407079646017699115044247787610619469026548672566371681415929203539823008849557522123893805309734513274336练习请修改程序，要求达到如下输出要求：1、 如果除得尽，则有多少位就输出多少位，如输入：1 2，则输出0.5；2、 如果输出的是循环小数，则要用小括号将循环部分括起来输出，如输入：3 7，则输出：0.(428571)；思考其实，如果是一个真正的高精度数除以另一个高精度数（或单精度数），方法虽然比较明确，但是编程工作还是比较复杂的，留给有兴趣的同学自己去思考完成。一般程序中，我们都应该尽量避免做高精度除法。例5、麦森数【问题描述】形如2P-1的素数称为麦森数，这时P一定也是个素数。但反过来不一定，即如果P是个素数，2P-1不一定也是素数。到1998年底，人们已找到了37个麦森数。最大的一个是P=3021377，它有909526位。麦森数有许多重要应用，它与完全数密切相关。 任务：输入P（P1000），计算2P-1的位数和最后500位数字（用十进制高精度数表示）输入：只包含一个整数P（P1000）输出：第一行：十进制高精度数2P-1的位数。 第2-11行：十进制高精度数2P-1的最后500位数字。（每行输出50位，共输出10行，不足500位时高位补0）不必验证2P-1与P是否为素数。 样例输入 1279样例输出 38600000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000010407932194664399081925240327364085538615262247266704805319112350403608059673360298012239441732324184842421613954281007791383566248323464908139906605677320762924129509389220345773183349661583550472959420547689811211693677147548478866962501384438260291732348885311160828538416585028255604666224831890918801847068222203140521026698435488732958028878050869736186900714720710555703168729087问题分析第一问是很简单的，只需要求一个对数而已，数学原理：十进制正整数n的位数为int(log10(n)+1。所以2P-1的位数int(p*log102)+1 。 第二问的关键是高精度乘法和指数幂的运算，而且由于题目要求最后500位数字，所以在计算乘法的时候我们只要求计算乘数的低500位就好了。指数幂的运算不能硬乘，而要采用分治算法，否则就超时了。分治递归算法求指数幂是非常经典的，其数学原理是an = a(n div 2)*a(n div 2)*f(a)，其中f(a) = 1(n MOD 2=0)或f(a) = a(n MOD 2=1)。参考程序Program ex5(input,output);Var a,b:array 0.1001 of Longint; n,i,j:Longint; Procedure quick(k:Longint); Var i,c,j,l:Longint; Begin If k=0 Then Exit; quick(k div 2); If k mod 20 Then c:=2 Else c:=1; For l:=1 to 500 Do For i:=1 to 500 Do Inc(bi+l-1,ai*al*c); For i:=1 to 500 DoBegin bi+1:=bi+1+bi div 10; bi:=bi mod 10; End; a:=b; Fillchar(b,sizeof(b),0); End; Begin Readln(n); Writeln(trunc(ln(2)/ln(10)*n)+1); a1:=1; quick(n); Dec(a1); For i:=500 downto 2 Do Begin Write(ai); If i mod 50=1 Then Writeln; End; Writeln(a1); End.例6、组合数的末十位问题描述求出C(n,r)的最后十位，其中00 do t:=t+s s:=s div p end u:=send方法之三是利用公式C（n，r）= C（n-1，r）+ C（n-1，r-1），将计算C（n，r）的过程化为加法来做，该方法即为求杨辉三角中的第n行中的第r列上的数，行列从0开始。以下是一个参考程序，请大家研究，它用的是哪种方法？参考程序program ex6(input,output);const max=20;var a:array0.30000 of integer; t,i,j,n,r,s:longint;b:array1.max of byte;procedure divide(n:integer);var m:integer;begin m:=2; while m=trunc(sqrt(n) do begin while (n mod m=0) and (n0) do begin n:=n div m; am:=am-1; end; m:=m+1; end; an:=an-1;end;procedure times(n:integer);var m:integer;begin m:=2; while m=trunc(sqrt(n) do beginwhile (n mod m=0) and (n0) do begin n:=n div m; am:=am+1; end; m:=m+1; end; an:=an+1;end;procedure add(d:integer);var dd,i,t:longint;begin i:=1;j:=0; while i=max do begin t:=bi*d+j; bi:=t mod 10; j:=t div 10; i:=i+1; end; i:=1;j:=0; while i=max do begin t:=bi+j; bi:=t mod 10; j:=t div 10; i:=i+1;end;end;begin main write(Input n,r:); readln(n,r); fillchar(a,sizeof(a),0); i:=n-r+1;j:=r; repeat if i1 then begin divide(j);j:=j-1;end; until (in) and (j0 do begin add(i); ai:=ai-1; end; end; for j:=10 downto 1 do write(bj);writelnend.部分测试数据Input n,r:10 40000000210Input n,r:1067 6224592928000Input n,r:29999 199995854060848Input n,r:29999 149979718772800Input n,r:29999 273814531330240三、作业（只要交源程序)1、蜜蜂路线（bee.pas,bee.in,bee.out,时限1秒）问题描述一只蜜蜂在下图所示的数字蜂房上爬动,已知它只能从标号小的蜂房爬到标号大的相邻蜂房,现在问你：蜜蜂从蜂房M开始爬到蜂房N，MN=10000，有多少种爬行路线