// Hopfield.C

#include "nsl_include.h"
#include "Hopfield.h"

int nslRand(int max)
{
	int tmp;
#ifdef NSL_PC
	tmp = rand()%max; 
#else
	tmp = nsl_rint(drand48()*(max-1)); // from 0 to 9
#endif
    return tmp;
}

HopfieldModel::HopfieldModel() : NslModel("HopfieldModel"), 
	in("in",this,10), hp("hopfield",this,10)
{
	makeConn();
}

void HopfieldModel::makeConn()
{
	nslConnect(in.out,hp.pat);
}
HopfieldInput::HopfieldInput(char* str, NslModule* parent,int size) : 
	NslModel(str,parent), out("out",this,size,size)
{
}

Hopfield::Hopfield(char* str, NslModule* parent,int size) : 
	NslModule(str,parent),
	pat("pat",this,size,size),	
	mp("mp",this,size,size),
	mf("mf",this,size,size),
	pmf("pmf",this,size,size),
	w("w",this,size,size,size,size),
	energy("energy",this),
	inOutError("inOutError",this)
{
	_size = size;
}

void Hopfield::simTrain()
{
	for (int k=0; k<_size; k++)
	  for (int l=0; l<_size; l++)
		for (int i=0; i<_size; i++)
	       for (int j=0; j<_size; j++)
			if (k!=l || i!=j) 
			   w[k][l][i][j]+=pat[k][l]*pat[i][j];
}

void Hopfield::initRun()
{
	mf = pat;
	pmf = 0;  
}
 
void Hopfield::simRun() 
{
	int i,j,k,l;

	k = nslRand(_size);
	l = nslRand(_size);
				
/*	mp[k][l] = 0;
	for (i=0; i<_size; i++)
	   for (j=0; j<_size; j++) 
		mp[k][l] += w[k][l][i][j]*mf[i][j];
*/
	mp[k][l] = (w[k][l]^mf).sum();	
	mf[k][l] = NSLstep(mp[k][l],0,-1,1);

	inOutError=0;
	for (i=0; i<_size; i++)
	   for (j=0; j<_size; j++) 	   
	 	inOutError = inOutError + fabs(pmf[i][j]-mf[i][j]);
	
//	inOutError.nslPrint();
	if (inOutError == 0) {
		cmd_out("Convergence");
		nslBreakSim();
		return;
	}

	pmf[k][l] = mf[k][l]; 

	computeEnergy();
}

void Hopfield::computeEnergy()
{
	energy = 0;
	for (int k=0; k<_size; k++)
	    for (int l=0; l<_size; l++)
		for (int i = 0; i < _size; i++)
	    	    for (int j = 0; j < _size; j++)
			energy = energy + w[k][l][i][j]*mp[k][l]*mp[i][j];
	
	energy = -0.5*energy;	
}

AslSchemaModel<HopfieldModel> _HopfieldModel("HopfieldModel");