output error "gsl: bisection.c:55: ERROR: function value is not finite" does not match with the line 55 in bisection.c file

Question

我想使用 GSL 二分例程求函数 f(x)=1/(2sin(x)^2)+log(sin(x))-k 的根，k 常数，在区间 [0;Pi/2].

     int main(void *params){
  
     struct func_params *part= (struct func_params*)params;
    
     int status;
     int iter = 0, max_iter = 10;
     const gsl_root_fsolver_type *R;
     gsl_root_fsolver *s;
     double x_lo = 1e-4, x_hi = M_PI/2.;

     gsl_function F;
     F.function = &my_func;
     F.params = &params;
     R = gsl_root_fsolver_bisection;
     s = gsl_root_fsolver_alloc(R);
     gsl_root_fsolver_set(s,&F,x_lo,x_hi);

     printf ("using %s method\n",gsl_root_fsolver_name(s));
     printf ("%5s [%9s, %9s] %9s %10s %9s\n","iter", "lower", "upper", "root","err", "err(est)");

    return(0);

}

This code returns a gsl: bisection.c:55: ERROR: function value is not finite Default GSL error handler invoked. 问题是，当我查看 bisection.c 文件时，第 55 行与之前的错误不对应......知道是什么发生了什么？

Answer 1

您必须提供一个区间，使函数只跨越 y = 0 一次。

例如，[0.1, 2.0]这里：

#include <stdio.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_roots.h>

double my_func(double x,  void* unused)
{
    double k = 1;
    return 1 / (2 * pow(sin(x), 2)) + log(sin(x)) - k;
}

int main()
{
    gsl_function F;
    F.function = &my_func;
    F.params = 0;

    gsl_root_fsolver *s = gsl_root_fsolver_alloc(gsl_root_fsolver_bisection);

    double x_lo = +0.1, x_hi = +2.0;
    gsl_root_fsolver_set(s, &F, x_lo, x_hi);

    for(int i = 0; i < 1000; i++) {

        int status = gsl_root_fsolver_iterate(s);

        double left_int  = gsl_root_fsolver_x_lower(s);
        double right_int = gsl_root_fsolver_x_upper(s);

        printf("iteration %03d: [%.04lf, %.04lf]\n", i, left_int, right_int);

        status = gsl_root_test_interval(left_int, right_int, 1.0e-5, 1.0e-15);
        if(status != GSL_CONTINUE) {

            printf("status: %s\n", gsl_strerror(status));
            printf("\nRoot interval = [%.010lf, %.010lf]\n", left_int, right_int);
            break;
        }
    }
}

输出：

$ gcc main.c -lgsl -lm && ./a.out
iteration 000: [0.1000, 1.0500]
iteration 001: [0.5750, 1.0500]
iteration 002: [0.5750, 0.8125]
iteration 003: [0.5750, 0.6938]
iteration 004: [0.5750, 0.6344]
iteration 005: [0.5750, 0.6047]
iteration 006: [0.5898, 0.6047]
iteration 007: [0.5973, 0.6047]
iteration 008: [0.5973, 0.6010]
iteration 009: [0.5973, 0.5991]
iteration 010: [0.5982, 0.5991]
iteration 011: [0.5987, 0.5991]
iteration 012: [0.5989, 0.5991]
iteration 013: [0.5989, 0.5990]
iteration 014: [0.5989, 0.5990]
iteration 015: [0.5989, 0.5990]
iteration 016: [0.5989, 0.5990]
iteration 017: [0.5989, 0.5990]
status: success

Root interval = [0.5989471436, 0.5989543915]