Algorythms/Orthogonal Polynomials/TEST/ExplicitExpressions.cpp

#include "ExplicitExpressions.h"
// Functions from this file are implemented using the book "Handbook of Mathematical Function With Formulas, Graphs and Mathematical Tables" 
// by Milton Abramowitz and Irene A. Stegun, tenth printing, 1972.
// Table 22.3 "Explicit Expressions"


double ExplicitExpressionCheb(const int& type, const int& n, const double& x)
{
    if (n == 0)
    {
        return 1;
    }
    if (type == 1)
    {
        if (n == 1)
        {
            return x;
        }
        int N = n / 2;
        double d = static_cast<double>(n) / 2;
        double c = pow(-1, 0) * std::tgamma(n) / (std::tgamma(1) * std::tgamma(n + 1));
        double g = pow(2*x, n);
        double res = 0;
        for (int m=0; m<=N; ++m)
        {
            c = pow(-1, m) * std::tgamma(n - m) / (std::tgamma(m + 1) * std::tgamma(n - 2*m + 1));
            g = pow(2*x, n - 2*m);
            res += c * g;
        }
    return d * res;
    }
    else if (type == 2)
    {
        if (n == 1)
        {
            return 2*x;
        }
        int N = n / 2;
        double d = 1;
        double c = 1;
        double g = pow(2*x, n);
        double res = 0;
        for (int m=0; m<=N; ++m)
        {
            c = pow(-1, m) * std::tgamma(n - m + 1) / (std::tgamma(m + 1) * std::tgamma(n - 2*m + 1));
            g = pow(2*x, n - 2*m);
            res += c * g;
        }
    return d * res;
    }
    else
    {
        printf("The order must be 1 or 2!");
        std::exit(9);
    }
}

double ExplicitExpressionLaguerre(const int& n, const double& x)
{
    if (n == 0)
    {
        return 1;
    }
    else if (n == 1)
    {
        return 1 - x;
    }
    int N = n;
    double d = 1;
    double c = pow(-1, 0) * std::tgamma(n) / (std::tgamma(1) * std::tgamma(n + 1));
    double g = pow(2*x, n);
    double res = 0;
    for (int m=0; m<=N; ++m)
    {
        c = pow(-1, m) * std::tgamma(n + 1) / (pow(std::tgamma(m + 1), 2) * std::tgamma(n - m + 1));
        g = pow(x, m);
        res += c * g;
    }
    return d * res;
}

double ExplicitExpressionHermite(const int& n, const double& x)
{
    if (n == 0)
    {
        return 1;
    }
    else if (n == 1)
    {
        return 2 * x;
    }
    int N = n / 2;
    double d = std::tgamma(n + 1);
    double c = pow(-1, 0) * std::tgamma(n) / (std::tgamma(1) * std::tgamma(n + 1));
    double g = pow(2*x, n);
    double res = 0;
    for (int m=0; m<=N; ++m)
    {
        c = pow(-1, m) / (std::tgamma(m + 1) * std::tgamma(n - 2 * m + 1));
        g = pow(2 * x, n - 2 * m);
        res += c * g;
    }
    return d * res;
}