我正在编写一个程序,该程序将以分段线性方式利用拉格朗日标准格式来插值n次多项式。我的代码在第一个子时间间隔以及第三个和第四个子时间间隔都能正常工作,但是由于某种原因,我无法弄清楚,我收到NaN作为第二个子时间间隔的输出。第二个子间隔在// P2注释行下计算。我已经绞尽脑汁,尝试了所有我能想到的改变来解决问题,但是没有运气。如果有人能提供一些见识,我将不胜感激。请注意,由于第三和第四个插值多项式遵循类似的方式,因此我仅包含了直到第二个插值多项式为止的代码。对于代码的残酷性,我事先表示歉意。我对C ++相对较新,还没有时间获得这种问题可能带来的优雅。再次感谢。
ofstream Outfile;
Outfile.open ("PiecewiseLagrange_D.dat");
double *P1 = new double[201]; //Polynomial 1
double *P2 = new double[201]; //Polynomial 2
double *P3 = new double[201]; //Polynomial 3
double *P4 = new double[201]; //Polynomial 4
double *x = new double[201]; //Interpolating points/x's
double *x1 = new double[(int)n+1]; //First subinterval/mesh/xi's
double *x2 = new double[(int)n+1]; //Second subinterval
double *x3 = new double[(int)n+1]; //Third subinterval
double *x4 = new double[(int)n+1]; //Fourth subinterval
double a, b; //interval end points
char func; //function selection
double xDifference1;
double xDifference2;
cout << "Enter an interval with integer end points (lesser value first)";
cin >> a >> b;
for (int i=0; i<=n; i++) //Initialize
{
P1[i] = 0;
P2[i] = 0;
P3[i] = 0;
P4[i] = 0;
x[i] = 0;
x1[i] = 0;
x2[i] = 0;
x3[i] = 0;
x4[i] = 0;
}
x1[0] = a;
for (int i=0; i<=n; i++)
{
x1[i] = x1[0] + i*(((b-a)/4)/n);
cout << x1[i] << endl;
}
cout << endl;
x2[0] = x1[(int) n];
for (int i=0; i<=n; i++)
{
x2[i] = x2[0] + i*(((b-a)/4)/n);
cout << x2[i] << endl;
}
cout << endl;
x3[0] = x2[(int) n];
for (int i=0; i<=n; i++)
{
x3[i] = x3[0] + i*(((b-a)/4)/n);
cout << x3[i] << endl;
}
cout << endl;
x4[0] = x3[(int) n];
for (int i=0; i<=n; i++)
{
x4[i] = x4[0] + i*(((b-a)/4)/n);
cout << x4[i] << endl;
}
cout << "Enter a function to evaluate (1,2, or 3):";
cin >> func;
//cout << "Polynomial is g1(x) on [" << a << "," << b << "]" << endl;
if (func == '1')
{
//P1
x[0] = a;
for (int i=0; i<=200; i++)
{
x[i] = x[0] + i*((x1[(int) n] - x1[0])/200);
}
for (int j=0; j<=200; j++)
{
xDifference1 = 0;
xDifference2 = 0;
for (int i=0; i<=n; i++)
{
xDifference1 = (x1[i] - x1[i+1]);
xDifference2 = (x1[i+1] - x1[i]);
P1[j] = F1(x1[i])*((x[j] - x1[i+1])/xDifference1) + F1(x1[i+1])*((x[j] - x1[i])/xDifference2);
}
Outfile << x[j] << " " << P1[j] << " " << F1(x[j]) << endl;
cout << setw(8) << x[j] << setw(12) << P1[j] << endl;
}
cout << endl;
//P2
x[0] = x1[(int) n];
for (int i=0; i<=200; i++)
{
x[i] = x[0] + i*((x2[(int) n] - x2[0])/200);
}
for (int j=0; j<=200; j++)
{
xDifference1 = 0;
xDifference2 = 0;
for (int i=0; i<=n; i++)
{
xDifference1 = (x2[i] - x2[i+1]);
xDifference2 = (x2[i+1] - x2[i]);
P2[j] = F1(x2[i])*((x[j] - x2[i+1])/xDifference1) + F1(x2[i+1])*((x[j] - x2[i])/xDifference2);
}
Outfile << x[j] << " " << P2[j] << " " << F1(x[j]) << endl;
cout << setw(8) << x[j] << setw(12) << P2[j] << " " << F1(x[j]) << endl;
}
cout << endl;
最佳答案
由于您没有调用任何其他函数,因此当您将零除以零(0.0 / 0.0
)时会得到NaN。在某些时候,您的xDifference1
和/或xDifference2
为零。
用非零除以零得到无穷大。
编辑但是,显然不是这样,进一步的研究表明,各种x
数组(包括x2
)中都有n+1
元素,它们通过0
索引到n
。在循环期间,您访问x2[i+1]
。由于i
在上一次迭代中将等于n
,因此您访问元素x2[n+1]
,该元素超出数组的范围,并导致未定义的行为。在这种情况下,该数组之后的随机内存会为x2
生成一个NaN,但不会为其他数组生成一个NaN。
在不相关的注释中,您为每个迭代分配给i
的内部P2[j]
循环,因此从循环中获得的唯一值是来自上一次迭代。您是要使用P2[j] += ...
吗?
关于c++ - 为什么我得到NaN作为计算结果? (带拉格朗日标准格式的分段线性插值),我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/33049199/