2014-04-13 48 views
0

我想通过Golden Section方法找到Gamma函数的最小点。但是当我执行程序时,我得到了分段错误错误。我认为,因为我是新手C用户,问题可能是由于调用功能Min_Search_Golden_Section错误。这是我的完整代码。我找不到我的错误。提前致谢。黄金分割例程分割错误

#include <stdio.h> 
#include <stdlib.h> 
#include <math.h> 
#include <float.h> 

#ifndef M_PI 
#define M_PI 3.14159265358979323846 
#endif 

#define A 12 
#define sqrt5 2.236067977499789696 

static int Stopping_Rule(double x0, double x1, double tolerance); 


double sp_gamma(double z) 
{ 
    const int a = A; 
    static double c_space[A]; 
    static double *c = NULL; 
    int k; 
    double accm; 

    if (c == NULL) { 
    double k1_factrl = 1.0; /* (k - 1)!*(-1)^k with 0!==1*/ 
    c = c_space; 
    c[0] = sqrt(2.0*M_PI); 
    for(k=1; k < a; k++) { 
     c[k] = exp(a-k) * pow(a-k, k-0.5)/k1_factrl; 
     k1_factrl *= -k; 
    } 
    } 
    accm = c[0]; 
    for(k=1; k < a; k++) { 
    accm += c[k]/(z + k); 
    } 
    accm *= exp(-(z+a)) * pow(z+a, z+0.5); /* Gamma(z+1) */ 
    return accm/z; 
} 

void Min_Search_Golden_Section(double (*f)(double), double* a, double *fa, 
            double* b, double* fb, double tolerance) 
{ 
    static const double lambda = 0.5 * (sqrt5 - 1.0); 
    static const double mu = 0.5 * (3.0 - sqrt5);   // = 1 - lambda 
    double x1; 
    double x2; 
    double fx1; 
    double fx2; 


       // Find first two internal points and evaluate 
       // the function at the two internal points. 

    x1 = *b - lambda * (*b - *a);        
    x2 = *a + lambda * (*b - *a);       
    fx1 = f(x1); 
    fx2 = f(x2); 

      // Verify that the tolerance is an acceptable number 

    if (tolerance <= 0.0) tolerance = sqrt(DBL_EPSILON) * (*b - *a); 

      // Loop by exluding segments from current endpoints a, b 
      // to current internal points x1, x2 and then calculating 
      // a new internal point until the length of the interval 
      // is less than or equal to the tolerance. 

    while (! Stopping_Rule(*a, *b, tolerance)) { 
     if (fx1 > fx2) { 
     *a = x1; 
     *fa = fx1; 
     if (Stopping_Rule(*a, *b, tolerance)) break; 
     x1 = x2; 
     fx1 = fx2; 
     x2 = *b - mu * (*b - *a); 
     fx2 = f(x2); 
     } else { 
     *b = x2; 
     *fb = fx2; 
     if (Stopping_Rule(*a, *b, tolerance)) break; 
     x2 = x1; 
     fx2 = fx1; 
     x1 = *a + mu * (*b - *a); 
     fx1 = f(x1); 
     } 
    } 
    return; 
} 


int main() 
{ 
    double x; 
    double a = 0.0, b = 4.0, fa = 0.00001, fb = 6.0; 
    double fx = sp_gamma(x); 
    Min_Search_Golden_Section(&fx, &a, &fa, &b, &fb, 0.0000001); 
    return 0; 
} 


static int Stopping_Rule(double x0, double x1, double tolerance) 
{ 
    double xm = 0.5 * fabs(x1 + x0); 

    if (xm <= 1.0) return (fabs(x1 - x0) < tolerance) ? 1 : 0; 
    return (fabs(x1 - x0) < tolerance * xm) ? 1 : 0; 
} 

回答

3

您应该会收到编译器错误。 Min_Search_Golden_Section的第一个参数应该是一个函数指针,但是您需要传递一个变量的地址。

当你遇到编译器错误时,修复它们 - 不要运行程序和希望。 :)

我猜你只是为了写:

Min_Search_Golden_Section(&sp_gamma, &a, &fa, &b, &fb, 0.0000001); 
+0

它可能是一个警告,而不是一个错误,在这种情况下,OP应该把他的警告级别。 – ooga

+0

它工作!我认为在使用指针时应该更加小心。谢谢马特:) – ChipotleHS