//-------------------------------------------------------------------
//	minassn.cpp
//
//	This program implements a variant of M.L. Balinski's primal 
//	algorithm for the mimimum-cost personnel assignment problem
//	for the number of jobs and workers specified by the user as
//	a command-line argument and using randomly generated costs.
//
//		to compile:  $ g++ minassn.cpp -o minassn
//		to execute:  $ ./ minassn <n>
//
//	Reference: 
//	  M. L. Balinski and R. E. Gomory, "A primal method for the 
//	assignment and transportation problems," Management Science
//	volume 10 (1964), pp. 578-593.
//
//	programmer: ALLAN CRUSE
//	date begun: 02 MAY 2011
//	completion: 08 MAY 2011
//-------------------------------------------------------------------

#include <stdio.h>	// for printf(), perror() 
#include <fcntl.h>	// for open() 
#include <stdlib.h>	// for exit() 
#include <unistd.h>	// for read(), write(), close() 
#include <time.h>	// for time()

#define LIMIT 	16	// maximum number of workers and jobs

int	l, k, M = 6;
int	cost[ LIMIT ][ LIMIT ];
int	r[ LIMIT ][ LIMIT ];
int	s[ LIMIT ][ LIMIT ];
int	assign[ LIMIT ];
int	rowlabel[ LIMIT ];
int	collabel[ LIMIT ];
int	u[ LIMIT ], v[ LIMIT ];
int	epsilon;

void inputdata( void )
{
	srand( time( NULL ) );

	for (int i = 0; i < M; i++)
	for (int j = 0; j < M; j++)
		cost[ i ][ j ] = rand() % 100 ;
}

void show_data( void )
{
	printf( "\n\n                 " ); 
	printf( "MINIMUM-COST ASSIGNMENT PROBLEM: " );
	printf( "%d WORKERS \n", M );
	printf( "\nCOST-MATRIX: \n" );
	for (int i = 0; i < M; i++)
		{
		printf( "\n      " );
		for (int j = 0; j < M; j++)
			printf( " %8d ", cost[i][j] );
		}
	printf( "\n\n" );
}

void initialize( void )
{
	// initialize 'assign[] to the diagonal-assignment 
	for (int j = 0; j < M; j++) assign[ j ] = j;

	// initialize 'v[]' to the current assignment costs
	for (int j = 0; j < M; j++) v[ j ] = cost[j][ assign[j] ];
	 
	// initialize 'u[]' to zeros
	for (int i = 0; i < M; i++) u[ i ] = 0; 

	// initialize rowlabel[] and collabel[]
	for (int i = 0; i < M; i++) rowlabel[i] = -1;
	for (int j = 0; j < M; j++) collabel[j] = -1;
}

void computeslack( void )
{
	for (int i = 0; i < M; i++)
	for (int j = 0; j < M; j++)
		s[i][j] = cost[i][j] - u[i] - v[j];
	epsilon = 0;
}

void labelpass( void )
{
	// reinitialize the label-vectors to -1 
	for (int i = 0; i < M; i++) rowlabel[i] = -1;
	for (int j = 0; j < M; j++) collabel[j] = -1;

	// do no labeling unless k and l are within bounds
	if (( k == M )||( l == M )) return;

	// convention: k = row-number, l = column-number
	// start by labeling column l with row-number k
	collabel[l] = k;

	// whenever a column gets labeled with a row-number, 
	// then label the row which is currently assigned to 
	// that column with this column-number
	rowlabel[ assign[ l ] ] = l;

	// now see if any further labeling is possible
	int	moretodo = 1;
	while ( moretodo )
		{
		moretodo = 0;
		//
		// whenever row i  has received a label, and
		// row i is different from  row k, if there 
		// are unlabeled columns whose slack-values 
		// are zeros in that row, then label those
		// columns with the row-number i.
		for (int i = 0; i < M; i++)
		    if (( rowlabel[ i ] != -1 )&&( i != k ))
		        for (int j = 0; j < M; j++)
			    if ( collabel[ j ] == -1 )
				{
				if ( s[i][j] == 0 )
				    {
				    collabel[j] = i;
				    rowlabel[ assign[j] ] = j;
				    moretodo = 1;
				    }
				}
		}


	// epsilon should be the minimum of the values s[i][j]
	// over all cells [i,j] such that row i is labeled and
	// column j is unlabeled and s[i][j] is nonnegative,  
	// if any such entries exist; and should be -s[k][l] 
	// if no such entries exist or in case row k receivd
	// a label.
	epsilon = -s[k][l];	// this will be positive

	// return if epsilon is positive and row k received a label
	if (( epsilon > 0 )&&( rowlabel[k] != -1 )) return;

	// otherwise maybe modify epsilon if row k remains unlabeled
	for (int i = 0; i < M; i++) 
		{
		if ( rowlabel[i] == -1 ) continue;
		// look at cells in labeled row 
		for (int j = 0; j < M; j++) 
			{
			if ( collabel[j] != -1 ) continue;
			// look at cells in unlabeled column
			if (( s[i][j] >= 0 )&&( s[i][j] < epsilon ))
					epsilon = s[i][j];
			}
		}

	if (( epsilon == 0 )||( rowlabel[k] != -1 )) 
			epsilon = -s[k][l];
}

