/*
 * CreateNodes.c:
 */

#include "../../toolkits/toolkits.h"
#include "../../classes/classes.h"
#include "../../shared/shared.h"
#include "../MeshPartitionx/MeshPartitionx.h"
#include "./ModelProcessorx.h"

void CreateNodes(Nodes** pnodes, IoModel* iomodel,int analysis,int finite_element,int approximation){

	/*Intermediaries*/
	int   i,j,counter,lid=0;
	bool *my_edges = NULL;
	bool *my_nodes = NULL;
	Node *node     = NULL;

	/*Recover pointer: */
	Nodes* nodes=*pnodes;

	/*First create nodes*/
	if(!nodes) nodes = new Nodes();

	switch(finite_element){
		case P1Enum:
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}
			break;

		case P1DGEnum:
			NodesPartitioning(&my_nodes,iomodel->my_elements,iomodel->my_vertices,iomodel,false);
			for(i=0;i<iomodel->numberofelements;i++){
				for(j=0;j<3;j++){
					if(my_nodes[3*i+j]){ 
						nodes->AddObject(new Node(iomodel->nodecounter+3*i+j+1,iomodel->nodecounter+3*i+j,lid++,iomodel->elements[+3*i+j]-1,iomodel,analysis,approximation));

					}
				}
			}
			break;

		case P1bubbleEnum:
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}
			for(i=0;i<iomodel->numberofelements;i++){
				if(iomodel->my_elements[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,approximation));
				}
			}
			break;

		case P1bubblecondensedEnum:
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}
			for(i=0;i<iomodel->numberofelements;i++){
				if(iomodel->my_elements[i]){
					node = new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,approximation);
					node->Deactivate();
					nodes->AddObject(node);
				}
			}
			break;

		case P1xP2Enum:
			EdgesPartitioning(&my_edges,iomodel);
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}

			counter = iomodel->numberofvertices;
			for(i=0;i<iomodel->numberofedges;i++){
				if(iomodel->edges[i*3+2]==2){
					if(my_edges[i]){
						node = new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,counter+1,lid++,0,iomodel,analysis,approximation);
						nodes->AddObject(node);
					}
					counter++;
				}
			}
			break;
		case P2xP1Enum:
			EdgesPartitioning(&my_edges,iomodel);
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}

			counter = iomodel->numberofvertices;
			for(i=0;i<iomodel->numberofedges;i++){
				if(iomodel->edges[i*3+2]!=2){
					if(my_edges[i]){
						node = new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,counter+1,lid++,0,iomodel,analysis,approximation);
						nodes->AddObject(node);
					}
					counter++;
				}
			}
			break;

		case P2Enum:
			EdgesPartitioning(&my_edges,iomodel);
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,approximation));
				}
			}
			for(i=0;i<iomodel->numberofedges;i++){
				if(my_edges[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,approximation));
				}
			}
			break;

		/*Stokes elements*/
		case P1P1Enum:
			_assert_(approximation==FSApproximationEnum);
			/*P1 velocity*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,FSvelocityEnum));
				}
			}
			/*P1 pressure*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,i,iomodel,analysis,FSpressureEnum));
				}
			}
			break;
		case P1P1GLSEnum:
			_assert_(approximation==FSApproximationEnum);
			/*P1 velocity*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,FSvelocityEnum));
				}
			}
			/*P1 pressure*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,i,iomodel,analysis,FSpressureEnum));
				}
			}
			break;
		case MINIcondensedEnum:
			_assert_(approximation==FSApproximationEnum);
			/*P1+ velocity (bubble statically condensed)*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,FSvelocityEnum));
				}
			}
			for(i=0;i<iomodel->numberofelements;i++){
				if(iomodel->my_elements[i]){
					node = new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,FSvelocityEnum);
					node->Deactivate();
					nodes->AddObject(node);
				}
			}
			/*P1 pressure*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+iomodel->numberofelements+i+1,iomodel->numberofvertices+iomodel->numberofelements+i,lid++,i,iomodel,analysis,FSpressureEnum));
				}
			}
			break;
		case MINIEnum:
			_assert_(approximation==FSApproximationEnum);
			/*P1+ velocity*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,FSvelocityEnum));
				}
			}
			for(i=0;i<iomodel->numberofelements;i++){
				if(iomodel->my_elements[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,FSvelocityEnum));
				}
			}
			/*P1 pressure*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+iomodel->numberofelements+i+1,iomodel->numberofvertices+iomodel->numberofelements+i,lid++,i,iomodel,analysis,FSpressureEnum));
				}
			}
			break;
		case TaylorHoodEnum:
			_assert_(approximation==FSApproximationEnum);
			/*P2 velocity*/
			EdgesPartitioning(&my_edges,iomodel);
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+i+1,i,lid++,i,iomodel,analysis,FSvelocityEnum));
				}
			}
			for(i=0;i<iomodel->numberofedges;i++){
				if(my_edges[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+i+1,iomodel->numberofvertices+i,lid++,0,iomodel,analysis,FSvelocityEnum));
				}
			}

			/*P1 pressure*/
			for(i=0;i<iomodel->numberofvertices;i++){
				if(iomodel->my_vertices[i]){
					nodes->AddObject(new Node(iomodel->nodecounter+iomodel->numberofvertices+iomodel->numberofedges+i+1,iomodel->numberofvertices+iomodel->numberofedges+i,lid++,i,iomodel,analysis,FSpressureEnum));
				}
			}
			break;

		default:
			_error_("Finite element "<<EnumToStringx(finite_element)<<" not supported yet");
	}

	/*Clean up*/
	xDelete<bool>(my_edges);
	xDelete<bool>(my_nodes);

	/*Assign output pointer: */
	*pnodes=nodes;
}
