b********r 发帖数: 1080 | 1 【 以下文字转载自 Computation 讨论区 】
发信人: bankbuster (恭喜发财), 信区: Computation
标 题: 怎么提高C++计算精度? C++ vs Matlab
发信站: BBS 未名空间站 (Thu Feb 23 02:00:48 2012, 美东)
做一个关于拉盖尔函数,L(n,x),的计算。
先用matlab。发现对任意固定的x,在n比较大,比如从n=50变成51(或是x比较大的时
候),函数的行为突然变了。变得和n<=50时完全不同。我怀疑是精度不够了。但
matlab并没报错。
于是改用c++。用了long double。结果n稍大就有Floating point exception. 这个错
误是什么意思?
比较了一下c++和matlab生成的拉盖尔函数。仅仅是拉盖尔函数本身。发现对相同的x,
在n<=13的时候结果相同。大于13,结果就不同了,这时候函数值还不到10。 而且c++对n
稍大就出Floating point exception.而matlab貌似还可以处理。
现在问题是该如何得到可靠的结果。我既不相信matlab,也不相信c++了。 |
r****t 发帖数: 10904 | |
h*******s 发帖数: 8454 | 3 http://www.mathworks.com/help/toolbox/symbolic/f1-5556.html
【在 b********r 的大作中提到】 : 【 以下文字转载自 Computation 讨论区 】 : 发信人: bankbuster (恭喜发财), 信区: Computation : 标 题: 怎么提高C++计算精度? C++ vs Matlab : 发信站: BBS 未名空间站 (Thu Feb 23 02:00:48 2012, 美东) : 做一个关于拉盖尔函数,L(n,x),的计算。 : 先用matlab。发现对任意固定的x,在n比较大,比如从n=50变成51(或是x比较大的时 : 候),函数的行为突然变了。变得和n<=50时完全不同。我怀疑是精度不够了。但 : matlab并没报错。 : 于是改用c++。用了long double。结果n稍大就有Floating point exception. 这个错 : 误是什么意思?
|
a****l 发帖数: 8211 | 4 matlab的程序能让你随便写一个什么程序就超过的话,你让mathwork里面那么多的编程
高手和数学phd的脸放到哪里?
【在 b********r 的大作中提到】 : 【 以下文字转载自 Computation 讨论区 】 : 发信人: bankbuster (恭喜发财), 信区: Computation : 标 题: 怎么提高C++计算精度? C++ vs Matlab : 发信站: BBS 未名空间站 (Thu Feb 23 02:00:48 2012, 美东) : 做一个关于拉盖尔函数,L(n,x),的计算。 : 先用matlab。发现对任意固定的x,在n比较大,比如从n=50变成51(或是x比较大的时 : 候),函数的行为突然变了。变得和n<=50时完全不同。我怀疑是精度不够了。但 : matlab并没报错。 : 于是改用c++。用了long double。结果n稍大就有Floating point exception. 这个错 : 误是什么意思?
|
d****n 发帖数: 1637 | 5 try "arbitrary precision" based C/C++programming |
g*****y 发帖数: 7271 | 6 gmp library
【在 d****n 的大作中提到】 : try "arbitrary precision" based C/C++programming
|
l********a 发帖数: 1154 | |
b********r 发帖数: 1080 | 8 谢谢上面几位。我在看gmp等C++处理精度的办法。不过还是要感慨一下,C++处理数值
比matlab差这么多? |
b********r 发帖数: 1080 | 9 还有那个Floating point exception是怎么引起的? |
y**b 发帖数: 10166 | 10 matlab专门搞数值的,这方面一般而言比较强。
可是你要特别高精度计算的话,matlab也未必能满足要求,
而C++有很多高精度库可以调用,这方面C++还是很强。
有很多,quad-precision是其中之一。自然这些库一般
都比较慢,毕竟是用软件去模拟硬件精度。
你的代码里面fact返回类型显然不能用int,而要用long。
另外浮点数一般不要用long double而要用double,
long double除了存储字节更多(我测试过这导致速度比double慢一倍),
而精度没有任何提高。唯一的好处是表示范围大了,能超过E+308/E-308很多。
【在 b********r 的大作中提到】 : 谢谢上面几位。我在看gmp等C++处理精度的办法。不过还是要感慨一下,C++处理数值 : 比matlab差这么多?
|
|
|
r****t 发帖数: 10904 | 11 有没有搞错,matlab 也是 c++ 啊。
【在 b********r 的大作中提到】 : 谢谢上面几位。我在看gmp等C++处理精度的办法。不过还是要感慨一下,C++处理数值 : 比matlab差这么多?
|
h****g 发帖数: 71 | 12 为何要用阶乘?
在能避免的情况下就不要用。
我想最好的方式应该是用递推公式
L_{k+1}(x) = [(2k+1-x)/(k+1)]L_k(x) + [k/(k+1)]L_{k-1}(x)
【在 b********r 的大作中提到】 : 还有那个Floating point exception是怎么引起的?
|
m*******1 发帖数: 58 | 13 试一试下面的程序。 和matlab的结果几乎一样。
----------------------------------------------------
double lagfun(int n, double x){
double v, b;
double p;
int k;
v = 1.0l;
b = 1.0l;
p = 1.0l;
for (k=n; k>=1; k--){
b = b*k/(n+p-k);
v = b - x*v/k;
}
return v;
} |
l***l 发帖数: 22 | 14 There is a better way to calculate Laguerre Polynomial. http://www.mathworks.com/matlabcentral/fileexchange/4912-laguerrepoly-m
As far as your code is concerned, I have 2 suggestions:
1. When you use sum = sum + ..., the rounding error accumulates. This issue
is pointed out in the first chapter of most numerical analysis textbooks. In
your case, try "for (int i=n;i>=0;i--)", it should help, but will not
completely solve the problem.
2. Another trick is, rather than using raw factorial n! which increases very
quickly, use log(n!). |
g******n 发帖数: 253 | 15 c 是通用的。matlab是专用的。不具可比性。可以考虑找个专用的数学库试试。。
【在 b********r 的大作中提到】 : 谢谢上面几位。我在看gmp等C++处理精度的办法。不过还是要感慨一下,C++处理数值 : 比matlab差这么多?
|
w********h 发帖数: 48 | 16 //32位机器上n>=13时溢出
int fact(int n) |
g*****y 发帖数: 7271 | 17 ft,你对阶乘到底有没有概念,居然拿int来存,我真是不知道该说啥了。
【在 b********r 的大作中提到】 : 还有那个Floating point exception是怎么引起的?
|
t****t 发帖数: 6806 | 18 提高计算精度不是一句话的事情, 很多学校专门有课讲数值计算. 你这什么概念也没有
的(比如说, 用int计算阶乘), 还是找本入门书看看吧. 不愿意看也行, 那你就靠
matlab好了, 就别想着自己写了, 这跟c++还是别的什么语言毫无关系.
【在 b********r 的大作中提到】 : 还有那个Floating point exception是怎么引起的?
|
t****t 发帖数: 6806 | 19 在x86-64上, long double用80位, 但是数学运算缺省是用SSE的(which is 64-bit).
你要注意一下.
【在 y**b 的大作中提到】 : matlab专门搞数值的,这方面一般而言比较强。 : 可是你要特别高精度计算的话,matlab也未必能满足要求, : 而C++有很多高精度库可以调用,这方面C++还是很强。 : 有很多,quad-precision是其中之一。自然这些库一般 : 都比较慢,毕竟是用软件去模拟硬件精度。 : 你的代码里面fact返回类型显然不能用int,而要用long。 : 另外浮点数一般不要用long double而要用double, : long double除了存储字节更多(我测试过这导致速度比double慢一倍), : 而精度没有任何提高。唯一的好处是表示范围大了,能超过E+308/E-308很多。
|
h*******s 发帖数: 8454 | 20 搭车问一下
lz的问题是不是直接用matlab
vpa(mfun('L', 20, 0.6), 100)
就行啦?
【在 y**b 的大作中提到】 : matlab专门搞数值的,这方面一般而言比较强。 : 可是你要特别高精度计算的话,matlab也未必能满足要求, : 而C++有很多高精度库可以调用,这方面C++还是很强。 : 有很多,quad-precision是其中之一。自然这些库一般 : 都比较慢,毕竟是用软件去模拟硬件精度。 : 你的代码里面fact返回类型显然不能用int,而要用long。 : 另外浮点数一般不要用long double而要用double, : long double除了存储字节更多(我测试过这导致速度比double慢一倍), : 而精度没有任何提高。唯一的好处是表示范围大了,能超过E+308/E-308很多。
|
|
|
h*********n 发帖数: 11319 | 21 可以试试softfloat 库
【在 b********r 的大作中提到】 : 还有那个Floating point exception是怎么引起的?
|
y**b 发帖数: 10166 | 22 我曾经用intel的编译器(好像也就它支持吧)来利用这个80bit,可是没看到什么效果。
好像是放内存里面的时候才80bit,取出来就64bit了,怎么回事有点记不清了。
但是用quad-precision就不同了,效果非常明显,控制并行计算的roundoff误差很有效。
【在 t****t 的大作中提到】 : 在x86-64上, long double用80位, 但是数学运算缺省是用SSE的(which is 64-bit). : 你要注意一下.
|
t****t 发帖数: 6806 | 23 use x87 for math.
效。
【在 y**b 的大作中提到】 : 我曾经用intel的编译器(好像也就它支持吧)来利用这个80bit,可是没看到什么效果。 : 好像是放内存里面的时候才80bit,取出来就64bit了,怎么回事有点记不清了。 : 但是用quad-precision就不同了,效果非常明显,控制并行计算的roundoff误差很有效。
|
t****t 发帖数: 6806 | 24 also, gcc support that too. for example:
long double foo(long double a, long double b, long double c)
{
return a+b*c;
}
compile with gcc -S -O2:
foo:
.LFB0:
.cfi_startproc
fldt 24(%rsp)
fldt 40(%rsp)
fmulp %st, %st(1)
fldt 8(%rsp)
faddp %st, %st(1)
ret
.cfi_endproc
it's clearly supported.
效。
【在 y**b 的大作中提到】 : 我曾经用intel的编译器(好像也就它支持吧)来利用这个80bit,可是没看到什么效果。 : 好像是放内存里面的时候才80bit,取出来就64bit了,怎么回事有点记不清了。 : 但是用quad-precision就不同了,效果非常明显,控制并行计算的roundoff误差很有效。
|
y**b 发帖数: 10166 | 25 good to know this, thanks.
【在 t****t 的大作中提到】 : also, gcc support that too. for example: : long double foo(long double a, long double b, long double c) : { : return a+b*c; : } : compile with gcc -S -O2: : foo: : .LFB0: : .cfi_startproc : fldt 24(%rsp)
|
y**b 发帖数: 10166 | 26 哦,x87浮点寄存器里面80bit,这是硬件特性,gcc应该支持。
【在 t****t 的大作中提到】 : use x87 for math. : : 效。
|
E*******1 发帖数: 3464 | 27 c++ is the core of matlab, if you are good enough at c++, no difference at
all
【在 b********r 的大作中提到】 : 谢谢上面几位。我在看gmp等C++处理精度的办法。不过还是要感慨一下,C++处理数值 : 比matlab差这么多?
|
m****s 发帖数: 1481 | 28 matlab也不够的时候就用 mathematica 吧,楼主说的那个matlab 50到51位的行为不
同问题我在matlab也遇到过,在mathematica上发现没问题,可以到几百都行为一致。 |
N***m 发帖数: 4460 | 29 很大的数做运算是不能死算的(除非结果就是很大的数。。。但是这种情况很少见),
比如算1000!/999!,你总不会真的去算阶乘吧?
【在 b********r 的大作中提到】 : 【 以下文字转载自 Computation 讨论区 】 : 发信人: bankbuster (恭喜发财), 信区: Computation : 标 题: 怎么提高C++计算精度? C++ vs Matlab : 发信站: BBS 未名空间站 (Thu Feb 23 02:00:48 2012, 美东) : 做一个关于拉盖尔函数,L(n,x),的计算。 : 先用matlab。发现对任意固定的x,在n比较大,比如从n=50变成51(或是x比较大的时 : 候),函数的行为突然变了。变得和n<=50时完全不同。我怀疑是精度不够了。但 : matlab并没报错。 : 于是改用c++。用了long double。结果n稍大就有Floating point exception. 这个错 : 误是什么意思?
|
p**o 发帖数: 3409 | 30 Matlab links against MKL (for dense matrices),
which is written in Fortran.
【在 E*******1 的大作中提到】 : c++ is the core of matlab, if you are good enough at c++, no difference at : all
|