高精度计算是一种程序设计的算法。由于中央处理器的字长限制,如32位CPU中一个整数最大只能取值4,294,967,295(=2^32-1),因此在超范围数值计算中,往往要采用模拟手段。通常使用分离字符的方法来处理数字数组。 维基百科【高精度计算】
在上一辑中,给大家介绍了如何使用matlab自带工具箱实现高精度计算(详见:如何用matlab做高精度计算?【第一辑】)。本期给大家带来两款来自File Exchange源代码共享资源库的宝贝,它们都是出自大神John D'Errico之手。前者是专门用于处理超大值整数运算的 —— Variable Precision Integer Arithmetic,对应数据类型为vpi,后者是用于处理浮点数计算的 —— HPF (a big decimal class),对应数据类型为hpf。
接下来咱们就一起来揭开这两款工具箱的神秘面纱。
一、Variable Precision Integer Arithmetic
从名称不难看出,VPI工具箱是专门用于处理整数相关计算的,而大神John D'Errico开发VPI工具箱的目的也很简单,因为他不能用MATLAB符号处理工具箱,而他的这种款工具箱意义已经超越了其初衷。正所谓前人栽树,后人乘凉。有了VPI这款工具箱,原则任意的大整数计算在MATLAB都可以实现。接下来就一起来看看如何使用这款工具吧!
1.1 创建VPI类型数据
vpi类型数据创建有以下四种方法:
① num = vpi; 创建一个值为0的变量;
② num = vpi([]); 创建一个空的vpi变量
③ num = vpi(1415926); 将整数转换成vpi类型数字
④ num = vpi('12345678900987654321'); 将符号型整数转换成vpi类型数字
注:为了防止因MATLAB精度限制造成超大整数无法完整准确表达,建议对于超大整数使用符号型整数来处理。
1.2 类型数据类型转换
使用double或single可以将vpi型数据转换成MATLAB的double或single型,但此操作会造成数据精度损失,如下图所示:
1.3 VPI类型数据的运算
vpi类型的数据除支持四则运算、比较运算外,还支持指数、开方等复杂运算规则,因vpi型数据是整数型数据,开方时只会得到最近的整数结果。
1.4 VPI类型数据进制转换
vpi型数据还有特殊功能,支持数据进制转换,相关的函数有:
vpi2base:vpi型数据到任意进制转换;
vpi2bin:vpi型数据到二进制转换;
vpi2english:vpi型数据的英语读法;
base2vpi:任意进制到vpi型数据;
bin2vpi:二进制到vpi型数据;
1.5 vpi类所支持的方法
代码语言:javascript复制Methods for class vpi:
abs double gcd isprime lt ne rdivide sum vpi2english
ceil eq ge isunit max nthroot rem times
comparemagnitudes exp gt iszero min num2str round trailingdigit
conv factor isequal lcm minus order shiftdec uminus
cumprod factorial iseven le mod plus sign unique
cumsum find isfinite leadingdigit mpower power single uplus
digits fix isinf log mrdivide prod sort vpi
disp floor isnan log10 mtimes quotient sortrows vpi2base
display full isnumeric log2 nchoosek randint sqrt vpi2bin
除了上面介绍的内容外,VPI工具箱更多功能函数小伙伴自己去探究吧,因篇幅不宜过长,咱就介绍到此。总而言之,VPI工具箱非常强大,能满足在matlab中的绝大部分整数计算需求。
二、High Precision Floating
大神John D'Errico在分享了超大整数计算工具箱后又分享了任意精度浮点数计算工具箱 —— High Precision Floating (HPF)。与上一期介绍的如何用matlab做高精度计算?【第一辑】一样,作为浮点数计算,首先需要人为自定义计算所需要的精度,如果不设置,则会使用默认的精度进行计算。
1.1 hpf类型数据数据精度定义
hpf型数据精度设置涉及到两个函数 —— DefaultNumberOfDigits与
DefaultDecimalBase。前者用于设置浮点数位数,后者用于设置存储数字的字段长度。
DefaultNumberOfDigits使用方式如下:
DefaultNumberOfDigits(nDigits,sDigits),其中nDigits是实际显示浮点数位数,sDigits是保护位数,因为浮点数计算比纯整数计算误差更不易掌控,因此数据实际显示位数后再加sDigits位数据来确保前面数据的准确性。
代码语言:javascript复制使用方式:
① DefaultNumberOfDigits(nDigits)
② DefaultNumberOfDigits(nDigits,sDigits)
③ DefaultNumberOfDigits nDigits
④ DefaultNumberOfDigits nDigits sDigits
DefaultDecimalBase使用方式如下:
DefaultDecimalBase(DBase),其中,DBase只能取1,2,3,4,5或6,DBase取不同值,能表示最大浮点数位数各不相同,如下:
代码语言:javascript复制DBase = 1 --> 3.6e14 decimal digits
DBase = 2 --> 3.6e12 decimal digits
DBase = 3 --> 3.6e10 decimal digits
DBase = 4 --> 3.6e8 decimal digits
DBase = 5 --> 3.6e6 decimal digits
DBase = 6 --> 36000 decimal digits
DBase越大,计算速度越快,反正越慢。故默认DBase为6,对于计算不是特别大时,建议不要修改Dbase。
2.2 hpf类型数据的创建
hpf类型数据创建有以下五种方法:
① num = hpf; 创建一个值为0的变量;
② num = hpf([]); 创建一个空的hpf变量
③ num = hpf(3.141592653589793238462); 将整数转换成hpf类型数字
④ num = hpf('3.141592653589793238462'); 将符号型整数转换成hpf类型数字
⑤ num = hpf('3.141592653589793238462',num); 将符号型整数转换成hpf类型数字,并指定位数为num。
和vpa函数一样,如果需要创建高精度hpf型数据,建议不要采用方法③,因为会受matlab的浮点数有效位数的影响导致被截断。
2.3 hpf类所支持的方法
和vpi类型数据一样,hpf类型数据也是一种数据类,在命令窗口输入methods('hpf')即可查看hpf所支持的所有方法。
代码语言:javascript复制Methods for class hpf:
abs asinh csc erf int32 lt plus sign uint32
acos atan cscd exp int64 lu power sin uint64
acosd atand csch factorial int8 mantissa prod sind uint8
acosh atanh cubrt fix isfinite max rank single uminus
acot augmentdigits cumprod floor isinf min rat sinh uplus
acotd ceil cumsum fminsearch isnan minus rdivide sort vpi
acsc chol det fractionalpart isnumeric mod reciprocal sqrt
acscd cos diag full iszero mpower round sum
adjustdecimalbase cosd disp ge le mrdivide roundn tan
asec cosh display gt linspace mtimes rref tand
asecd cot double hpf log ne sec tanh
asin cotd eps int16 log10 nthroot secd times
asind coth eq int2str log2 num2str sech uint16
Static methods:
eye ones ten zeros
2.4 hpf型数据应用案例 (来源于作者提供的帮助文档)
使用Chudnovsky公式计算π:
代码语言:javascript复制DefaultNumberOfDigits 5000 10
DefaultDecimalBase 5
a1 = hpf(13591409); % small integers, so exact
a2 = hpf(545140134);
a3 = hpf(640320);
a4 = a3.*a3.*a3;
num = hpf(1);
den = sqrt(a4);
tol = eps(num);
piinv = hpf(0);
k = 0;
while 1
piterm = num.*(a1 a2.*k)./den;
piinv = piinv piterm;
if abs(piterm) < tol
break
end
k = k 1;
num = num.*prod(hpf(6*k (-5:0)));
den = den.*a4.*(3*k-2).*(3*k-1).*(3*k).*hpf(k).^3;
num = -num;
end
piest = reciprocal(12.*piinv)
Ref: en.wikipedia.org/wiki/Chudnovsky_algorithm
以上就是matlab做高精度计算第二辑的全部内容,更多地应用伙伴们自己去开发,本文仅是作为引子的作用。下一辑将会为大家介绍一款收费的高精度计算工具箱 —— Multiprecision Computing Toolbox,执行效率远超matlab自带的vpa工具箱。欲知后事如何,且看下回分解!
参考资料:
[1] www.mathworks.com/help/symbolic/vpa.html
[2] https://zh.m.wikipedia.org/zh-hans/高精度计算