Проблемы с GSL с использованием полиномов Лежандра
Я пытаюсь обновить старый код, в котором использовалась версия GSL с устаревшими функциями, но мне не удается найти, как использовать новую версию функций нормализованных полиномов Лежандра. Вот фрагмент, который суммирует проблему:
#include <iostream>
#include <gsl/gsl_sf_legendre.h>
#include <cmath>
#define GSL_NEW
using namespace std;
int main() {
int order = 17;
int ntheta = 36;
double theta_step = M_PI / ntheta;
double c, theta;
double legendre[ntheta][order+1];
for( int m = 0; m <= order; m += 2) {
for(int l = m; l <= ntheta; l += 2 ) {
for( int t = 0; t < ntheta; t++ ) {
theta = ( ntheta + 0.5 ) * theta_step;
c = cos(theta);
if ( l == m ) {
#ifdef GSL_NEW
gsl_sf_legendre_array( GSL_SF_LEGENDRE_SPHARM, order, c, &legendre[t][l] );
cout << legendre[t][l] << endl;
#else
gsl_sf_legendre_sphPlm_array(order, m, c, &legendre[t][l] );
cout << legendre[t][l] << endl;
#endif
}
}
}
}
}
Когда я компилирую с использованием GSL 1.9, я использую устаревшую функцию gsl_sf_legendre_sphPlm_array, а когда я вычисляю с помощью GSL 2.5, я использую новую функцию gsl_sf_legendre_array, которая явно требует ключевого слова для нормализации (GSL_SF_LEGENDRE_SPHARM) и не запрашивает параметр m. Старая версия дает мне стабильные результаты, а новая возвращает мне ошибку сегментации после 25 циклов. Может ли кто-нибудь из вас помочь мне исправить код и объяснить, что я делаю неправильно?
Ответы 1
Я попытался скомпилировать с отладочными символами (флаг -g в приведенной ниже команде, который предполагает, что ваша программа называется «main.cpp»)
usr@cmptr $ g++ -g main.cpp -lgsl -lgslcblas -lm -o main
И запустил программу с помощью отладчика, например. gdb (отладчик GNU):
usr@cmptr $ gdb main
(gdb) run main
Starting program: /home/usr/Desktop/main
0, 0, 0.282095
1, 0, 0.282095
2, 0, 0.282095
3, 0, 0.282095
4, 0, 0.282095
5, 0, 0.282095
6, 0, 0.282095
7, 0, 0.282095
8, 0, 0.282095
9, 0, 0.282095
10, 0, 0.282095
11, 0, 0.282095
12, 0, 0.282095
13, 0, 0.282095
14, 0, 0.282095
15, 0, 0.282095
16, 0, 0.282095
17, 0, 0.282095
18, 0, 0.282095
19, 0, 0.282095
20, 0, 0.282095
21, 0, 0.282095
22, 0, 0.282095
23, 0, 0.282095
24, 0, 0.282095
Program received signal SIGSEGV, Segmentation fault.
0x0000555555554ded in main () at main.cpp:26
26 cout << t << ", " << l << ", " << legendre[t][l] << endl;
(gdb) quit
Ошибка в строке 26, где вы пишете на экран. Причина, вероятно, в том, что вы пытаетесь получить доступ к массиву за его пределами. Может быть, посмотреть, как «правильно» выделить память в руководство по полиномам Лежандра gsl. Я специально думаю о функциях gsl_sf_legendre_array_n и gsl_sf_legendre_array_index.
Note: I have not used this part of GSL myself, but I usually find that when using GSL it is useful to use temporary variables and arrays that you pack/unpack before/after calling functions. Not pretty but less error prone, and since you are using the "plusspluss" version of C you can always wrap the implementation in a
class. Hopefully it is helpful.
Редактировать:
Я попытался увеличить размер массива legendre, и программа завершается, когда размер установлен на:
double legendre[ntheta + 10][order + 4];
Может быть, кто-то, кто знает, как работают функции, может ответить, почему...