// PROGRAM graph_sol
// ロジスティック写像の軌跡をグラフを使って求める

#include "TrueBASIC.h"
int main();
void initial(double *x0, double *r, int *iterate);
void draw_function(double r, int iterate);
void trajectory(double *x, double r, int iterate);
double f(double x, double r, int iterate);

int main()
{
    int iterate;
    double r, x;

    GWopen(0);

    initial(&x, &r, &iterate);
    draw_function(r, iterate);
    trajectory(&x, r, iterate);             // 中止するには任意のキーを押す

    GWquit();
    return 0;
}

void initial(double *x0, double *r, int *iterate)
{
    char Stmp1_[_LBUFF_];

    printf("control parameter r = ");
    fgets(Stmp1_, _LBUFF_, stdin);
    sscanf(Stmp1_, "%lg", r);
    printf("initial value of x = ");
    fgets(Stmp1_, _LBUFF_, stdin);
    sscanf(Stmp1_, "%lg", x0);
    printf("iterate of f(x) = ");
    fgets(Stmp1_, _LBUFF_, stdin);
    sscanf(Stmp1_, "%d", iterate);
}

void draw_function(double r, int iterate)
{
    int i, nplot;
    double delta, margin, x, y;
    char Stmp1_[_LBUFF_];

    nplot = 200;                            // 関数の計算される点の数
    delta = 1.0/nplot;
    margin = 0.1;
    GWindow((float)(-margin), (float)(-margin), (float)(1 + margin), 
            (float)(1 + margin));
    sprintf(Stmp1_, "r = %f", r);
    GWputtxt(0.0f, 1.0f, Stmp1_);
    GWline(0.0f, 0.0f, 1.0f, 1.0f);         // 対角線 y = x を描く
    GWmove2(0.0f, 1.0f);        // 軸を描く
    GWline2(0.0f, 0.0f);
    GWline2(1.0f, 0.0f);
    GWsetpen(-1, -1, -1, -1);               // ペンを離す
    GWsetpen(RED, -1, -1, -1);
    x = 0;
    for(i = 0; i <= nplot; ++i)
    {
        y = f(x, r, iterate);
        GWline2((float)(x), (float)(y));
        x += delta;
    }
}

void trajectory(double *x, double r, int iterate)
{
    int k;
    double x0, y, y0;

    y0 = 0;
    x0 = (*x);
    GWsetpen(-1, -1, -1, -1);
    GWsetpen(BLUE, -1, -1, -1);
    do
    {
        y = f(*x, r, iterate);
        GWline2((float)(x0), (float)(y0));
        GWline2((float)(x0), (float)(y));
        GWline2((float)(y), (float)(y));
        x0 = y;
        y0 = y;
        (*x) = y;
    } while(!(TBkeyinput()));
    k = TBgetkey();
}

double f(double x, double r, int iterate)
// f は再帰呼び出しにより定義されている
{
    double f_, y;

    if(iterate > 1)
    {
        y = f(x, r, iterate - 1);
        f_ = 4*r*y*(1 - y);
    }
    else
        f_ = 4*r*x*(1 - x);
    return f_;
}