/***************************************************************************
                          nuk.cpp  -  description
                             -------------------
    begin                : Die Feb 19 12:20:21 CET 2002
    copyright            : (C) 2002 by Sven Klauke
    email                : sklauke@wiwi.uni-bielefeld.de
    Institute of Mathematical Economics (IMW)
    University of Bielefeld, Germany
    http://www.wiwi.uni-bielefeld.de/~imw/   I
 ***************************************************************************/

/***************************************************************************
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/

// Modified by Alex Kontorovich & John F. Nash Jr. 6/11/2002

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

void calc_nuc(double *v, double min_schrittweite, double **imputation_new)
{
	int i,j;
	int weiter;
	double x_s,v_s;
	double summe=0.0;
	int angenommen;
	int erfolglos;
	int t1,t2;
	int s_hyp=0;
	double schrittweite;
	
	int versuche;
//	int kontrolle[coalitions];
	double *imputation_old;

	double *excess_old;
	double *excess_new;

	imputation_old=(double *) malloc( n * sizeof(double));
	excess_old=(double *) malloc( coalitions * sizeof(double));
	excess_new=(double *) malloc( coalitions * sizeof(double));

	for (i = 0; i < n ; i++)
	{
		if(v[coalitions-1]==1)
			imputation_old[i] = v[coalitions - 1] / n;
		else
			imputation_old[i]=(int)(v[coalitions-1]/n);
	}

	for (i=0;i<n;i++)
		(*imputation_new)[i] = imputation_old[i];

	for (i = 0; i < coalitions; i++)
	{
		x_s = 0;
		for (j = 0; j < n; j++)
			x_s += j_in_i (j, i) * imputation_old[j];
		
		v_s = v[i];
		excess_old[i] = v_s - x_s;
	}
	excess_old[0]=0;
	
	if(n<=4)
		excess_old[coalitions-1]=fabs(excess_old[coalitions-1]/min_schrittweite);
	else
		excess_old[coalitions-1]=0.0;

	// now need to sort excess in reverse

	reverse_sort(&excess_old);

	for(schrittweite=0.1;schrittweite>min_schrittweite;schrittweite/=10.0)
	{
		versuche=10000;
		weiter = 1;
		angenommen=0;
		erfolglos=0;
		t1=t2=0;
		do {
			weiter = 1;
			summe=0.0;
			if (angenommen)
			{
				for (i = 0; i < n; i++)
				{
					(*imputation_new)[i] = imputation_old[i] + (double)(j_in_i(i,t1)-j_in_i(i,t2)) 
						* schrittweite;
					//imputation_new[i] = imputation_old[i] 
					//  + (double)(j_in_i(i,t1)-j_in_i(i,t2)) * schrittweite * (1.0+double(zufall(8)));
					summe+=(*imputation_new)[i];
				}
			}
			else
			{
				t1=zufall(coalitions-1);
				t2=zufall(coalitions-1);
				s_hyp=zufall(n-1);
				summe=0.0;
//				kontrolle[t1]=kontrolle[t1]+1;
//				kontrolle[t2]=kontrolle[t2]+1;
				
				for (i = 0; i < n; i++)
				{
					(*imputation_new)[i] = imputation_old[i] + (double)(j_in_i(i,t1)-j_in_i(i,t2)) 
						* schrittweite;
					
					summe+=(*imputation_new)[i];
				}
			}

			if(n>4)
			  (*imputation_new)[s_hyp]+=(v[coalitions-1]-summe);

			for (i = 1; i < coalitions; i++)
			{
				x_s = 0;
				
				for (j = 0; j < n; j++)
					x_s += (double)(j_in_i (j, i)) * (*imputation_new)[j];
				
				v_s = v[i];
				//excess_new[i] = (double)((int)(((v_s - x_s)/min_schrittweite+0.5)*min_schrittweite));
				excess_new[i]=v_s-x_s;
				
				if(fabs(excess_new[i])<min_schrittweite) 
					excess_new[i]=0.0;
			}
			excess_new[0]=0.0;

			if(n<=4)
				excess_new[coalitions-1]=fabs(excess_new[coalitions-1])/min_schrittweite;

			reverse_sort(&excess_new);

			if (less_than(excess_new,excess_old)) // is new excess better?
			// here we have found a better excess, so start over
			{
				weiter = 1;
				angenommen = 1;
				erfolglos=0;
				
				for (i=0;i<n;i++)
					imputation_old[i] = (*imputation_new)[i];
				for (i=0;i<coalitions;i++)
					excess_old[i] = excess_new[i];
			}
			else
			// imputation was no better, so if done enough times, quit out
			{
				erfolglos++;
				if (erfolglos >= versuche)
					weiter = 0;
				else
					weiter = 1;
				for (i=0;i<n;i++)
					(*imputation_new)[i] = imputation_old[i];
				angenommen = 0;
			}
		} while(weiter);
		
		for(i=0;i<n;i++)
			printf("%1.7f ", (*imputation_new)[i]);
		printf("\n");
	}

	free(excess_new);
	free(excess_old);
	free(imputation_old);
	return;
}