/*
   This is an example of a user MetaModel. It calculates function values
   using the average of the known points, weighted with the inverse of the
   geometric distance to each point.
 
   This code is heavily documented to serve as an example on how to write
   user MetaModels.

   Please note that this algorithm may not be numerically stable, or
   even useful for actual optimization. It is here to serve as an
   example, and shouldn't be used in a production environmet.


   Livermore Software Technology Corporation (LSTC) holds the copyright
   for this code. LSTC hereby grants anyone the right to create and
   distribute derivative work based on this code under their own terms.
   
   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
   CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES
   INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
   MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
   DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR
   CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
   INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
   (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
   GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
   BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
   LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
   OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
   POSSIBILITY OF SUCH DAMAGE.

*/

#include "UserMetaModel.h"

#include <stdio.h>
#include <fstream>
using namespace std;


PyObject *pName, *pModule, *pFunc;
PyObject *pFuncBuild, *pFuncVal;
int built = 0;

/* Our initialization routine, called by LS-OPT. */
void UserMetaModel_Init(UserMetaModel* mm)
{
	/* This is a requirement of the protocol! */
	mm->size_user = sizeof(UserMetaModel);

	/* Present ourself to LS-OPT */
	mm->description = "Interface for Python-metamodels";
	/* Here, we can allocate memory that isn't really needed for persistant
     state.
     
	This memory must be freed in  UserMetaModel_CleanUp.
	*/
}

void UserMetaModel_CleanUp(UserMetaModel* mm)
{
  /* If we would have allocated memory in UserMetaModel_Init, here is where
     we would free it */
	if (mm->userstate) free(mm->userstate);
		printf("Cleaned");

	/*Py_XDECREF(pFuncVal);
	Py_XDECREF(pFuncBuild);
	Py_DECREF(pModule);
	Py_Finalize();*/
}

int UserMetaModel_GetConstantCount(UserMetaModel* mm, int numpoints_good)
{
	/* As constants, we actually just use our points! */
	return (mm->numvar + 1) * numpoints_good;
}

void UserMetaModel_Build(UserMetaModel* mm, double ** x_case_good,
			 double * funcvals_good,int numpoints_good)
{
	int i,j;
	built=1;
	/*Conversion of the datatypes and arrangement as arguments for Python*/
	InitPython(mm);
	PyObject *coo,*coos, *vals, *pArgs;

	pArgs = PyTuple_New(2);
	coos = PyList_New(numpoints_good);
	vals = PyList_New(numpoints_good);
	
	//Arguments are saved into lists
	for (i=0; i<numpoints_good;i++) {
		coo = PyList_New(mm->numvar);
		for (j=0; j<mm->numvar;j++) {
			PyList_SET_ITEM(coo,j,PyFloat_FromDouble(x_case_good[i][j]));
		}
		PyList_SET_ITEM(coos,i,coo);
		PyList_SET_ITEM(vals,i,PyFloat_FromDouble(funcvals_good[i]));
	}
	PyTuple_SetItem(pArgs, 0, coos);
	PyTuple_SetItem(pArgs, 1, vals);

	double result = RunPython(mm,pArgs,pFuncBuild);

	//Saving the data
	int c=0;
	for (i=0; i<numpoints_good; i++) {
		for (j=0; j<mm->numvar; j++)
			mm->constants[c++] = x_case_good[i][j];
		mm->constants[c++] = funcvals_good[i];
	}
}

double UserMetaModel_FuncVal(UserMetaModel* mm, double *coords)
{
	if (built==0) printf("Metamodel has not been built - no calculation of points possible.");

	int i;
	int numvar = mm->numvar;
	int numpoints_good = mm->num_constants / (mm->numvar+1);
	PyObject *coo, *pArgs;

	//Collecting and converting the given coordinates as arguments in lists for Python
	pArgs = PyTuple_New(1);
	coo = PyList_New(numvar);
	for (i=0; i<numvar; i++) 
		PyList_SET_ITEM(coo,i,PyFloat_FromDouble(coords[i]));
	PyTuple_SetItem(pArgs, 0, coo);

	return RunPython(mm,pArgs,pFuncVal);
}

