/*!\file:  NodesPartitioning.cpp
 * \brief: partition elements and nodes and vertices
 */ 

#ifdef HAVE_CONFIG_H
	#include <config.h>
#else
#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
#endif

#include <string.h>
#include "../../objects/objects.h"
#include "../../shared/shared.h"
#include "../../EnumDefinitions/EnumDefinitions.h"
#include "../../io/io.h"
#include "../../include/include.h"
#include "../MeshPartitionx/MeshPartitionx.h"
#include "../ModelProcessorx/ModelProcessorx.h"

void  DiscontinuousGalerkinNodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel);
void  ContinuousGalerkinNodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel);

void  NodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel, bool continuous){
	
	/*First thing, this is a new partition for a new analysis_type, therefore, to avoid a leak, erase the nodes partition that might come through pmy_nodes: */
	xDelete<bool>(*pmy_nodes);

	/*Now, depending on whether we are running galerkin discontinous or continuous elements, carry out a different partition of the nodes: */
	if(continuous==true)
		ContinuousGalerkinNodesPartitioning(pmy_nodes,my_elements, my_vertices, iomodel);
	else
		DiscontinuousGalerkinNodesPartitioning(pmy_nodes,my_elements, my_vertices, iomodel);
}

void  ContinuousGalerkinNodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel){

	/*as many nodes as there are vertices */
	int    numberofvertices;

	/*output: */
	bool* my_nodes=NULL;

	/*Fetch parameters: */
	iomodel->Constant(&numberofvertices,MeshNumberofverticesEnum);

	my_nodes=xNew<bool>(numberofvertices);
	for(int i=0;i<numberofvertices;i++) my_nodes[i]=(bool)my_vertices[i];

	/*Assign output pointers:*/
	*pmy_nodes=my_nodes;
}


void  DiscontinuousGalerkinNodesPartitioning(bool** pmy_nodes,bool* my_elements, int* my_vertices, IoModel* iomodel){

	int    numberofelements;

	/*Fetch parameters: */
	iomodel->Constant(&numberofelements,MeshNumberofelementsEnum);

	/*each element has it own nodes (as many as vertices) + additional nodes from neighbouring elements for each edge. This yields to a very different partition for 
	 * the nodes and the vertices. The vertices are similar to continuous galerkin, but the nodes partitioning involves edges, which mess up sorting of 
	 * ids. */
	
	int i,j;
	int    dim;

	/*output: */
	bool*   my_nodes=NULL;

	int  i1,i2;
	int  cols;
	int  e1,e2;
	int  pos;
	int  numberofedges;
	int *edges         = NULL;
	int *elements      = NULL;

	/*Fetch parameters: */
	iomodel->Constant(&dim,MeshDimensionEnum);

	/*Build discontinuous node partitioning
	 *  - there are three nodes per element (discontinous)
	 *  - for each element present of each partition, its three nodes will be in this partition
	 *  - the edges require the dofs of the 2 nodes of each elements sharing the edge.
	 *    if the 2 elements sharing the edge are on 2 different cpus, we must duplicate
	 *    the two nodes that are not on the cpus so that the edge can access the dofs of
	 *    all its 4 nodes
	 */

	/*Allocate*/
	my_nodes=xNewZeroInit<bool>(3*numberofelements);

	/*First: add all the nodes of all the elements belonging to this cpu*/
	if (dim==2){
		for (i=0;i<numberofelements;i++){
			if (my_elements[i]){
				my_nodes[3*i+0]=true;
				my_nodes[3*i+1]=true;
				my_nodes[3*i+2]=true;
			}
		}
	}
	else{
		_error2_("not implemented yet");
	}

	/*Second: add all missing nodes*/

	/*Get edges and elements*/
	iomodel->FetchData(&edges,&numberofedges,&cols,MeshEdgesEnum);
	iomodel->FetchData(&elements,NULL,NULL,MeshElementsEnum);
	if (cols!=4) _error2_("field edges should have 4 columns");

	/*!All elements have been partitioned above, only create elements for this CPU: */
	for (i=0;i<numberofedges;i++){

		/*Get left and right elements*/
		e1=edges[4*i+2]-1; //edges are [node1 node2 elem1 elem2]
		e2=edges[4*i+3]-1; //edges are [node1 node2 elem1 elem2]

		/* 1) If the element e1 is in the current partition
		 * 2) and if the edge of the element is shared by another element (internal edge)
		 * 3) and if this element is not in the same partition:
		 * we must clone the nodes on this partition so that the loads (Numericalflux)
		 * will have access to their properties (dofs,...)*/
		if(my_elements[e1] && e2!=-2 && !my_elements[e2]){

			/*1: Get vertices ids*/
			i1=edges[4*i+0];
			i2=edges[4*i+1];

			/*2: Get the column where these ids are located in the index*/
			pos=UNDEF;
			for(j=0;j<3;j++){
				if (elements[3*e2+j]==i1) pos=j;
			}

			/*3: We have the id of the elements and the position of the vertices in the index
			 * we can now create the corresponding nodes:*/
			if (pos==0){
				my_nodes[e2*3+0]=true;
				my_nodes[e2*3+2]=true;
			}
			else if(pos==1){
				my_nodes[e2*3+1]=true;
				my_nodes[e2*3+0]=true;
			}
			else if (pos==2){
				my_nodes[e2*3+2]=true;
				my_nodes[e2*3+1]=true;
			}
			else{
				_error2_("Problem in edges creation");
			}
		}
	}

	/*Free data: */
	xDelete<int>(elements);
	xDelete<int>(edges);

	/*Assign output pointers:*/
	*pmy_nodes=my_nodes;
}
