b********r
发帖数: 1080 | 1 做一个关于拉盖尔函数，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++了。 | l********a
发帖数: 1154 | | m*******1
发帖数: 58 | 3 Matlab 还是比较可靠的。应该采用比较简洁的计算方法。 看下面的程序， 与Matlab
的结果基本一样( 差别< 10^(-11) ).
------------------------------------------------------------
#include
#include
#include
double lagfun(int n, double x);
int main(int argc, char **argv){
double x, y;
int k, n;
x = 0.6l;
n = 100;
for (k=0; k<=n; k++){
y = lagfun(k, x);
printf("k=%3d, \t %26.14f\n", k, y);
}
return 0;
}
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;
} | b********r
发帖数: 1080 | 4 做一个关于拉盖尔函数，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++了。
--------更新---------
C++程序摘要如下。我知道子程序lag里那几个阶乘可以化简一下。但还是不相信C++
居然连n=20都算不了。这程序里有什么错吗？
int main()
{
.....
int n=20;
long double x=0.6;
cout<
.....
}
// factorial function
int fact(int n)
{
if (n>1)
{
return n*fact(n-1);
}
return 1;
}
// Laguerre(n,x)
long double lag(int n , long double x) //
{
long double sum=0.0;
for (int i=0;i
{
sum=sum+fact(n)/fact(n-i)/pow(fact(i),2)*pow(-x,i);
}
return sum;
} | l********a
发帖数: 1154 | | m*******1
发帖数: 58 | 6 Matlab 还是比较可靠的。应该采用比较简洁的计算方法。 看下面的程序， 与Matlab
的结果基本一样( 差别< 10^(-11) ).
------------------------------------------------------------
#include
#include
#include
double lagfun(int n, double x);
int main(int argc, char **argv){
double x, y;
int k, n;
x = 0.6l;
n = 100;
for (k=0; k<=n; k++){
y = lagfun(k, x);
printf("k=%3d, \t %26.14f\n", k, y);
}
return 0;
}
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;
} | m*****e
发帖数: 15 | 7 第一个程序要算阶乘，用的是整型，13!的时候就溢出了。 | g**********t
发帖数: 475 | 8 Yes, overflow.
Lz might try to calculate in log scale and then transform the final result
back to linear scale.
【在 m*****e 的大作中提到】
: 第一个程序要算阶乘，用的是整型，13!的时候就溢出了。
| g**********t
发帖数: 475 | |
|
