[162b93]: src / double_int_to_double.cc.template  Maximize  Restore  History

Download this file

113 lines (107 with data), 3.3 kB

DEFUN_DLD(GSL_OCTAVE_NAME, args, nargout, "\
  -*- texinfo -*-\n\
@deftypefn {Loadable Function} {@var{y} =} GSL_OCTAVE_NAME (@var{x}, @var{n})\n\
@deftypefnx {Loadable Function} {[@var{y}, @var{err}] =} GSL_OCTAVE_NAME (@dots{})\n\
\n\
GSL_FUNC_DOCSTRING
\n\
@var{err} contains an estimate of the absolute error in the value @var{y}.\n\
\n\
This function is from the GNU Scientific Library,\n\
see @url{http://www.gnu.org/software/gsl/} for documentation.\n\
@end deftypefn\n\
")
{
    int i;
    dim_vector dv;
    
    gsl_set_error_handler (octave_gsl_errorhandler);
    
    if(args.length() != 2) {
    print_usage ();
    return octave_value();
    }
    if(!args(0).is_real_type() || !args(1).is_real_type()) {
        error("The arguments must be real.");
    print_usage ();
    return octave_value();
    }

    // Nice combinatorial explosion here
    NDArray x = args(0).array_value();
    NDArray n = args(1).array_value();
    if(n.length() == x.length()) {
    dv = x.dims();
        NDArray y(dv);
    int len = x.length();
    if(nargout < 2) {
        for(i = 0; i < len; i++) {
        y.xelem(i) = GSL_FUNC_NAME (x.xelem(i),
                                            static_cast<int>(n.xelem(i)));
        }
        return octave_value(y);
    } else {
        NDArray err(dv);
        gsl_sf_result result;
        octave_value_list retval;
        for(i = 0; i < len; i++) {
        GSL_FUNC_NAME_e (x.xelem(i),
                                 static_cast<int>(n.xelem(i)), &result);
        y.xelem(i) = result.val;
        err.xelem(i) = result.err;
        }
        retval(1) = octave_value(err);
        retval(0) = octave_value(y);
        return retval;
    }
    } else if(n.length() == 1) {
    dv = x.dims();
        NDArray y(dv);
    int len = x.length();
    int nint = static_cast<int>(n.xelem(0));
    if(nargout < 2) {
        for(i = 0; i < len; i++) {
        y.xelem(i) = GSL_FUNC_NAME (x.xelem(i), nint);
        }
        return octave_value(y);
    } else {
        NDArray err(dv);
        gsl_sf_result result;
        octave_value_list retval;
        for(i = 0; i < len; i++) {
        GSL_FUNC_NAME_e (x.xelem(i), nint, &result);
        y.xelem(i) = result.val;
        err.xelem(i) = result.err;
        }
        retval(1) = octave_value(err);
        retval(0) = octave_value(y);
        return retval;
    }
    } else if(x.length() == 1) {
    dv = n.dims();
        NDArray y(dv);
    int len = n.length();
    double xdouble = x.xelem(0);
    if(nargout < 2) {
        for(i = 0; i < len; i++) {
        y.xelem(i) = GSL_FUNC_NAME (xdouble,
                                            static_cast<int>(n.xelem(i)));
        }
        return octave_value(y);
    } else {
        NDArray err(dv);
        gsl_sf_result result;
        octave_value_list retval;
        for(i = 0; i < len; i++) {
        GSL_FUNC_NAME_e (xdouble,
                                 static_cast<int>(n.xelem(i)), &result);
        y.xelem(i) = result.val;
        err.xelem(i) = result.err;
        }
        retval(1) = octave_value(err);
        retval(0) = octave_value(y);
        return retval;
    }
    } else {
    error("First and second argument must either have the same size, or one of them must be scalar.");
    print_usage ();
    }

    return octave_value();

}