#include <iostream>
#include <armadillo>
#include "dG_parallel.h"
#include "dG_library.h"

using namespace arma;

namespace dG_global
{
  int n_partitions;
  urowvec E_partitions; // hned pri nacteni ze souboru odecist jednicku (pro pripad, ze nebude master procesor, se budou procesory cislovat od 0)
  field<urowvec> partitions;
  
  field<uvec> mapM_local, vmapM_local, mapP_local, vmapP_local;
  field<uvec> mapB_local, vmapB_local;
  field<uvec> mapM_inner_local, vmapM_inner_local, mapP_inner_local, vmapP_inner_local, mapM_interfaces_local, vmapM_interfaces_local;
  field<uvec> neighbouring_partitions_IDs_local, neighbouring_partitions_IDs_unique_local;
  field<uvec> permutation_vector_local;
  
  uvec mapM_inner, vmapM_inner, mapP_inner, vmapP_inner, mapM_interfaces, vmapM_interfaces;
  uvec neighbouring_partitions_IDs, neighbouring_partitions_IDs_unique;
  uvec permutation_vector;
  
  uvec GlobalToLocal_map(uword global_nr, urowvec &vector_of_elems)
  {
    uword elem_global = floor(global_nr / (N_faces * N_fp));
	uvec elem_local = find(vector_of_elems == elem_global);
	uvec local_nr = zeros<uvec>(1);
	if(!elem_local.is_empty())
	{
	  uword elem_nr = global_nr - elem_global * N_faces * N_fp; // ID uzlu
	  local_nr(0) = elem_local(0) * N_faces * N_fp + elem_nr;
	}
	if(elem_local.is_empty())
	{
	  local_nr.reset();
	}
	return local_nr;
  }
  
  uvec GlobalToLocal_vmap(uword global_nr, urowvec &vector_of_elems)
  {
	uword elem_global = floor(global_nr / N_p);
	uvec elem_local = find(vector_of_elems == elem_global);
	uvec local_nr = zeros<uvec>(1);
	if(!elem_local.is_empty())
	{
	  uword elem_nr = global_nr - elem_global * N_p; // ID uzlu
	  local_nr(0) = elem_local(0) * N_p + elem_nr;
	}
	if(elem_local.is_empty())
	{
	  local_nr.reset();
	}
	return local_nr;
  }
  
  uvec LocalToGlobal_map(uword local_nr, urowvec &vector_of_elems)
  {
    uword elem_local = floor(local_nr / (N_faces * N_fp));
	uword elem_global = vector_of_elems(elem_local);
	uword elem_nr = local_nr - elem_local * N_faces * N_fp;
	uword global_nr = elem_global * N_faces * N_fp + elem_nr;
	uword neigh_face = floor(elem_nr / N_fp);
	uvec localtoglobal_temp = zeros<uvec>(3,1);
	localtoglobal_temp(0) = global_nr; localtoglobal_temp(1) = elem_global; localtoglobal_temp(2) = neigh_face;
	return localtoglobal_temp;
  }
  
