PROGRAM bifurcate
! plot values of x for different values of r
CALL initial(x,r,rmax,nvalues,dr,ntransient,nplot)
FOR ir = 0 to nvalues
    CALL output(x,r,ntransient,nplot)
    LET r = r + dr
NEXT ir
! maximum value of r done separately to avoid r > 1
CALL output(x,rmax,ntransient,nplot)
END

SUB initial(x0,r0,rmax,nvalues,dr,ntransient,nplot)
    INPUT prompt "initial value of control parameter r = ": r0
    ! important that r not be greater than 1
    INPUT prompt "maximum value of r = ": rmax
    ! suggest dr <= 0.01
    INPUT prompt "incremental change of r = ": dr
    INPUT prompt "number of iterations not plotted = ": ntransient
    INPUT prompt "number of iterations plotted = ": nplot
    LET nvalues = (rmax - r0)/dr  ! number of r values plotted
    LET nvalues = int(nvalues)
    LET x0 = 0.5                  ! initial value
    CLEAR
    LET xmax = 1                  ! maximum value of x
    LET mx = 0.05*xmax            ! margin
    SET WINDOW r0-dr,rmax+dr,-mx,xmax + mx
    BOX LINES r0,rmax,0,1
END SUB

SUB output(x,r,ntransient,nplot)
    DECLARE DEF f
    SET COLOR "black/white"
    SET CURSOR 1,1
    PRINT "          ";           ! erase previous output 
    SET CURSOR 1,1
    PRINT "r ="; r
    FOR i = 1 to ntransient       ! x values not plotted
        LET x = f(x,r)
    NEXT i
    SET COLOR "red"
    FOR i = 1 to 0.5*nplot
        LET x = f(x,r)
        ! show different x-values for given value of r
        PLOT r,x
    NEXT i
    ! change color to see if values of x have converged
    SET COLOR "blue"
    FOR i = (0.5*nplot + 1) to nplot
        LET x = f(x,r)
        PLOT r,x
    NEXT i
END SUB

DEF f(x,r) = 4*r*x*(1 - x)