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
|
