Next: , Previous: bilinear interpolation, Up: ncap2 netCDF Arithmetic Processor


4.1.19 GSL special functions

As of version 3.9.6 (released January, 2009), NCO can link to the GNU Scientific Library (GSL). ncap can access most GSL special functions including Airy, Bessel, error, gamma, beta, hypergeometric, and Legendre functions and elliptical integrals. GSL must be version 1.4 or later. To list the GSL functions available with your NCO build, use ncap2 -f | grep ^gsl.

The function names used by ncap2 mirror their GSL names. The NCO wrappers for GSL functions automatically call the error-handling version of the GSL function when available 1. This allows NCO to return a missing value when the GSL library encounters a domain error or a floating point exception. The slow-down due to calling the error-handling version of the GSL numerical functions was found to be negligible (please let us know if you find otherwise).

Consider the gamma function.
The GSL function prototype is
int gsl_sf_gamma_e(const double x, gsl_sf_result * result) The ncap script would be:

     lon_in[lon]={-1,0.1,0,2,0.3};
     lon_out=gsl_sf_gamma(lon_in);
     lon_out= _, 9.5135, 4.5908, 2.9915

The first value is set to _FillValue since the gamma function is undefined for negative integers. If the input variable has a missing value then this value is used. Otherwise, the default double fill value is used (defined in the netCDF header netcdf.h as NC_FILL_DOUBLE = 9.969e+36).

Consider a call to a Bessel function with GSL prototype
int gsl_sf_bessel_Jn_e(int n, double x, gsl_sf_result * result)

An ncap script would be

     lon_out=gsl_sf_bessel_Jn(2,lon_in);
     lon_out=0.11490, 0.0012, 0.00498, 0.011165

This computes the Bessel function of order n=2 for every value in lon_in. The Bessel order argument, an integer, can also be a non-scalar variable, i.e., an array.

     n_in[lon]={0,1,2,3};
     lon_out=gsl_sf_bessel_Jn(n_in,0.5);
     lon_out= 0.93846, 0.24226, 0.03060, 0.00256

Arguments to GSL wrapper functions in ncap must conform to one another, i.e., they must share the same sub-set of dimensions. For example: three_out=gsl_sf_bessel_Jn(n_in,three_dmn_var_dbl) is valid because the variable three_dmn_var_dbl has a lon dimension, so n_in in can be broadcast to conform to three_dmn_var_dbl. However time_out=gsl_sf_bessel_Jn(n_in,time) is invalid.

Consider the elliptical integral with prototype int gsl_sf_ellint_RD_e(double x, double y, double z, gsl_mode_t mode, gsl_sf_result * result)

     three_out=gsl_sf_ellint_RD(0.5,time,three_dmn_var_dbl);

The three arguments are all conformable so the above ncap call is valid. The mode argument in the function prototype controls the convergence of the algorithm. It also appears in the Airy Function prototypes. It can be set by defining the environment variable GSL_PREC_MODE. If unset it defaults to the value GSL_PREC_DOUBLE. See the GSL manual for more details.

     export GSL_PREC_MODE=0 // GSL_PREC_DOUBLE
     export GSL_PREC_MODE=1 // GSL_PREC_SINGLE
     export GSL_PREC_MODE=2 // GSL_PREC_APPROX

The ncap wrappers to the array functions are slightly different. Lets consider the following gsl prototype
int gsl_sf_bessel_Jn_array(int nmin, int nmax, double x, double *result_array)

     b1=lon.double();
     x=0.5;
     status=gsl_sf_bessel_Jn_array(1,4,x,&b1);
     print(status);
     b1=0.24226, 0.0306, 0.00256, 0.00016 ;
     

This calculates the bessel function of x=0.5 for n=1 to 4. The first three arguments are scalar values. if a non-scalar variable is supplied as an argument then only the first value is used.The final argument is the variable where the results go ( note the '&' this indicates a call by reference). This final argument must be of type double and must be of least size (nmax-nmin+1). If either of these conditions are not met then then the function will blow out with an error message. The function/wrapper returns a status flag. Zero indicates success.

