/*
 *                  THE APPROACH TO EQUILIBRIUM
 */

#include"lib351.h"

void main(){

    // output initialization
	ofstream dataOutStr;
    outFileMaker("m357p1c.out", &dataOutStr);
    
    // program variable definitions
    int nCoins, nFlips, nStep, nRanCoin, nFlip, nDot, nRun, nTails, nStatus;
	double dblRanNum, dblNormedRanNum, dFracTails;
	const double m_pi = 3.141592654;
    const double dblSeed = 13.;
    const double dblConst = 899.;
    const double dblMod = 32768.;
    const double dblMaxNum = dblMod - 1.;
    const int nDots = 71;
    const int nLines = 50;
	
    char cLine[72];
    int nStatuses[1001];

	// initialize the random number generator
    dblRanNum = dblSeed;

		// cycle through the runs
	nRun = 0;
	bool goout = false;
	while (!goout){

    	nRun++;
		
		// data input
		cout <<"Number of Coins ="<< endl;
		cin >> nCoins;
		cout <<"\nNumber of Flips ="<< endl;
		cin >> nFlips;


    	// set the step in number of turns from one line to the next
    	nStep = int(nFlips/nLines); if(nStep == 0){nStep = 1;}

    	// initialize the status of each coin to 'up'
    	for (int nCoin = 1; nCoin <= nCoins; nCoin++){
	 		nStatuses[nCoin] = 1;
    	}
    	dFracTails = 0;

		// output the header for this run
		dataOutStr<<"\f"<<endl;
		dataOutStr<<"THE APPROACH TO EQUILIBRIUM"<<endl;
    	dataOutStr<<"nCoins = "<<nCoins<<";   nFlips = "<<nFlips<<endl<<endl;
    	dataOutStr<<"flip  fraction"<<endl<<" no.  tails"<<endl;
		// output top chart boundary
		dataOutStr<<"            ";
		for (nDot = 0; nDot <= nDots; nDot++) {
			dataOutStr<<".";
		}
		dataOutStr<<endl;
        // load the first chart line with blanks
        for(nDot = 1; nDot <= nDots-1; nDot++){
        	cLine[nDot] = ' ';
        }
		// set the non-blank line elements
		cLine[0] = '*';
		cLine[50] = '|';
		cLine[nDots] = '.';
		// output the first chart line
		dataOutStr<<setiosflags(ios::fixed)<<setw(4)<<0<<setprecision(4)<<setw(7)<<0.0<<" ";
		for (nDot = 0; nDot <= nDots; nDot++) {
			dataOutStr<<cLine[nDot];
		}
		dataOutStr<<endl;

		// flipping
	 	nTails = 0;

    	for (nFlip=1; nFlip<=nFlips; nFlip++) {
    		
    		// get the next random number in (0,1)
    		dblRanNum = fmod(dblConst*dblRanNum,dblMod);
	    	dblNormedRanNum = dblRanNum/dblMaxNum;

        	// determine the coin to be flipped
        	nRanCoin = int(dblNormedRanNum*(nCoins - 1.0) + 1.0);
        	
        	// if it is 'heads', make it 'tails'
        	nStatus = nStatuses[nRanCoin];
			if (nStatus == 1) {
        		nStatuses[nRanCoin] = 0;
        		nTails++;
        	}
        	// if it is 'tails', make it 'heads'
        	if (nStatus == 0) {
        		nStatuses[nRanCoin] = 1;
        		nTails--;
        	}
        	
        	dFracTails = double(nTails)/double(nCoins);

        	// if this is an output line, do it
        	if ((nFlip/nStep)*nStep == nFlip) {
        	
        		// load the chart line with blanks
        		for(nDot = 1; nDot <= nDots-1; nDot++){
        		    cLine[nDot] = ' ';
        		}
		
				// set the non-blank line elements
				cLine[0] = '.';
				cLine[50] = '|';
				cLine[nDots] = '.';
				cLine[int(dFracTails*100.)] = '*';
		
				// output the line
				dataOutStr<<setiosflags(ios::fixed)<<setw(4)<<nFlip<<setprecision(4)<<setw(7)<<dFracTails<<" ";
				for (nDot = 0; nDot <= nDots; nDot++) {
					dataOutStr<<cLine[nDot];
				}
				dataOutStr<<endl;
			}
    	}
		// output bottom chart boundary
		dataOutStr<<"            ";
		for (nDot = 0; nDot <= nDots; nDot++) {
			dataOutStr<<".";
		}
		dataOutStr<<endl;
		
		// another run?
		goout = query();
	}
    // close the output file
	fileClose("m357p1c.out", &dataOutStr);
};