  void ParallelStartUp()
  {
    int debug_of_maps = 0;
    int debug_of_lengths = 0;
    int debug_unique_IDs = 0;
    
    filename = directoryname + "E_partitions_for_h=" + mesh_size_parameter_string + ".txt";
    E_partitions.load(filename);
    n_partitions = max(E_partitions);
    E_partitions = E_partitions - 1; // kvuli cislovani od 0 v C++
  
    umat partitions_matrix = zeros<umat>(n_partitions,K);
	uvec pointers = zeros<uvec>(n_partitions,1);
	uword part_ID;
	for(int i = 0; i < K; i++)
	{
      part_ID = E_partitions(i);
      partitions_matrix(part_ID,pointers(part_ID)) = i;
      pointers(part_ID) = pointers(part_ID) + 1;
	}
	partitions = field<urowvec>(n_partitions,1);
	urowvec partitions_vector;
	for(int i = 0; i < n_partitions; i++)
	{
      partitions_vector = partitions_matrix.row(i);
      partitions_vector.resize(pointers(i));
      partitions(i,0) = partitions_vector;
	}
	if(debug_of_maps){cout << partitions(0,0) << endl;}
	
	uword interface_node_ID = K * N_faces * N_fp;
	// mapM, vmapM, mapP, vmapP - v globalni siti
	mapM_local = field<uvec>(n_partitions,1);
	vmapM_local = field<uvec>(n_partitions,1);
	mapP_local = field<uvec>(n_partitions,1);
	vmapP_local = field<uvec>(n_partitions,1);
	uvec mapM_local_temp, vmapM_local_temp, mapP_local_temp, vmapP_local_temp;
	int K_local, n_local;
	urowvec elems_in_partition;
	uvec indeces_local, indeces_global;
	for(int i = 0; i < n_partitions; i++)
	{
	  elems_in_partition = partitions(i,0);
      K_local = elems_in_partition.n_cols;
      n_local = K_local * N_faces * N_fp;
      mapM_local_temp = zeros<uvec>(n_local,1);
      vmapM_local_temp = zeros<uvec>(n_local,1);
      mapP_local_temp = zeros<uvec>(n_local,1);
      vmapP_local_temp = zeros<uvec>(n_local,1);
      for(int j = 0; j < K_local; j++)
      {
        indeces_local = linspace<uvec>(j * N_faces * N_fp,(j + 1) * N_faces * N_fp - 1,N_faces * N_fp);
        indeces_global = linspace<uvec>(elems_in_partition(j) * N_faces * N_fp,(elems_in_partition(j) + 1) * N_faces * N_fp - 1,N_faces * N_fp);
        mapM_local_temp(indeces_local) = mapM(indeces_global);
        vmapM_local_temp(indeces_local) = vmapM(indeces_global);
        mapP_local_temp(indeces_local) = mapP(indeces_global);
        vmapP_local_temp(indeces_local) = vmapP(indeces_global);
      }
      mapM_local(i,0) = mapM_local_temp;
      vmapM_local(i,0) = vmapM_local_temp;
      mapP_local(i,0) = mapP_local_temp;
      vmapP_local(i,0) = vmapP_local_temp;
	}
	if(debug_of_maps){cout << "mapM_local(0,0) = \n" << mapM_local(0,0) << endl;}
	if(debug_of_maps){cout << "mapP_local(0,0) = \n" << mapP_local(0,0) << endl;}
	// mapM, vmapM, mapP, vmapP - v lokalni siti
	uvec mapM_partition, vmapM_partition, mapP_partition, vmapP_partition;
	uword mapM_global_ID, vmapM_global_ID, mapP_global_ID, vmapP_global_ID;
	uvec mapM_local_ID, vmapM_local_ID, mapP_local_ID, vmapP_local_ID;
	for(int i = 0; i < n_partitions; i++)
	{
	  elems_in_partition = partitions(i,0);
	  mapM_partition = mapM_local(i,0);
      vmapM_partition = vmapM_local(i,0);
      mapP_partition = mapP_local(i,0);
      vmapP_partition = vmapP_local(i,0);
      n_local = mapM_partition.n_rows;
      mapM_local_temp = zeros<uvec>(n_local,1);
      vmapM_local_temp = zeros<uvec>(n_local,1);
      mapP_local_temp = zeros<uvec>(n_local,1);
      vmapP_local_temp = zeros<uvec>(n_local,1);
      for(int j = 0; j < n_local; j++)
      {
        mapM_global_ID = mapM_partition(j);
        mapM_local_ID = GlobalToLocal_map(mapM_global_ID,elems_in_partition);
        vmapM_global_ID = vmapM_partition(j);
        vmapM_local_ID = GlobalToLocal_vmap(vmapM_global_ID,elems_in_partition);
        mapP_global_ID = mapP_partition(j);
        mapP_local_ID = GlobalToLocal_map(mapP_global_ID,elems_in_partition);
        vmapP_global_ID = vmapP_partition(j);
        vmapP_local_ID = GlobalToLocal_vmap(vmapP_global_ID,elems_in_partition);
        
        mapM_local_temp(j) = mapM_local_ID(0);
        vmapM_local_temp(j) = vmapM_local_ID(0);
        if(mapP_local_ID.is_empty())
        {
          mapP_local_temp(j) = interface_node_ID;
        }
        if(!mapP_local_ID.is_empty())
        {
          mapP_local_temp(j) = mapP_local_ID(0);
        }
        if(vmapP_local_ID.is_empty())
        {
          vmapP_local_temp(j) = interface_node_ID;
        }
        if(!vmapP_local_ID.is_empty())
        {
          vmapP_local_temp(j) = vmapP_local_ID(0);
        }
      }
      mapM_local(i,0) = mapM_local_temp;
      vmapM_local(i,0) = vmapM_local_temp;
      mapP_local(i,0) = mapP_local_temp;
      vmapP_local(i,0) = vmapP_local_temp;
	}
	if(debug_of_lengths){cout << "mapM_local(0,0) ma " << mapM_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "mapM_local(0,0) = \n" << mapM_local(0,0) << endl;}
	if(debug_of_lengths){cout << "mapP_local(0,0) ma " << mapP_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "mapP_local(0,0) = \n" << mapP_local(0,0) << endl;}
	
	// urceni bodu, ktere lezi na fyzikalni hranici
	mapB_local = field<uvec>(n_partitions,1);
	vmapB_local = field<uvec>(n_partitions,1);
	for(int i = 0; i < n_partitions; i++)
	{
      mapB_local(i,0) = mapM_local(i,0)(find(mapM_local(i,0) == mapP_local(i,0)));
      vmapB_local(i,0) = vmapM_local(i,0)(find(mapM_local(i,0) == mapP_local(i,0)));
	}
	if(debug_of_lengths){cout << "mapB_local(0,0) ma " << mapB_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "maB_local(0,0) = \n" << mapB_local(0,0) << endl;}

	// urceni bodu, ktere lezi uvnitr dane partition (nesousedi s jinou partition)
	mapM_inner_local = field<uvec>(n_partitions,1);
	vmapM_inner_local = field<uvec>(n_partitions,1);
	mapP_inner_local = field<uvec>(n_partitions,1);
	vmapP_inner_local = field<uvec>(n_partitions,1);
	for(int i = 0; i < n_partitions; i++)
	{
      mapM_inner_local(i,0) = mapM_local(i,0)(find(mapP_local(i,0) != interface_node_ID));
      vmapM_inner_local(i,0) = vmapM_local(i,0)(find(mapP_local(i,0) != interface_node_ID));
      mapP_inner_local(i,0) = mapP_local(i,0)(find(mapP_local(i,0) != interface_node_ID));
      vmapP_inner_local(i,0) = vmapP_local(i,0)(find(mapP_local(i,0) != interface_node_ID));
	}
	if(debug_of_lengths){cout << "mapM_inner_local(0,0) ma " << mapM_inner_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "mapM_inner_local(0,0) = \n" << mapM_inner_local(0,0) << endl;}
	if(debug_of_lengths){cout << "mapP_inner_local(0,0) ma " << mapP_inner_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "mapP_inner_local(0,0) = \n" << mapP_inner_local(0,0) << endl;}

	// urceni bodu na hranicich mezi castmi site
	mapM_interfaces_local = field<uvec>(n_partitions,1);
	vmapM_interfaces_local = field<uvec>(n_partitions,1);
	for(int i = 0; i < n_partitions; i++)
	{
      mapM_interfaces_local(i,0) = mapM_local(i,0)(find(mapP_local(i,0) == interface_node_ID));
      vmapM_interfaces_local(i,0) = vmapM_local(i,0)(find(mapP_local(i,0) == interface_node_ID));
	}
	if(debug_of_lengths){cout << "mapM_interfaces_local(0,0) ma " << mapM_interfaces_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "mapM_interfaces_local(0,0) = \n" << mapM_interfaces_local(0,0) << endl;}

	// urceni casti site, ve ktere lezi sousedni body k bodum z mapM_interfaces_local
	neighbouring_partitions_IDs_local = field<uvec>(n_partitions,1);
	uvec mapM_interfaces_partition;
	uvec neigh_parts_temp;
	uword mapM_interface_local_ID;
	uvec localtoglobal;
	uword neigh_elem_global, neigh_part_temp;
	for(int i = 0; i < n_partitions; i++)
	{
      elems_in_partition = partitions(i,0);
      mapM_interfaces_partition = mapM_interfaces_local(i,0);
      n_local = mapM_interfaces_partition.n_rows;
      neigh_parts_temp = zeros<uvec>(n_local,1);
      for(int j = 0; j < n_local; j++)
      {     
        mapM_interface_local_ID = mapM_interfaces_partition(j);
        localtoglobal = LocalToGlobal_map(mapM_interface_local_ID,elems_in_partition);
        
        neigh_elem_global = EtoE(localtoglobal(1),localtoglobal(2));
        neigh_part_temp = E_partitions(neigh_elem_global);
       
        neigh_parts_temp(j) = neigh_part_temp;
      }
      neighbouring_partitions_IDs_local(i,0) = neigh_parts_temp;
	}
	if(debug_of_lengths){cout << "neighbouring_partitions_IDs_local(0,0) ma " << neighbouring_partitions_IDs_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "neighbouring_partitions_IDs_local(0,0) = \n" << neighbouring_partitions_IDs_local(0,0) << endl;}
	neighbouring_partitions_IDs_unique_local = field<uvec>(n_partitions,1);
	for(int i = 0; i < n_partitions; i++)
	{
      neighbouring_partitions_IDs_unique_local(i,0) = unique(neighbouring_partitions_IDs_local(i,0));
	}
	if(debug_unique_IDs){cout << "neighbouring_partitions_IDs_unique_local(0,0) ma " << neighbouring_partitions_IDs_unique_local(0,0).n_rows << " prvku" << endl;}
	if(debug_unique_IDs){cout << "neighbouring_partitions_IDs_unique_local(0,0) = \n" << neighbouring_partitions_IDs_unique_local(0,0) << endl;}
	
	// napocitani permutacnich vektoru
	permutation_vector_local = field<uvec>(n_partitions,1);
	uvec mapM_interfaces_global;
	uword mapM_local_ID_temp;
	uvec mapM_interface_global_IDs, global_neighbours;
	uvec current_neighbouring_partitions, current_neighbouring_partitions_unique;
	int n_neighbouring_partitions;
	uword partition_ID;
	urowvec elems_in_neighbouring_partition;
	uvec indeces, indeces_temp;
	int n2_local;
	uvec mapM_interfaces_neighbouring_partition, mapM_interfaces_neighbouring_partition_global;
	uvec permutation_vector_temp;
	uword temp_ID;
	uvec position;
	for(int i = 0; i < n_partitions; i++)
	{
      elems_in_partition = partitions(i,0);
      mapM_interfaces_partition = mapM_interfaces_local(i,0);
      n_local = mapM_interfaces_partition.n_rows;
      mapM_interfaces_global = zeros<uvec>(n_local,1);
      for(int j = 0; j < n_local; j++)
      {
        mapM_local_ID_temp = mapM_interfaces_partition(j);
        localtoglobal = LocalToGlobal_map(mapM_local_ID_temp,elems_in_partition);
        mapM_interfaces_global(j) = localtoglobal(0);
      }
      mapM_interface_global_IDs = mapP(mapM_interfaces_global);
      global_neighbours = zeros<uvec>(n_local,1);
      current_neighbouring_partitions = neighbouring_partitions_IDs_local(i,0);
      current_neighbouring_partitions_unique = neighbouring_partitions_IDs_unique_local(i,0);
      n_neighbouring_partitions = current_neighbouring_partitions_unique.n_rows;
      //cout << "kod dospel az zde" << endl;
      for(int j = 0; j < n_neighbouring_partitions; j++)
      {
        partition_ID = current_neighbouring_partitions_unique(j);
        elems_in_neighbouring_partition = partitions(partition_ID,0);
        indeces = find(current_neighbouring_partitions == partition_ID);
        mapM_interfaces_neighbouring_partition = mapM_interfaces_local(partition_ID,0)(find(neighbouring_partitions_IDs_local(partition_ID,0) == i));
        n2_local = mapM_interfaces_neighbouring_partition.n_rows;
        mapM_interfaces_neighbouring_partition_global = zeros<uvec>(n2_local,1);
        for(int k = 0; k < n2_local; k++)
        {
            mapM_local_ID_temp = mapM_interfaces_neighbouring_partition(k);
            localtoglobal = LocalToGlobal_map(mapM_local_ID_temp,elems_in_neighbouring_partition);
            mapM_interfaces_neighbouring_partition_global(k) = localtoglobal(0);
        }
        global_neighbours(indeces) = mapM_interfaces_neighbouring_partition_global;
      }
      permutation_vector_temp = zeros<uvec>(n_local,1);
      for(int j = 0; j < n_local; j++)
      {
        temp_ID = mapM_interface_global_IDs(j);
        position = find(global_neighbours == temp_ID);
        permutation_vector_temp(j) = position(0);
      }
      permutation_vector_local(i,0) = permutation_vector_temp;
	}
	if(debug_of_lengths){cout << "permutation_vector_local(0,0) ma " << permutation_vector_local(0,0).n_rows << " prvku" << endl;}
	if(debug_of_maps){cout << "permutation_vector_local(0,0) = \n" << permutation_vector_local(0,0) << endl;}
  }
}
