Golden section search

Question

I am trying to find the maximum value in function zt between tl=0 and tu=10 but the code is not working to find the maximum value. What is wrong with my code?

#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

My suspicion is that the incorrect results are caused by the integer calculations going on. Let's look at the function 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);
}

Note that most of the numbers there are integer literals (that is, there is no fractional part). This means that integer operations occur until a double value gets involved.

One notable difference between floating-point arithmetic and integer arithmetic in with the division operator /:

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

To fix this issue in your code, you could change the return statement to:

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);

Another possible issue is that some round brackets are misplaced, causing the results to evaluate to something different. This can be better detected (if it is present) by using variables to hold each notable sub-expression, and then doing calculations with those.