double UserMetaModel_FuncGrad(UserMetaModel* mm, int var, double *coords)
{
	double *mycoords;
	if (!mm->userstate)
		mm->userstate = malloc(sizeof(double) * mm->numvar);
	mycoords = (double*)mm->userstate;
	mm->error = UMM_WARNING;
	mm->error_message = mm->stringparam;
	/* Quick and dirty numerical Eulerian differentiator */ 
	memcpy(mycoords,coords,sizeof(double) * mm->numvar);
	mycoords[var] += 0.0001;
	return (UserMetaModel_FuncVal(mm,mycoords) - UserMetaModel_FuncVal(mm,coords)) / 0.0001;
}

void InitPython(UserMetaModel* mm){
	//Used to load the Python script "pyInterface.py" and it's contained functions "build()" and "funcval()"
	Py_Initialize();
	pName = PyString_FromString("pyInterface");
	pModule = PyImport_Import(pName);
	Py_DECREF(pName);
	if (pModule != NULL) {
		pFuncBuild = PyObject_GetAttrString(pModule, "build");
		if (pFuncBuild == NULL) {
			PyErr_Print();
			fprintf(stderr, "Failed to load function \"build()\"\n");
			mm->error = UMM_ERROR;
			mm->error_message = "Failed to load function \"build()\"";
		}
		pFuncVal = PyObject_GetAttrString(pModule, "funcval");
		if (pFuncVal == NULL) {
			PyErr_Print();
			fprintf(stderr, "Failed to load function \"funcval()\"\n");
			mm->error = UMM_ERROR;
			mm->error_message = "Failed to load function \"funcval()\"";
		}
	}
	else {
		PyErr_Print();
		fprintf(stderr, "Failed to load module \"pyInterface.py\"\n");
		mm->error = UMM_ERROR;
		mm->error_message = "Failed to load module \"pyInterface.py\"";
	}
}
double RunPython(UserMetaModel* mm,PyObject *pArgs,PyObject *pFunc)
{
	PyObject *pValue;
	if (pModule != NULL) {
		if (pFunc && PyCallable_Check(pFunc)) {
			pValue = PyObject_CallObject(pFunc, pArgs);
			Py_DECREF(pArgs);
			if (pValue == NULL) {
				PyErr_Print();
				fprintf(stderr,"Call failed\n");
				mm->error = UMM_ERROR;
				mm->error_message = "Call has failed, no result has been given. Critical Error.";
				return 0;
			}
		}
		else {
			if (PyErr_Occurred())
				PyErr_Print();
			fprintf(stderr, "Cannot find function\n");
			mm->error = UMM_ERROR;
			mm->error_message = "Python function is not loaded. Critical error.";
			return 0;
		}
	}
    else {
		mm->error = UMM_ERROR;
		mm->error_message = "Python module is not loaded. Critical error.";
		return 0;
    }

	double value;
	if ((pValue != NULL) && (pFunc != pFuncBuild)) {
		if (PyFloat_Check(pValue)) {
			value = PyFloat_AsDouble(pValue);
		}else if(PyLong_Check(pValue)) {
			value = PyLong_AsDouble(pValue);
			mm->error = UMM_WARNING;
			mm->error_message = "Result has been given as LONG. Please return float.";
			if (value = -1.0) {
				mm->error = UMM_WARNING;
				mm->error_message = "PyLong_AsDouble returned -1.0. May be datatype conversion error. Please return float.";
			}
		}else if(PyInt_Check(pValue)) {
			value = (double) PyInt_AsLong(pValue);
			mm->error = UMM_WARNING;
			mm->error_message = "Result has been given as INT. Please return float.";
		}else{
			mm->error = UMM_ERROR;
			mm->error_message = "Result has not been given as numeric. Please return float.";
			value=0;
		}
		Py_DECREF(pValue);
	} else {
		value=0;
	}
	return value;
}