2014-04-13 67 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