Рекурсивная программа для получения полиномов Лежандра

Question

Я пытаюсь написать на C ++ функцию для конструктора полиномов Лежандра, которая выводит коэффициенты полинома степени m . Полиномы следуют за простым рекурсивным отношением .

Теперь я пытаюсь реализовать это отношение, и оно отлично работает для каждого n до 8, начиная с 9,вокруг 7-й итерации он выбирает значение мусора, которого нет ни в одном из предыдущих векторов коэффициентов. Интересно, как я могу решить это. Я представляю вам мой код:

#include <vector>
#include <cmath>
std::vector<double> set_coeffs(int m){
std::vector<double> coeffs;

if (m == 0) //Casos iniciales para empezar la recursión
{               
    coeffs.push_back(1);
} else if (m == 1)
{
    coeffs.push_back(1);
    coeffs.push_back(0);
} else if (m == 2) //Puse también el caso 2 porque de otro modo ocurre el mismo problema pero desde n=5
{
    coeffs.push_back(1.5);
    coeffs.push_back(0);
    coeffs.push_back(-0.5);
} 
 else
{
    std::vector<double> v = set_coeffs(m-1);
    std::vector<double> u = set_coeffs(m-2);
    std::cout << "inicia cicle" << std::endl;

    double a = (2* ((double)m) -1)/((double)m);
    double b = (((double)m)-1)/((double)m);

    coeffs.push_back(a*v[0]);
    coeffs.push_back(a*v[1]);

    for (int i = 0; i < m-1; i++)
    {
        double c = a*v[i+2] - b*u[i];
        std::cout << m << " " << v[i+2] << " " << u[i] << " " << c <<std::endl;
        coeffs.push_back(c);
    }

    std::cout << "termina cicle" << std::endl;
} 
return coeffs;
}

Answer 1

Вызов set_coeffs(n) возвращает вектор с n + 1 элементами.

Следовательно, std::vector<double> v = set_coeffs(m-1); имеет m элементов.

В вашем цикле for, i начинается с 0 до m-2, и вы получаете доступ к v[i+2]. На последней итерации будет получен доступ к v[m], который находится за пределами.

Answer 2

Как отмечено в interjay * answer , отправленный код неправильно учитывает различия в размерах вектора.

Чтобы сделатьсравнение индексов легче (по крайней мере, ИМХО, но это также должно помочь в дальнейших вычислениях с использованием этих коэффициентов), я бы сохранял коэффициенты в обратном порядке. Я бы также воспользовался их конкретным шаблоном (только нечетные или , четные отличны от нуля).

Возможно, стоит упомянуть, что рекурсивный подход такого типабудет неэффективным, когда применяется к большим векторам.

Тем не менее, это возможная реализация:

std::vector<double> Legendre_coefficients(int m)
{
    if (m == 0) 
    {               
        return {1};
    }
    if (m == 1)
    {
        return {0, 1};
    }

    // Initialize with zero, only at most (half + 1) of the terms will be changed later
    std::vector<double> coeffs(m + 1);

    // Consider some form of memoization instead of this recursion
    std::vector<double> v = Legendre_coefficients(m - 1);
    std::vector<double> u = Legendre_coefficients(m - 2);

    // using literals of floating point type, 'm' is promoted by the compiler
    double a = (2.0 * m - 1.0) / m;
    double b = (m - 1.0) / m;

    int first = 1;
    // If 'm' is even, so is (m - 2) and its first element is zero. It can be skipped.
    // It also avoids annoying '-0' in the output 
    if ( m % 2 == 0 )
    {
        coeffs[0] = -b * u[0];
        first = 2;
    }
    for (int i = first; i < m - 1; i += 2)
    {
        coeffs[i] = (a * v[i - 1] - b * u[i]);
    }
    coeffs[m] = a * v[m - 1];

    return coeffs;
}

Testable здесь .