import std.stdio: printf;
import std.math: sqrt, pow;
import std.random;


struct Rating
{
	size_t user;
	size_t object;
	double value;
}


class ReputationAlgorithm
{
	this() {}
	this(ref Rating[] ratings,
	     ref double[] reputationUser,
	     ref double[] reputationObject) {}
}


class AvgWeighted : ReputationAlgorithm
{
	double[] weightSum;

	uint opCall(ref Rating[] ratings,
	            ref double[] reputationUser,
	            ref double[] reputationObject)
	{
		weightSum.length = reputationObject.length;
		weightSum[] = 0;

		reputationObject[] = 0;

		foreach(ref const(Rating) r; ratings) {
			reputationObject[r.object] += reputationUser[r.user]*r.value;
			weightSum[r.object] += reputationUser[r.user];
		}

		foreach(size_t o, ref double r; reputationObject)
			r /= weightSum[o];

		return 0;
	}

	this() {}

	this(ref Rating[] ratings,
	     ref double[] reputationUser,
	     ref double[] reputationObject)
	{
		opCall(ratings,reputationUser,reputationObject);
	}
}


final class AvgArithmetic : AvgWeighted
{
	uint opCall(ref Rating[] ratings,
	            ref double[] reputationUser,
	            ref double[] reputationObject)
	{
		reputationUser[] = 1;
		return super.opCall(ratings,reputationUser,reputationObject);
	}

	this() {}

	this(ref Rating[] ratings,
	     ref double[] reputationUser,
	     ref double[] reputationObject)
	{
		opCall(ratings,reputationUser,reputationObject);
	}
}


final class Yzlm : AvgWeighted
{
	double beta;
	double convergenceRequirement;
	double errorMin;
	size_t[] userLinks;
	double[] weightSum;
	double[] oldReputationObject;

	double objectReputationUpdate(ref Rating[] ratings,
	                              ref double[] reputationUser,
	                              ref double[] reputationObject)
	{
		double diff = 0;
		double[] tmp = oldReputationObject;
		oldReputationObject = reputationObject;
		reputationObject = tmp;
		reputationObject[] = 0;

		super.opCall(ratings,reputationUser,reputationObject);

		foreach(size_t object, ref const(double) r; reputationObject) {
			auto aux = r - oldReputationObject[object];
			diff += aux*aux;
		}

		return sqrt(diff);
	}

	void userReputationUpdate(ref Rating[] ratings,
	                          ref double[] reputationUser,
	                          ref double[] reputationObject)
	{
		reputationUser[] = 0;

		foreach(ref const(Rating) r; ratings) {
			auto aux = r.value - reputationObject[r.object];
			reputationUser[r.user] += aux*aux;
		}

		foreach(size_t user, ref double r; reputationUser) {
			if(userLinks[user]>0)
				r = pow( (r/userLinks[user]) + errorMin, -beta);
		}
	}

	uint opCall(ref Rating[] ratings,
	            ref double[] reputationUser,
	            ref double[] reputationObject)
	{
		userLinks.length = reputationUser.length;
		userLinks[] = 0;

		weightSum.length = reputationObject.length;
		oldReputationObject.length = reputationObject.length;

		foreach(ref const(Rating) r; ratings)
			userLinks[r.user]++;

		double diff;
		uint iterations=0;
		do {
			userReputationUpdate(ratings,
			                     reputationUser,
			                     reputationObject);
			diff = objectReputationUpdate(ratings,
			                              reputationUser,
			                              reputationObject);
			++iterations;
		} while(diff > convergenceRequirement);
		printf("Exited in %u iterations with diff = %g < %g\n",
		         iterations,diff,convergenceRequirement);
		return iterations;
	}

	this() {}

	this(double b,
	     double c,
	     double e)
	{
		beta = b;
		convergenceRequirement = c;
		errorMin = e;
		assert(beta>=0);
		assert(convergenceRequirement>0);
		assert(errorMin>=0);
	}

	this(ref Rating[] ratings,
	     ref double[] reputationUser,
	     ref double[] reputationObject,
	     double b,
	     double c,
	     double e)
	{
		this(b,c,e);

		opCall(ratings,
		       reputationUser,
		       reputationObject);
	}
}


void main()
{
	Rating[] ratings;
	double[] reputationObject, reputationUser;
	double[] objectQuality, userError;
	uint iterTotal = 0;
	Mt19937 rng;
	auto aa = new AvgArithmetic;
	auto yzlm = new Yzlm(0.8,1e-12,1e-36);

	rng.seed(1001);

	reputationObject.length = 1000;
	reputationUser.length = 1000;

	objectQuality.length = reputationObject.length;
	userError.length = reputationUser.length;

	ratings.length = reputationObject.length*reputationUser.length;

	for(uint i=0;i<100;++i) {
		foreach(ref double Q; objectQuality)
			Q = uniform(0.0,10.0,rng);

		foreach(ref double sigma2; userError)
			sigma2 = uniform(0.0,1.0,rng);

		size_t pos = 0;
		foreach(size_t object, ref double Q; objectQuality) {
			foreach(size_t user, ref double sigma2; userError) {
				ratings[pos++] = Rating(user,
				                        object,
				                        uniform((Q-sigma2), (Q+sigma2), rng));
			}
		}
		ratings.length = pos;

		printf("We now have %u ratings.\n",ratings.length);

		aa(ratings,reputationUser,reputationObject);
		iterTotal += yzlm(ratings,reputationUser,reputationObject);

		double deltaQ = 0;
		foreach(size_t object, ref const(double) r; reputationObject)
			deltaQ += pow( (r-objectQuality[object]), 2.0);
		deltaQ = sqrt(deltaQ/reputationObject.length);
		printf("[%u] Error in quality estimate: %g\n",i,deltaQ);
	}
	printf("Total iterations: %u\n",iterTotal);
}