我有这个方程式

然后从下面找到多项式

我尝试这样实现它:
for (int n=0;n<order;n++){
df[n][0]=y[n];
for (int i=0;i<N;i++){ //N number of points
df[n][i]+=factorial(n,i)*y[i+n-1];
}
}
for (int i=0;i<N;i++){
term=factorial(s,i);
result*=df[0][i]*term;
sum+=result;
}
return sum;1)我不确定如何实现function.As中每个参数的符号,你可以看到它变成了‘正’,‘负’,‘正’……
2)我不确定有没有什么错误……
谢谢!
int fact(int n){
//3!=1*2*3
if (n==0) return 1;
else
return n*fact(n-1);
}
double factorial(double s,int n){
//(s 3)=s*(s-1)*(s-2)/6
if ((n==0) &&(s==0)) return 1;
else
return fact(s)/fact(n);
}发布于 2013-01-15 09:58:14
嗯,我知道你想要近似计算给定x=X的值f(x),使用等距点的牛顿插值多项式(更具体地说,牛顿-格雷戈里正向差分插值多项式)。假设s=(X- x0 )/h,其中x0是第一个x,h是获得x的其余部分的步骤,对于这些步骤,您知道f的确切值:
double coef (double s, int k)
{
double c(1);
for (int i=1; i<=k ; ++i)
c *= (s-i+1)/i ;
return c;
}
double P_interp_value(double s, int Num_of_intervals , double f[] /* values of f in these points */) // P_n_s
{
int N=Num_of_intervals ;
double *df0= new double[N+1]; // calculing df only for point 0
for (int n=0 ; n<=N ; ++n) // n here is the order
{
df0[n]=0;
for (int k=0, sig=-1; k<=n; ++k, sig=-sig) // k here is the "x point"
{
df0[n] += sig * coef(n,k) * f[n-k];
}
}
double P_n_s = 0;
for (int k=0; k<=N ; ++k ) // here k is the order
{
P_n_s += coef(s,k)* df0[k];
}
delete []df0;
return P_n_s;
}
int main()
{
double s=0.415, f[]={0.0 , 1.0986 , 1.6094 , 1.9459 , 2.1972 };
int n=1; // Num of interval to use during aproximacion. Max = 4 in these example
while (true)
{
std::cin >> n;
std::cout << std::endl << "P(n=" << n <<", s=" << s << ")= " << P_interp_value(s, n, f) << std::endl ;
}
}它打印:
1
P(n=1,s=0.415)= 0.455919
2
P(n=2,s=0.415)= 0.527271
3.
P(n=3,s=0.415)= 0.55379
4.
P(n=4,s=0.415)= 0.567235
比较:http://ecourses.vtu.ac.in/nptel/courses/Webcourse-contents/IIT-KANPUR/Numerical%20Analysis/numerical-analysis/Rathish-kumar/rathish-oct31/fratnode8.html
它起作用了。现在我们可以开始优化这些代码了。
发布于 2013-01-15 01:42:01
最简单的解决方案可能是将符号保存在一个变量中,然后在每次循环中将其相乘。类似于:
sign = 1.0;
for ( int i = 0; i < N; ++ i ) {
term = factorial( s, i );
result *= df[0][i] * term;
sum += sign * result;
sign = - sign;
}发布于 2013-01-15 02:21:38
您不能执行pow( -1, m )。
您可以编写自己的代码:
inline int minusOnePower( unsigned int m )
{
return (m & 1) ? -1 : 1;
}您可能想要构建一些计算值的表。
https://stackoverflow.com/questions/14323276
复制相似问题