void showtable( void )
{
	printf( "\n" );
	printf( "     V" );
	for (int j = 0; j < M; j++) printf( " %8d ", v[j] );
	printf( "\n" );
	printf( "   U " );
	printf( "\n" );
	for (int i = 0; i < M; i++)
		{
		printf( "%4d  ", u[i] );
		for (int j = 0; j < M; j++)
			{
			if ( i == assign[j] ) 
				printf( "       [%d]", s[i][j] );
			else 	{ 
			
				printf( " %8d", s[i][j] );
				if (( i == k )&&( j == l )) printf( "*" );
				else	printf( " " );
				}
			}

		if ( rowlabel[i] != -1 ) 
			printf( "  { %d }", 1 + rowlabel[i] );
		printf( "\n" );
		}
	printf( "\n" );
	printf( "         " );
	for (int j = 0; j < M; j++) 
		{
		if ( collabel[j] != -1 ) 	
			printf( "  { %d }   ", 1 + collabel[j] );
		else	printf( "          " );
		}
	printf( "\n" );	

	printf( "\n" );
	printf( "epsilon = %d ", epsilon );
	printf( "\n" );

	int	total = 0;
	for (int j = 0; j < M; j++)
		{
		int	i = assign[j];
		total += cost[i][j]; 
		}
	
	int	dualsum = 0;
	for (int i = 0; i < M; i++) dualsum += u[i] + v[i];

	printf( "\n" );
	printf( "current assignment: " );
	for (int i = 0; i < M; i++) printf( " %d ", 1 + assign[i] );
	printf( "\n" );
	printf( "current cost: %d \n", total );
	printf( "sum of duals = %d \n", dualsum );
	printf( "\n" );
}

void revise( void )
{
	// if row i is labeled, increate u[i] by epsilon
	for (int i = 0; i < M; i++) 
		if ( rowlabel[i] != -1 ) u[i] += epsilon;
	
	// if column j is labeled, decrease v[j] by epsilon
	for (int j = 0; j < M; j++) 
		if ( collabel[j] != -1 ) v[j] -= epsilon;

	printf( "REVISED DUAL VARIABLES\n" );
	computeslack();
}

void change( void )
{

	// show the current assignment:
	printf( "OLD ASSIGNMENT: " );
	for (int j = 0; j < M; j++) printf( " %2d ", 1+assign[j] );
	printf( "\n" );


	// change 'assign[]' entries in the unique labeled circuit
	int	i = k;
	int	j = l;
	printf( "CYCLE: (%d,%d) ", 1+k, 1+l );

	do	{
		j = rowlabel[i];
		printf( "--> (%d,%d) ", 1+i, 1+j );

		i = collabel[j];
		printf( "--> (%d,%d) ", 1+i, 1+j );

		// change row-number assigned to column j
		assign[j] = i;
		}
	while (( i != k )&&( j != l ));
	printf( "\n" );
	printf( "CHANGED PRIMAL VARIABLES \n" );

	// show the changed assignment:
	printf( "NEW ASSIGNMENT: " );
	for (int j = 0; j < M; j++) printf( " %2d ", 1+assign[j] );
	printf( "\n" );

	
	// reduce the value of v[l] by epsilon
	v[l] -= epsilon;
	printf( "REVISED DUAL VARIABLES \n" );

	computeslack();
	printf( "RECOMPUTED SLACK VALUES\n" );
}

void showanswer( void )
{
	printf( "The minimum-cost assignment has been obtained \n" );
	printf( "\n" );
}

int main( int argc, char **argv )
{
	// allow a user to override the default problem-size
	if ( argc > 1 ) M = atoi( argv[ 1 ] );
	if (( M < 2 )||( M >= LIMIT )) 
		{ fprintf( stderr, "argument out-of-range\n" ); exit(1); }

	inputdata();
	show_data();
	initialize();
	computeslack();

	for (l = 0; l < M; l++)
	for (k = 0; k < M; k++)
		{
		while ( s[k][l] < 0 )
			{
			labelpass();
			showtable();
			if ( rowlabel[k] == -1 ) revise();
			else	change();
			}
		}

	labelpass();
	showtable();
	showanswer();
}