[C++] Kwadratury Gausa-Hermite

pag15
Użytkownik
Użytkownik
Posty: 70
Rejestracja: 2 gru 2014, o 01:26
Płeć: Mężczyzna
Lokalizacja: Koszalin
Podziękował: 10 razy

[C++] Kwadratury Gausa-Hermite

Post autor: pag15 »

Mam problem, ponieważ muszę zaimplementować w C++ metodę całkowania numerycznego Gausa-Hermite , posiadam surowy kod z książki ,, Numerical Recipesin C' lecz nie jestem w stanie zrozumieć co on zwraca i jak go poprawnie wywołać, mógłby mi ktoś pomóc, chciałbym zwracać wynik całki (poniżej zamieszczam procedurę):

Kod: Zaznacz cały

#include <cmath>
#define EPS 3.0e-14
#define PIM4 0.7511255444649425
#define MAXIT 10

void gauher(double *x, double *w, int n)
{

    int its,j,m;
    double p1,p2,p3,pp,z,z1;

    m =(n+1)/2;

    for (int i = 1 ; i <= m ; i++)
    {
        if (i == 1)
        {
            z = sqrt((double)(2*n+1))-1.85575*pow((double)(2*n+1),-0.16667);
        }
        else if (i == 2)
        {
            z -= 1.14*pow((double)n,0.426)/z;
        }
        else if (i == 3)
        {
            z = 1.86*z-0.86*x[1];
        }
        else if (i == 4)
        {
            z=1.91*z-0.91*x[2];
        }
        else
        {
            z=2.0*z-x[i-2];
        }
        for (its = 1; its <= MAXIT ; its++)
        {
            p1=PIM4;
            p2=0.0;
            for (j=1; j<=n; j++)
            {
                p3=p2;
                p2=p1;
                p1=z*sqrt(2.0/j)*p2-sqrt(((double)(j-1))/j)*p3;
            }

            pp = sqrt((double)2*n)*p2;
            z1 = z;
            z = z1-p1/pp;
            if (fabs(z-z1) <= EPS)
            {
                break;
            }
        }
        if (its > MAXIT)
        {
            cout << "Za duzo iteracji";
        }

        x[i]=z;
        x[n+1-i] = -z;
        w[i]=2.0/(pp*pp);
        w[n+1-i]=w[i];
    }

}
#undef EPS
#undef PIM4
#undef MAXIT
Ostatnio zmieniony 4 maja 2016, o 09:11 przez Afish, łącznie zmieniany 1 raz.
Powód: Poprawa wiadomości.
ODPOWIEDZ