C ++中的高斯传奇

Question

我正在尝试根据以下算法创建 gauss-legendre 代码：

对于 n 个点

也就是说，它创建了一个 2n 方程系统（如果我们要求对 2n-1 阶的多项式准确，

ti 是 n 阶legendre 多项式的根。legendre 多项式给出：

和无线：

我的代码是：

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <iomanip>
#include <cmath>


using namespace std;

const double pi=3.14;


//my function with limits (-1,1)
double f(double x){
double y;
y=(pi/4.0)*(log((pi*(x+1.0))/4.0 +1.0));
return y;

}

double legendre (int n){

    double *L,*w,*t;
    double x,sum1,sum2,result;
    L=new double [n];
    w=new double [n];
    t=new double [n];


        while(n<10){

         L[0]=1;
         L[1]=x;


        //legendre coef
        for (int i=1;i<=10;i++){
        L[i+1]=((2.0*i+1.0)*x*L[i] - i*L[i-1])/(i+1.0);


        }

        //weights w
        w=0;
        for (int i=1;i<=10;i++){
        w[i]+=(2.0*(1.0-x*x))/(i*i*(L[i-1]*L[i-1]));
        }


        //sums  w*t
        for (int i=1;i<=10;i++){
            sum1=0.0; //for k=1,3,5,2n-1
            for (int k=1;k<=2*n-1;k+=2){
                sum1+=w[i]*(pow(t[i],k));
            }
                sum1=0;
                sum2=0.0;//for k=0,2,4,2n-2
                for(int k=0;k<=2*n-2;k+=2){
                    sum2+=w[i]*(pow(t[i],k));
                }
                sum2=2.0/n;
        }


    }

    result=w*f(*t);

    return result;

}


int main()
{
    double eps=1e-8;//accuracy
    double exact=0.8565899396;//exact solution for the integral
    double error=1.0;
    double result;

    int n=1;//initial point


    while (fabs(result-exact)>eps) {
        result=legendre(n);
        cout <<"\nFor n = "<<n<<",error = "<<fabs(error-exact)<<",value = "<<result;

    n++;
    }

    return 0;
}

我的问题是：

1) 编译器在结果=w*f(*t); 处给我：错误：'double*' 和'double' 类型的无效操作数到二进制'operator*' -->

2) 我不确定我是否做对了整个事情。我的意思是，如果我将所有事情结合在一起并且我是否正确实施了算法。

Answer 1

我不知道算法，但你的代码是错误的。
第一：

        while(n<10)
        {
         L[0]=1;
         L[1]=x;
        //legendre coef
        for (int i=1;i<=10;i++){
        L[i+1]=((2.0*i+1.0)*x*L[i] - i*L[i-1])/(i+1.0);
        }

您必须更改n 的值（递增、递减等），否则这是一个无限循环。

第二：

//weights w
    w=0;
    for (int i=1;i<=10;i++){
    w[i]+=(2.0*(1.0-x*x))/(i*i*(L[i-1]*L[i-1]));
    }

w 是一个指针。如果你想重置它，使用memset(w,0,sizeof(double)*n);不要让它等于0。

最后：

result=w*f(*t);

由于您将w 和t 指针用作数组，因此您必须提供某种索引，例如result=w[ind] * f(t[ind]);。在这里，您只是将指针 w 的值，而不是 w 指向的值（顺便说一下，w 的值从您的代码中为 0）与 @ 指向的双精度数组的第一个值相乘987654333@.

我也无法得到您的代码与问题中的公式之间的任何关系。所以我的谦虚建议是不要为此使用 C 或 C++。如果必须，请不要使用指针，因为您似乎不熟悉它们。您可以轻松地使用 std::vector 代替指针。

Answer 2

w 是一个指针，你正试图将它与某个东西相乘......你必须使用索引

w[index] * f(*t)

*t 也是 t 数组的第一个元素。是这个意思吗？

Answer 3

关于您的算法，x（横坐标值）应该是勒让德多项式的零点。我没有看到你在任何地方定义它们。定义它们有点痛苦。我正在做类似的事情，发现这个（它是一个 Matlab 文件，而不是一个 C++ 文件）定义了 N xi 和 wi 值。该算法运行良好：http://www.mathworks.com/matlabcentral/fileexchange/4540-legendre-gauss-quadrature-weights-and-nodes

Answer 4

alglib 和 gsl 都有高斯正交的实现。不过，两者都是 GPL 许可的。