Lets look at another array function
int gsl_sf_legendre_Pl_array( int lmax, double x, double *result_array);

     a1=time.double();
     x=0.3;
     status=gsl_sf_legendre_Pl_array(a1.size()-1, x,&a1);
     print(status);

This call calculates P_l(0.3) for l=0..9. Note |x|<=1, otherwise there will be a domain error. See the GSL documentation for more details.

Below is table detailing what GSL functions have been implemented. This table is correct for GSL version 1.10. To see what functions are available on your build run the command ncap2 -f |grep ^gsl . To see this table along with the GSL C function prototypes look at the spreadsheet doc/nco_gsl.ods.

GSL NAME I NCAP FUNCTION CALL
gsl_sf_airy_Ai_e Y gsl_sf_airy_Ai(dbl_expr)
gsl_sf_airy_Bi_e Y gsl_sf_airy_Bi(dbl_expr)
gsl_sf_airy_Ai_scaled_e Y gsl_sf_airy_Ai_scaled(dbl_expr)
gsl_sf_airy_Bi_scaled_e Y gsl_sf_airy_Bi_scaled(dbl_expr)
gsl_sf_airy_Ai_deriv_e Y gsl_sf_airy_Ai_deriv(dbl_expr)
gsl_sf_airy_Bi_deriv_e Y gsl_sf_airy_Bi_deriv(dbl_expr)
gsl_sf_airy_Ai_deriv_scaled_e Y gsl_sf_airy_Ai_deriv_scaled(dbl_expr)
gsl_sf_airy_Bi_deriv_scaled_e Y gsl_sf_airy_Bi_deriv_scaled(dbl_expr)
gsl_sf_airy_zero_Ai_e Y gsl_sf_airy_zero_Ai(uint_expr)
gsl_sf_airy_zero_Bi_e Y gsl_sf_airy_zero_Bi(uint_expr)
gsl_sf_airy_zero_Ai_deriv_e Y gsl_sf_airy_zero_Ai_deriv(uint_expr)
gsl_sf_airy_zero_Bi_deriv_e Y gsl_sf_airy_zero_Bi_deriv(uint_expr)
gsl_sf_bessel_J0_e Y gsl_sf_bessel_J0(dbl_expr)
gsl_sf_bessel_J1_e Y gsl_sf_bessel_J1(dbl_expr)
gsl_sf_bessel_Jn_e Y gsl_sf_bessel_Jn(int_expr,dbl_expr)
gsl_sf_bessel_Jn_array Y status=gsl_sf_bessel_Jn_array(int,int,double,&var_out)
gsl_sf_bessel_Y0_e Y gsl_sf_bessel_Y0(dbl_expr)
gsl_sf_bessel_Y1_e Y gsl_sf_bessel_Y1(dbl_expr)
gsl_sf_bessel_Yn_e Y gsl_sf_bessel_Yn(int_expr,dbl_expr)
gsl_sf_bessel_Yn_array Y gsl_sf_bessel_Yn_array
gsl_sf_bessel_I0_e Y gsl_sf_bessel_I0(dbl_expr)
gsl_sf_bessel_I1_e Y gsl_sf_bessel_I1(dbl_expr)
gsl_sf_bessel_In_e Y gsl_sf_bessel_In(int_expr,dbl_expr)
gsl_sf_bessel_In_array Y status=gsl_sf_bessel_In_array(int,int,double,&var_out)
gsl_sf_bessel_I0_scaled_e Y gsl_sf_bessel_I0_scaled(dbl_expr)
gsl_sf_bessel_I1_scaled_e Y gsl_sf_bessel_I1_scaled(dbl_expr)
gsl_sf_bessel_In_scaled_e Y gsl_sf_bessel_In_scaled(int_expr,dbl_expr)
gsl_sf_bessel_In_scaled_array Y staus=gsl_sf_bessel_In_scaled_array(int,int,double,&var_out)
gsl_sf_bessel_K0_e Y gsl_sf_bessel_K0(dbl_expr)
gsl_sf_bessel_K1_e Y gsl_sf_bessel_K1(dbl_expr)
gsl_sf_bessel_Kn_e Y gsl_sf_bessel_Kn(int_expr,dbl_expr)
gsl_sf_bessel_Kn_array Y status=gsl_sf_bessel_Kn_array(int,int,double,&var_out)
gsl_sf_bessel_K0_scaled_e Y gsl_sf_bessel_K0_scaled(dbl_expr)
gsl_sf_bessel_K1_scaled_e Y gsl_sf_bessel_K1_scaled(dbl_expr)
gsl_sf_bessel_Kn_scaled_e Y gsl_sf_bessel_Kn_scaled(int_expr,dbl_expr)
gsl_sf_bessel_Kn_scaled_array Y status=gsl_sf_bessel_Kn_scaled_array(int,int,double,&var_out)
gsl_sf_bessel_j0_e Y gsl_sf_bessel_J0(dbl_expr)
gsl_sf_bessel_j1_e Y gsl_sf_bessel_J1(dbl_expr)
gsl_sf_bessel_j2_e Y gsl_sf_bessel_j2(dbl_expr)
gsl_sf_bessel_jl_e Y gsl_sf_bessel_jl(int_expr,dbl_expr)
gsl_sf_bessel_jl_array Y status=gsl_sf_bessel_jl_array(int,double,&var_out)
gsl_sf_bessel_jl_steed_array Y gsl_sf_bessel_jl_steed_array
gsl_sf_bessel_y0_e Y gsl_sf_bessel_Y0(dbl_expr)
gsl_sf_bessel_y1_e Y gsl_sf_bessel_Y1(dbl_expr)
gsl_sf_bessel_y2_e Y gsl_sf_bessel_y2(dbl_expr)
gsl_sf_bessel_yl_e Y gsl_sf_bessel_yl(int_expr,dbl_expr)
gsl_sf_bessel_yl_array Y status=gsl_sf_bessel_yl_array(int,double,&var_out)
gsl_sf_bessel_i0_scaled_e Y gsl_sf_bessel_I0_scaled(dbl_expr)
gsl_sf_bessel_i1_scaled_e Y gsl_sf_bessel_I1_scaled(dbl_expr)
gsl_sf_bessel_i2_scaled_e Y gsl_sf_bessel_i2_scaled(dbl_expr)
gsl_sf_bessel_il_scaled_e Y gsl_sf_bessel_il_scaled(int_expr,dbl_expr)
gsl_sf_bessel_il_scaled_array Y status=gsl_sf_bessel_il_scaled_array(int,double,&var_out)
gsl_sf_bessel_k0_scaled_e Y gsl_sf_bessel_K0_scaled(dbl_expr)
gsl_sf_bessel_k1_scaled_e Y gsl_sf_bessel_K1_scaled(dbl_expr)
gsl_sf_bessel_k2_scaled_e Y gsl_sf_bessel_k2_scaled(dbl_expr)
gsl_sf_bessel_kl_scaled_e Y gsl_sf_bessel_kl_scaled(int_expr,dbl_expr)
gsl_sf_bessel_kl_scaled_array Y status=gsl_sf_bessel_kl_scaled_array(int,double,&var_out)
gsl_sf_bessel_Jnu_e Y gsl_sf_bessel_Jnu(dbl_expr,dbl_expr)
gsl_sf_bessel_Ynu_e Y gsl_sf_bessel_Ynu(dbl_expr,dbl_expr)
gsl_sf_bessel_sequence_Jnu_e N gsl_sf_bessel_sequence_Jnu
gsl_sf_bessel_Inu_scaled_e Y gsl_sf_bessel_Inu_scaled(dbl_expr,dbl_expr)
gsl_sf_bessel_Inu_e Y gsl_sf_bessel_Inu(dbl_expr,dbl_expr)
gsl_sf_bessel_Knu_scaled_e Y gsl_sf_bessel_Knu_scaled(dbl_expr,dbl_expr)
gsl_sf_bessel_Knu_e Y gsl_sf_bessel_Knu(dbl_expr,dbl_expr)
gsl_sf_bessel_lnKnu_e Y gsl_sf_bessel_lnKnu(dbl_expr,dbl_expr)
gsl_sf_bessel_zero_J0_e Y gsl_sf_bessel_zero_J0(uint_expr)
gsl_sf_bessel_zero_J1_e Y gsl_sf_bessel_zero_J1(uint_expr)
gsl_sf_bessel_zero_Jnu_e N gsl_sf_bessel_zero_Jnu
gsl_sf_clausen_e Y gsl_sf_clausen(dbl_expr)
gsl_sf_hydrogenicR_1_e N gsl_sf_hydrogenicR_1
gsl_sf_hydrogenicR_e N gsl_sf_hydrogenicR
gsl_sf_coulomb_wave_FG_e N gsl_sf_coulomb_wave_FG
gsl_sf_coulomb_wave_F_array N gsl_sf_coulomb_wave_F_array
gsl_sf_coulomb_wave_FG_array N gsl_sf_coulomb_wave_FG_array
gsl_sf_coulomb_wave_FGp_array N gsl_sf_coulomb_wave_FGp_array
gsl_sf_coulomb_wave_sphF_array N gsl_sf_coulomb_wave_sphF_array
gsl_sf_coulomb_CL_e N gsl_sf_coulomb_CL
gsl_sf_coulomb_CL_array N gsl_sf_coulomb_CL_array
gsl_sf_coupling_3j_e N gsl_sf_coupling_3j
gsl_sf_coupling_6j_e N gsl_sf_coupling_6j
gsl_sf_coupling_RacahW_e N gsl_sf_coupling_RacahW
gsl_sf_coupling_9j_e N gsl_sf_coupling_9j
gsl_sf_coupling_6j_INCORRECT_e N gsl_sf_coupling_6j_INCORRECT
gsl_sf_dawson_e Y gsl_sf_dawson(dbl_expr)
gsl_sf_debye_1_e Y gsl_sf_debye_1(dbl_expr)
gsl_sf_debye_2_e Y gsl_sf_debye_2(dbl_expr)
gsl_sf_debye_3_e Y gsl_sf_debye_3(dbl_expr)
gsl_sf_debye_4_e Y gsl_sf_debye_4(dbl_expr)
gsl_sf_debye_5_e Y gsl_sf_debye_5(dbl_expr)
gsl_sf_debye_6_e Y gsl_sf_debye_6(dbl_expr)
gsl_sf_dilog_e N gsl_sf_dilog
gsl_sf_complex_dilog_xy_e N gsl_sf_complex_dilog_xy_e
gsl_sf_complex_dilog_e N gsl_sf_complex_dilog
gsl_sf_complex_spence_xy_e N gsl_sf_complex_spence_xy_e
gsl_sf_multiply_e N gsl_sf_multiply
gsl_sf_multiply_err_e N gsl_sf_multiply_err
gsl_sf_ellint_Kcomp_e Y gsl_sf_ellint_Kcomp(dbl_expr)
gsl_sf_ellint_Ecomp_e Y gsl_sf_ellint_Ecomp(dbl_expr)
gsl_sf_ellint_Pcomp_e Y gsl_sf_ellint_Pcomp(dbl_expr,dbl_expr)
gsl_sf_ellint_Dcomp_e Y gsl_sf_ellint_Dcomp(dbl_expr)
gsl_sf_ellint_F_e Y gsl_sf_ellint_F(dbl_expr,dbl_expr)
gsl_sf_ellint_E_e Y gsl_sf_ellint_E(dbl_expr,dbl_expr)
gsl_sf_ellint_P_e Y gsl_sf_ellint_P(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_ellint_D_e Y gsl_sf_ellint_D(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_ellint_RC_e Y gsl_sf_ellint_RC(dbl_expr,dbl_expr)
gsl_sf_ellint_RD_e Y gsl_sf_ellint_RD(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_ellint_RF_e Y gsl_sf_ellint_RF(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_ellint_RJ_e Y gsl_sf_ellint_RJ(dbl_expr,dbl_expr,dbl_expr,dbl_expr)
gsl_sf_elljac_e N gsl_sf_elljac
gsl_sf_erfc_e Y gsl_sf_erfc(dbl_expr)
gsl_sf_log_erfc_e Y gsl_sf_log_erfc(dbl_expr)
gsl_sf_erf_e Y gsl_sf_erf(dbl_expr)
gsl_sf_erf_Z_e Y gsl_sf_erf_Z(dbl_expr)
gsl_sf_erf_Q_e Y gsl_sf_erf_Q(dbl_expr)
gsl_sf_hazard_e Y gsl_sf_hazard(dbl_expr)
gsl_sf_exp_e Y gsl_sf_exp(dbl_expr)
gsl_sf_exp_e10_e N gsl_sf_exp_e10
gsl_sf_exp_mult_e Y gsl_sf_exp_mult(dbl_expr,dbl_expr)
gsl_sf_exp_mult_e10_e N gsl_sf_exp_mult_e10
gsl_sf_expm1_e Y gsl_sf_expm1(dbl_expr)
gsl_sf_exprel_e Y gsl_sf_exprel(dbl_expr)
gsl_sf_exprel_2_e Y gsl_sf_exprel_2(dbl_expr)
gsl_sf_exprel_n_e Y gsl_sf_exprel_n(int_expr,dbl_expr)
gsl_sf_exp_err_e Y gsl_sf_exp_err(dbl_expr,dbl_expr)
gsl_sf_exp_err_e10_e N gsl_sf_exp_err_e10
gsl_sf_exp_mult_err_e N gsl_sf_exp_mult_err
gsl_sf_exp_mult_err_e10_e N gsl_sf_exp_mult_err_e10
gsl_sf_expint_E1_e Y gsl_sf_expint_E1(dbl_expr)
gsl_sf_expint_E2_e Y gsl_sf_expint_E2(dbl_expr)
gsl_sf_expint_En_e Y gsl_sf_expint_En(int_expr,dbl_expr)
gsl_sf_expint_E1_scaled_e Y gsl_sf_expint_E1_scaled(dbl_expr)
gsl_sf_expint_E2_scaled_e Y gsl_sf_expint_E2_scaled(dbl_expr)
gsl_sf_expint_En_scaled_e Y gsl_sf_expint_En_scaled(int_expr,dbl_expr)
gsl_sf_expint_Ei_e Y gsl_sf_expint_Ei(dbl_expr)
gsl_sf_expint_Ei_scaled_e Y gsl_sf_expint_Ei_scaled(dbl_expr)
gsl_sf_Shi_e Y gsl_sf_Shi(dbl_expr)
gsl_sf_Chi_e Y gsl_sf_Chi(dbl_expr)
gsl_sf_expint_3_e Y gsl_sf_expint_3(dbl_expr)
gsl_sf_Si_e Y gsl_sf_Si(dbl_expr)
gsl_sf_Ci_e Y gsl_sf_Ci(dbl_expr)
gsl_sf_atanint_e Y gsl_sf_atanint(dbl_expr)
gsl_sf_fermi_dirac_m1_e Y gsl_sf_fermi_dirac_m1(dbl_expr)
gsl_sf_fermi_dirac_0_e Y gsl_sf_fermi_dirac_0(dbl_expr)
gsl_sf_fermi_dirac_1_e Y gsl_sf_fermi_dirac_1(dbl_expr)
gsl_sf_fermi_dirac_2_e Y gsl_sf_fermi_dirac_2(dbl_expr)
gsl_sf_fermi_dirac_int_e Y gsl_sf_fermi_dirac_int(int_expr,dbl_expr)
gsl_sf_fermi_dirac_mhalf_e Y gsl_sf_fermi_dirac_mhalf(dbl_expr)
gsl_sf_fermi_dirac_half_e Y gsl_sf_fermi_dirac_half(dbl_expr)
gsl_sf_fermi_dirac_3half_e Y gsl_sf_fermi_dirac_3half(dbl_expr)
gsl_sf_fermi_dirac_inc_0_e Y gsl_sf_fermi_dirac_inc_0(dbl_expr,dbl_expr)
gsl_sf_lngamma_e Y gsl_sf_lngamma(dbl_expr)
gsl_sf_lngamma_sgn_e N gsl_sf_lngamma_sgn
gsl_sf_gamma_e Y gsl_sf_gamma(dbl_expr)
gsl_sf_gammastar_e Y gsl_sf_gammastar(dbl_expr)
gsl_sf_gammainv_e Y gsl_sf_gammainv(dbl_expr)
gsl_sf_lngamma_complex_e N gsl_sf_lngamma_complex
gsl_sf_taylorcoeff_e Y gsl_sf_taylorcoeff(int_expr,dbl_expr)
gsl_sf_fact_e Y gsl_sf_fact(uint_expr)
gsl_sf_doublefact_e Y gsl_sf_doublefact(uint_expr)
gsl_sf_lnfact_e Y gsl_sf_lnfact(uint_expr)
gsl_sf_lndoublefact_e Y gsl_sf_lndoublefact(uint_expr)
gsl_sf_lnchoose_e N gsl_sf_lnchoose
gsl_sf_choose_e N gsl_sf_choose
gsl_sf_lnpoch_e Y gsl_sf_lnpoch(dbl_expr,dbl_expr)
gsl_sf_lnpoch_sgn_e N gsl_sf_lnpoch_sgn
gsl_sf_poch_e Y gsl_sf_poch(dbl_expr,dbl_expr)
gsl_sf_pochrel_e Y gsl_sf_pochrel(dbl_expr,dbl_expr)
gsl_sf_gamma_inc_Q_e Y gsl_sf_gamma_inc_Q(dbl_expr,dbl_expr)
gsl_sf_gamma_inc_P_e Y gsl_sf_gamma_inc_P(dbl_expr,dbl_expr)
gsl_sf_gamma_inc_e Y gsl_sf_gamma_inc(dbl_expr,dbl_expr)
gsl_sf_lnbeta_e Y gsl_sf_lnbeta(dbl_expr,dbl_expr)
gsl_sf_lnbeta_sgn_e N gsl_sf_lnbeta_sgn
gsl_sf_beta_e Y gsl_sf_beta(dbl_expr,dbl_expr)
gsl_sf_beta_inc_e N gsl_sf_beta_inc
gsl_sf_gegenpoly_1_e Y gsl_sf_gegenpoly_1(dbl_expr,dbl_expr)
gsl_sf_gegenpoly_2_e Y gsl_sf_gegenpoly_2(dbl_expr,dbl_expr)
gsl_sf_gegenpoly_3_e Y gsl_sf_gegenpoly_3(dbl_expr,dbl_expr)
gsl_sf_gegenpoly_n_e N gsl_sf_gegenpoly_n
gsl_sf_gegenpoly_array Y gsl_sf_gegenpoly_array
gsl_sf_hyperg_0F1_e Y gsl_sf_hyperg_0F1(dbl_expr,dbl_expr)
gsl_sf_hyperg_1F1_int_e Y gsl_sf_hyperg_1F1_int(int_expr,int_expr,dbl_expr)
gsl_sf_hyperg_1F1_e Y gsl_sf_hyperg_1F1(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_U_int_e Y gsl_sf_hyperg_U_int(int_expr,int_expr,dbl_expr)
gsl_sf_hyperg_U_int_e10_e N gsl_sf_hyperg_U_int_e10
gsl_sf_hyperg_U_e Y gsl_sf_hyperg_U(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_U_e10_e N gsl_sf_hyperg_U_e10
gsl_sf_hyperg_2F1_e Y gsl_sf_hyperg_2F1(dbl_expr,dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_2F1_conj_e Y gsl_sf_hyperg_2F1_conj(dbl_expr,dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_2F1_renorm_e Y gsl_sf_hyperg_2F1_renorm(dbl_expr,dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_2F1_conj_renorm_e Y gsl_sf_hyperg_2F1_conj_renorm(dbl_expr,dbl_expr,dbl_expr,dbl_expr)
gsl_sf_hyperg_2F0_e Y gsl_sf_hyperg_2F0(dbl_expr,dbl_expr,dbl_expr)
gsl_sf_laguerre_1_e Y gsl_sf_laguerre_1(dbl_expr,dbl_expr)
gsl_sf_laguerre_2_e Y gsl_sf_laguerre_2(dbl_expr,dbl_expr)
gsl_sf_laguerre_3_e Y gsl_sf_laguerre_3(dbl_expr,dbl_expr)
gsl_sf_laguerre_n_e Y gsl_sf_laguerre_n(int_expr,dbl_expr,dbl_expr)
gsl_sf_lambert_W0_e Y gsl_sf_lambert_W0(dbl_expr)
gsl_sf_lambert_Wm1_e Y gsl_sf_lambert_Wm1(dbl_expr)
gsl_sf_legendre_Pl_e Y gsl_sf_legendre_Pl(int_expr,dbl_expr)
gsl_sf_legendre_Pl_array Y status=gsl_sf_legendre_Pl_array(int,double,&var_out)
gsl_sf_legendre_Pl_deriv_array N gsl_sf_legendre_Pl_deriv_array
gsl_sf_legendre_P1_e Y gsl_sf_legendre_P1(dbl_expr)
gsl_sf_legendre_P2_e Y gsl_sf_legendre_P2(dbl_expr)
gsl_sf_legendre_P3_e Y gsl_sf_legendre_P3(dbl_expr)
gsl_sf_legendre_Q0_e Y gsl_sf_legendre_Q0(dbl_expr)
gsl_sf_legendre_Q1_e Y gsl_sf_legendre_Q1(dbl_expr)
gsl_sf_legendre_Ql_e Y gsl_sf_legendre_Ql(int_expr,dbl_expr)
gsl_sf_legendre_Plm_e Y gsl_sf_legendre_Plm(int_expr,int_expr,dbl_expr)
gsl_sf_legendre_Plm_array Y status=gsl_sf_legendre_Plm_array(int,int,double,&var_out)
gsl_sf_legendre_Plm_deriv_array N gsl_sf_legendre_Plm_deriv_array
gsl_sf_legendre_sphPlm_e Y gsl_sf_legendre_sphPlm(int_expr,int_expr,dbl_expr)
gsl_sf_legendre_sphPlm_array Y status=gsl_sf_legendre_sphPlm_array(int,int,double,&var_out)
gsl_sf_legendre_sphPlm_deriv_array N gsl_sf_legendre_sphPlm_deriv_array
gsl_sf_legendre_array_size N gsl_sf_legendre_array_size
gsl_sf_conicalP_half_e Y gsl_sf_conicalP_half(dbl_expr,dbl_expr)
gsl_sf_conicalP_mhalf_e Y gsl_sf_conicalP_mhalf(dbl_expr,dbl_expr)
gsl_sf_conicalP_0_e Y gsl_sf_conicalP_0(dbl_expr,dbl_expr)
gsl_sf_conicalP_1_e Y gsl_sf_conicalP_1(dbl_expr,dbl_expr)
gsl_sf_conicalP_sph_reg_e Y gsl_sf_conicalP_sph_reg(int_expr,dbl_expr,dbl_expr)
gsl_sf_conicalP_cyl_reg_e Y gsl_sf_conicalP_cyl_reg(int_expr,dbl_expr,dbl_expr)
gsl_sf_legendre_H3d_0_e Y gsl_sf_legendre_H3d_0(dbl_expr,dbl_expr)
gsl_sf_legendre_H3d_1_e Y gsl_sf_legendre_H3d_1(dbl_expr,dbl_expr)
gsl_sf_legendre_H3d_e Y gsl_sf_legendre_H3d(int_expr,dbl_expr,dbl_expr)
gsl_sf_legendre_H3d_array N gsl_sf_legendre_H3d_array
gsl_sf_legendre_array_size N gsl_sf_legendre_array_size
gsl_sf_log_e Y gsl_sf_log(dbl_expr)
gsl_sf_log_abs_e Y gsl_sf_log_abs(dbl_expr)
gsl_sf_complex_log_e N gsl_sf_complex_log
gsl_sf_log_1plusx_e Y gsl_sf_log_1plusx(dbl_expr)
gsl_sf_log_1plusx_mx_e Y gsl_sf_log_1plusx_mx(dbl_expr)
gsl_sf_mathieu_a_array N gsl_sf_mathieu_a_array
gsl_sf_mathieu_b_array N gsl_sf_mathieu_b_array
gsl_sf_mathieu_a N gsl_sf_mathieu_a
gsl_sf_mathieu_b N gsl_sf_mathieu_b
gsl_sf_mathieu_a_coeff N gsl_sf_mathieu_a_coeff
gsl_sf_mathieu_b_coeff N gsl_sf_mathieu_b_coeff
gsl_sf_mathieu_ce N gsl_sf_mathieu_ce
gsl_sf_mathieu_se N gsl_sf_mathieu_se
gsl_sf_mathieu_ce_array N gsl_sf_mathieu_ce_array
gsl_sf_mathieu_se_array N gsl_sf_mathieu_se_array
gsl_sf_mathieu_Mc N gsl_sf_mathieu_Mc
gsl_sf_mathieu_Ms N gsl_sf_mathieu_Ms
gsl_sf_mathieu_Mc_array N gsl_sf_mathieu_Mc_array
gsl_sf_mathieu_Ms_array N gsl_sf_mathieu_Ms_array
gsl_sf_pow_int_e N gsl_sf_pow_int
gsl_sf_psi_int_e Y gsl_sf_psi_int(int_expr)
gsl_sf_psi_e Y gsl_sf_psi(dbl_expr)
gsl_sf_psi_1piy_e Y gsl_sf_psi_1piy(dbl_expr)
gsl_sf_complex_psi_e N gsl_sf_complex_psi
gsl_sf_psi_1_int_e Y gsl_sf_psi_1_int(int_expr)
gsl_sf_psi_1_e Y gsl_sf_psi_1(dbl_expr)
gsl_sf_psi_n_e Y gsl_sf_psi_n(int_expr,dbl_expr)
gsl_sf_synchrotron_1_e Y gsl_sf_synchrotron_1(dbl_expr)
gsl_sf_synchrotron_2_e Y gsl_sf_synchrotron_2(dbl_expr)
gsl_sf_transport_2_e Y gsl_sf_transport_2(dbl_expr)
gsl_sf_transport_3_e Y gsl_sf_transport_3(dbl_expr)
gsl_sf_transport_4_e Y gsl_sf_transport_4(dbl_expr)
gsl_sf_transport_5_e Y gsl_sf_transport_5(dbl_expr)
gsl_sf_sin_e N gsl_sf_sin
gsl_sf_cos_e N gsl_sf_cos
gsl_sf_hypot_e N gsl_sf_hypot
gsl_sf_complex_sin_e N gsl_sf_complex_sin
gsl_sf_complex_cos_e N gsl_sf_complex_cos
gsl_sf_complex_logsin_e N gsl_sf_complex_logsin
gsl_sf_sinc_e N gsl_sf_sinc
gsl_sf_lnsinh_e N gsl_sf_lnsinh
gsl_sf_lncosh_e N gsl_sf_lncosh
gsl_sf_polar_to_rect N gsl_sf_polar_to_rect
gsl_sf_rect_to_polar N gsl_sf_rect_to_polar
gsl_sf_sin_err_e N gsl_sf_sin_err
gsl_sf_cos_err_e N gsl_sf_cos_err
gsl_sf_angle_restrict_symm_e N gsl_sf_angle_restrict_symm
gsl_sf_angle_restrict_pos_e N gsl_sf_angle_restrict_pos
gsl_sf_angle_restrict_symm_err_e N gsl_sf_angle_restrict_symm_err
gsl_sf_angle_restrict_pos_err_e N gsl_sf_angle_restrict_pos_err
gsl_sf_zeta_int_e Y gsl_sf_zeta_int(int_expr)
gsl_sf_zeta_e Y gsl_sf_zeta(dbl_expr)
gsl_sf_zetam1_e Y gsl_sf_zetam1(dbl_expr)
gsl_sf_zetam1_int_e Y gsl_sf_zetam1_int(int_expr)
gsl_sf_hzeta_e Y gsl_sf_hzeta(dbl_expr,dbl_expr)
gsl_sf_eta_int_e Y gsl_sf_eta_int(int_expr)
gsl_sf_eta_e Y gsl_sf_eta(dbl_expr)


Footnotes

[1] These are the GSL standard function names postfixed with _e. NCO calls these functions automatically, without the NCO command having to specifically indicate the _e function suffix.