/*
 * Calculator by method of brute force Monte-Carlo
 *
 * Problem proposed by Mato Nagel (2010) regarding the consequence of
 * a "Dunning / Kruger Effect" on elections in a Democracy. This program
 * by Telford Tendys, Copyright (C) 2012, gives a significantly different result
 * to the conclusion of Mato Nagel, suggesting that one of us is wrong.
 *
 * This program may be re-distributed without modification.
 *
 * See also:
 *
 *          A Mathematical Model of Democratic Elections
 *                           Mato Nagel
 * Center for Nephrology and Metabolic Disorders, A.-Schweitzer-Ring 32,
 *                   D-02943 Weisswasser, Germany
 *
 */

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

#define NUM_PEOPLE     1000000
#define NUM_ELECTIONS  1000000
#define REFIL_CQ       1000

double cq[ NUM_PEOPLE ];
int vote[ NUM_PEOPLE ];

/*
 * Just take one value for convenience, inefficient but never mind.
 * I'm very impressed with my quad-core "AMD Phenom(tm) II", especially the price!
 *
 * See also:
 *       http://en.wikipedia.org/wiki/Normal_distribution#Generating_values_from_normal_distribution
 */
double Box_Muller( void )
{
	for(;;)
	{
		double bm;

		double u = random();
		double v = random();
		u /= RAND_MAX;
		v /= RAND_MAX;

		bm = sqrt( -2 * log( u )) * cos( 2 * M_PI * v );
		/*
		 * Note that standard normal has mean of 0 and SD of 1,
		 * but Psychologist normal is mean 100, and SD of 15.
		 * Also for no particular reason, reject the possibility
		 * of anything being below 0 nor greater than 200.
		 */
		bm *= 15;
		bm += 100;
		if( bm >= 0 && bm <= 200 ) { return( bm ); }
	}
}

int cq_compar( const void *p1, const void *p2 )
{
	register const double *dp1;
	register const double *dp2;

	dp1 = p1;
	dp2 = p2;
	if( *dp1 > *dp2 ) { return( 1 ); }
	if( *dp1 < *dp2 ) { return( 0 ); }
	/* Not expecting any NaN, Inf, etc */
	return( 0 );
}

/*
 * Very simple array fill with normal random numbers.
 * Then sort the array from lowest to highest.
 */
int fill_cq()
{
	int i;

	/*
	 * Put a bunch of normally distributed random numbers together
	 */
	for( i = 0; i < NUM_PEOPLE; i++ )
	{
		cq[ i ] = Box_Muller();
	}

	/*
 	 * Sort the stack (standard library sort function)
	 */
	qsort( cq, NUM_PEOPLE, sizeof( cq[ 0 ]), &cq_compar );
	return( 0 );
}

/*
 * Note that vote is an integer -- we only accept whole votes.
 *
 * Thus, we need many trials to get a suitable spread of winners,
 * based on uniform random voting. We presume the random() function
 * is good to go on any modulus (people are allowed to vote for
 * themselves, and the top guy ALWAYS votes for himself). Actually,
 * I'm not sure GNU's random() can really hold up to this type of
 * usage but my results are so massively different to Nagel's
 * calculation that really the fine details of the random generator
 * don't freak me out.
 */
int fill_vote( void )
{
	int i;
	double x = 0;

	/* Empty out the ballot box */
	bzero( vote, sizeof( vote[ 0 ]) * NUM_PEOPLE );

	for( i = 0; i < ( NUM_PEOPLE - 1 ); i++ )
	{
		int j = random();
		/* Lowest CQ casts completely random vote */
		j %= ( NUM_PEOPLE - i );
		j += i;
		/* Put a ballot against the candidate */
		vote[ j ]++;
	}
	vote[ NUM_PEOPLE - 1 ]++; /* Always votes for self */
	return( 0 );
}

/*
 * The winner is the person with the most votes.
 * Yup, that's how Democracy works.
 *
 * However, the broad spread of this system makes a tie somewhat likely.
 * We could do a run-off election, that's a heap of extra code for not
 * much value, just print the two extreme candidates and plot both.
 *
 * The difference between tied candidates turns out to not be all that big.
 */
int declare_winner( FILE *out )
{
	int i;
	int best1 = 0;
	int best2 = 0;

	for( i = 1; i < NUM_PEOPLE; i++ )
	{
		/* best1 favours low CQ in case of a tie-breaker */
		if( vote[ i ] >  vote[ best1 ]) { best1 = i; }
		/* best2 favours high CQ in case of tie-breaker */
		if( vote[ i ] >= vote[ best2 ]) { best2 = i; }
	}
	fprintf( out, "%f\t%f\n", cq[ best1 ], cq[ best2 ]);
	return( 0 );
}

/*
 * ====== COMPILED AGAINST THE FOLLOWING ======
 *
 * GNU C Library stable release version 2.12.2, by Roland McGrath et al.
 * Copyright (C) 2010 Free Software Foundation, Inc.
 * This is free software; see the source for copying conditions.
 * There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A
 * PARTICULAR PURPOSE.
 * Compiled by GNU CC version 4.5.2.
 * Compiled on a Linux 2.6.38 system on 2011-07-08.
 * Available extensions:
 *         C stubs add-on version 2.1.2
 *         crypt add-on version 2.1 by Michael Glad and others
 *         Gentoo patchset 3
 *         GNU Libidn by Simon Josefsson
 *         Native POSIX Threads Library by Ulrich Drepper et al
 *         Support for some architectures added on, not maintained in glibc core.
 *         BIND-8.2.3-T5B
 * libc ABIs: UNIQUE IFUNC
 * For bug reporting instructions, please see:
 * <http://www.gnu.org/software/libc/bugs.html>.
 *
 * ------------------------------------------------------
 */

int main( void )
{
	int i, refil;
	FILE *fp;

	srand( 123456 );
	fill_cq();

	/* Dump out a basic set of CQ numbers to show population */
	fp = fopen( "population_cq.dat", "w" );
	fprintf( fp, "CQ\n" );
	for( i = 0; i < NUM_PEOPLE; i++ )
	{
		fprintf( fp, "%f\n", cq[ i ]);
	}
	fclose( fp );

	/* Output the results of multiple elections */
	fp = fopen( "winner_cq.dat", "w" );
	fprintf( fp, "WIN1\tWIN2\n" );
	refil = 0;
	for( i = 0; i < NUM_ELECTIONS; i++ )
	{
		fill_vote();
		declare_winner( fp );
		/* Every so often re-generate CQ numbers */
		if( refil++ >= REFIL_CQ )
		{
			fill_cq();
			refil = 0;
		}
		fflush( fp );
	}
	fclose( fp );
	return( 0 );
}