[58acde]: src / legendre_sphPlm_array.cc Maximize Restore History

Download this file

legendre_sphPlm_array.cc    79 lines (67 with data), 2.5 kB

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
/*
Copyright (C) 2004 Teemu Ikonen <tpikonen@pcu.helsinki.fi>
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program; If not, see <http://www.gnu.org/licenses/>.
Changes:
6/18/2008 R. Rogers Removed unsed variable i
*/
#include <octave/oct.h>
#include <gsl/gsl_sf.h>
#include <gsl/gsl_sf_result.h>
#include <gsl/gsl_errno.h>
void octave_gsl_errorhandler (const char * reason, const char * file,
int line, int gsl_errno)
{
error("GSL error %d at %s, line %d: %s\n", gsl_errno, file, line, reason);
}
DEFUN_DLD(legendre_sphPlm_array, args, nargout, "\
-*- texinfo -*-\n\
@deftypefn {Loadable Function} {@var{y} =} legendre_sphPlm_array (@var{lmax}, @var{m}, @var{x})\n\
\n\
This function computes an array of normalized associated Legendre functions\n\
@iftex\n\
@tex\n\
$\\sqrt{(2l+1)/(4\\pi)} \\sqrt{(l-m)!/(l+m)!} P_l^m(x)$\n\
for $m >= 0, l = |m|, ..., lmax, |x| <= 1.0$\n\
@end tex\n\
@end iftex\n\
@ifinfo\n\
sqrt((2l+1)/(4*pi)) * sqrt((l-m)!/(l+m)!) Plm (x)\n\
for m >= 0, l = |m|, ..., lmax, |x| <= 1.0\n\
@end ifinfo\n\
\n\
This function is from the GNU Scientific Library,\n\
see @url{http://www.gnu.org/software/gsl/} for documentation.\n\
@end deftypefn")
{
gsl_set_error_handler (octave_gsl_errorhandler);
if(args.length() != 3) {
print_usage ();
return octave_value_list();
}
if(!args(0).is_real_type()
|| !args(1).is_real_type()
|| !args(2).is_real_type()
|| !args(0).is_scalar_type()
|| !args(1).is_scalar_type()
|| !args(2).is_scalar_type()) {
error("The arguments must be real scalars.");
print_usage ();
return octave_value_list();
}
int lmax = static_cast<int>(args(0).scalar_value());
int m = static_cast<int>(args(1).scalar_value());
double x = args(2).scalar_value();
RowVector *y = new RowVector(gsl_sf_legendre_array_size (lmax, m));
double *yd = y->fortran_vec();
gsl_sf_legendre_sphPlm_array (lmax, m, x, yd);
return octave_value_list(octave_value(*y));
}