c - 黄金分割搜索

Question

我试图找到函数 zt 中 tl=0 和 tu=10 之间的最大值，但代码无法正常工作找到最大值。我的代码有什么问题吗？

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

#define e 2.718281828459

double zt(double t)
{
    return (100 + (80 / 15)* (55 + 80*9.81 / 15) * (1 - pow(e, -(15/80)*t)) - (80*9.81 / 15)*t);
}



void Golden(double tl, double tu)
{
    int i;                              //이터레이션 i
    int imax = 100;                 //i가 무한히 반복되지 않도록 max값 설정
    double es = 0.01;   //상대오차가 충분히 작아졌을때 중단하기 위한 es 설정
    long double t1 = tl + (tu - tl)*((pow(5,0.5) - 1) / 2);
    long double t2 = tu - (tu - tl)*((pow(5,0.5) - 1) / 2);
    long double pz;    
    long double pzold = 0;
    double pt;
    long double ea;                     // 상대오차

    printf("Golden Section Search \n");

    for (i = 1; i <= imax; i++)
    {

        if (zt(t1) > zt(t2))
        {
            pz = zt(t1);
            pt = t1;
        }

        else if (zt(t2)>zt(t1))
        {
            pz = zt(t2);
            pt = t2;
        }

        ea = fabs((tu - tl) / pt) * 100; //

        printf("iteration=%d Peak altitude=%.10f Peak Time=%.10f Error=%.15f\n", i, pz, pt, ea);

        pzold = pz;                 //xrold 값을 xr 로 설정하여 다음 이터레이션의 상대오차를 구하기위해 정의한다

        if (ea < es)
        {
            break;              //ea가 es보다 작아질 경우 xr이 근에 충분히 근접하였으므로 중단
        }


        if (zt(t1)>zt(t2))          //조건문으로 함수 f(xr)의 음수 양수 판단후 다음 iteration 에서의 xu, xl 값 설정
        {
            tl = t2;
            t2 = t1;
            t1= tl + fabs(tu - tl)*((pow(5,0.5) - 1) / 2);
        }  

        else if (zt(t2)>zt(t1))
        {
            tu = t1;
            t1 = t2;
            t2= tu - fabs(tu - tl)*((pow(5,0.5) - 1) / 2);
        }


    }

}

int main(void)
{
    Golden(0, 10);

    return 0;
}

Answer 1

我怀疑错误的结果是由正在进行的整数计算引起的。让我们看一下函数zt:

double zt(double t)
{
    return (100 + (80 / 15)* (55 + 80*9.81 / 15) * (1 - pow(e, -(15/80)*t)) - (80*9.81 / 15)*t);
}

请注意，大多数数字都是整数文字(即没有小数部分)。这意味着整数运算会发生，直到涉及 double 值为止。

使用除法运算符进行浮点算术和整数算术之间的一个显着区别/:

int val_1 = 10 / 4; // 2, because the fractional part is discarded
double val_2 = 10.0 / 4.0; // 2.5

要解决代码中的此问题，您可以将 return 语句更改为:

return (100.0 + (80.0 / 15.0) * (55.0 + 80.0*9.81 / 15.0)) * (1.0 - pow(e, -(15.0/80.0)*t)) - (80.0*9.81 / 15.0);

另一个可能的问题是一些圆括号放错了位置，导致结果评估不同。通过使用变量来保存每个值得注意的子表达式，然后用它们进行计算，可以更好地检测到这一点(如果存在)。