I have a code for doing integrations and in the last for loop starting on line 67 I have a for loop which accumulates the values of the function at randomly generated points to get the integral(Monte Carlo integration). Unfortunately after the loop finishes I get NAN as the result for "monte2" variable. I have written a printf statement inside the for loop to pinpoint the mistake only to notice that after 235.494781 the sum turns into -nan. What may be the reason behind this problem? I am running Ubuntu 12.04.3 LTS 32-bit and compile the plain C code with GCC version 4.6.3. I appreciate your help, the code is as follows:
我有一个用于进行积分的代码,并且在第67行开始的最后一个for循环中,我有一个for循环,它在随机生成的点累积函数的值以获得积分(蒙特卡洛积分)。不幸的是,在循环结束后,我得到NAN作为“monte2”变量的结果。我在for循环中写了一个printf语句来查明错误,只是注意到在235.494781之后总和变成了-nan。这个问题背后可能是什么原因?我正在运行Ubuntu 12.04.3 LTS 32位并使用GCC版本4.6.3编译普通C代码。感谢您的帮助,代码如下:
P.S: The code is originally written in Code Blocks on Windows 8 64-bit by a friend of mine if this makes a difference.
P.S:代码最初由我的一位朋友在Windows 8 64位的代码块中编写,如果这有所不同的话。
#include<stdio.h>
#include<math.h>
float I1(float x)
{
return exp(-x)*cos(x*x)*cos(x*x);
}
float I2(float t)
{
return cos(log(t)*log(t))*cos(log(t)*log(t));
}
float random()
{
float a;
a=rand()%1000;
a=a/1000*20;
// printf("%.15f\t%f\n",I1(a),a);
return a;
}
float random2()
{
float a;
a=rand()%1000;
a/=1000;
// printf("%.15f\t%f\n",I2(a),a);
return a;
}
int main()
{
FILE *data=fopen("data.txt","w");
FILE *data2=fopen("data2.txt","w");
float trap=0,monte=0,sum=0, monte2=0;
float a[1000],b[1000],dt=0.005;
int i;
/* Part 1 */
for(i=0;i<1000;i++)
a[i]=I1(i*dt);
for(i=0;i<1000;i++)
fprintf(data,"%f\t%f\n",i*dt,a[i]);
for(i=1;i<1000;i++)
trap+=(a[i]+a[i-1])/2*dt;
printf("The integral value of I1 is = %f with trapezoid rule\n",trap);
for(i=0;i<500;i++)
monte+=I1(random());
printf("The Monte Carlo Technique value for I1 is %f with 500 samples\n",monte/500*20);
/* Part 2 */
dt=0.001;
printf(" \n");
for(i=1;i<=1000;i++)
b[i]=I2(i*dt);
for(i=1;i<=1000;i++)
fprintf(data2,"%f\t%f\n",i*dt,b[i]);
for(i=2;i<=1000;i++)
trap+=(b[i]+b[i-1])/2*dt;
printf("The integral value of I2 is = %f with trapezoid rule\n",trap/2);
for(i=0;i<500;i++)
{
monte2+=I2(random2());
printf("%f \n", monte2);
}
printf("The Monte Carlo Technique value of I2 is %f with 500 samples\n",monte2/500);
printf("\n");
printf("Comment 1: Two values obtained with trapezoid rule is close to each other;however,they are not exactly same.\n");
printf("\n");
printf("Comment 2: The integral value and monte carlo value of I1 is closer than the integral value and monte carlo value of I2.This means that we have better expectation value of I1 with monte carlo technique with 500 samples.\n");
fclose(data2);
fclose(data);
return 0;
}
1 个解决方案
#1
5
Your function call
你的函数调用
monte2+=I2(random2());
may produce NaN
. This is because random2
may returns 0
. log 0
is infinity. This will cause cos(log(t)*log(t))*cos(log(t)*log(t))
to produce NaN
.
可能会产生NaN。这是因为random2可能返回0. log 0是无穷大。这将导致cos(log(t)* log(t))* cos(log(t)* log(t))产生NaN。
See the graph for log
function:
查看日志功能图:
Note that the graph gets arbitrarily close to the y axis, but does not meet or intersect it1.
请注意,图形任意接近y轴,但不与it1相交或相交。
1. Source Wikipedia
1.来源*
#1
5
Your function call
你的函数调用
monte2+=I2(random2());
may produce NaN
. This is because random2
may returns 0
. log 0
is infinity. This will cause cos(log(t)*log(t))*cos(log(t)*log(t))
to produce NaN
.
可能会产生NaN。这是因为random2可能返回0. log 0是无穷大。这将导致cos(log(t)* log(t))* cos(log(t)* log(t))产生NaN。
See the graph for log
function:
查看日志功能图:
Note that the graph gets arbitrarily close to the y axis, but does not meet or intersect it1.
请注意,图形任意接近y轴,但不与it1相交或相交。
1. Source Wikipedia
1.来源*