test17.cpp

#include "preprocessor_defines.dat"
#include "Class_Global.hpp"
#include "Class_Para_Tree.hpp"
#include "User_Data_Comm.hpp"
#include "User_Data_LB.hpp"

using namespace std;

// =================================================================================== //

int main(int argc, char *argv[]) {

#if NOMPI==0
    MPI::Init(argc, argv);

    {
#endif
        int iter = 0;
        Class_Para_Tree<2> pablo17a;
        Class_Para_Tree<2> pablo17b;

        int idx = 0;
        pablo17a.setBalance(idx,false);

        for (iter=1; iter<6; iter++){
            pablo17a.adaptGlobalRefine();
            pablo17b.adaptGlobalRefine();
        }

        double xc, yc;
        xc = yc = 0.5;
        double radius = 0.25;

        uint32_t nocts = pablo17a.getNumOctants();
        uint32_t nghosts = pablo17a.getNumGhosts();
        vector<double> oct_data(nocts, 0.0);
        vector<double> ghost_data(nghosts, 0.0);

        for (int i=0; i<nocts; i++){
            vector<double> center = pablo17a.getCenter(i);
            oct_data[i] = (pow((center[0]-xc),2.0)+pow((center[1]-yc),2.0));
            vector<vector<double> > nodes = pablo17a.getNodes(i);
            for (int j=0; j<global2D.nnodes; j++){
                double x = nodes[j][0];
                double y = nodes[j][1];
                if ((pow((x-xc),2.0)+pow((y-yc),2.0) <= pow(radius,2.0))){
                    if (center[0]<=xc){

                        pablo17a.setMarker(i,1);
                    }
                    else{

                        pablo17a.setMarker(i,-1);
                    }
                }
            }
        }

        iter = 0;
        pablo17a.updateConnectivity();
        pablo17b.updateConnectivity();
        pablo17a.write("Pablo17a_iter"+to_string(static_cast<unsigned long long>(iter)));
//      pablo7b.write("Pablo7b_iter"+to_string(static_cast<unsigned long long>((iter)));

        int start = 1;
        for (iter=start; iter<start+2; iter++){
            for (int i=0; i<nocts; i++){
                vector<double> center = pablo17a.getCenter(i);
                vector<vector<double> > nodes = pablo17a.getNodes(i);
                for (int j=0; j<global2D.nnodes; j++){
                    double x = nodes[j][0];
                    double y = nodes[j][1];
                    if ((pow((x-xc),2.0)+pow((y-yc),2.0) <= pow(radius,2.0))){
                        if (center[0]<xc){

                            pablo17a.setMarker(i,1);
                        }
                        else{

                            pablo17a.setMarker(i,-1);
                        }
                    }
                }
            }

            vector<double> oct_data_new;
            vector<uint32_t> mapper;
            vector<bool> isghost;
            pablo17a.adapt(true);
            nocts = pablo17a.getNumOctants();
            oct_data_new.resize(nocts, 0.0);

            for (uint32_t i=0; i<nocts; i++){
                pablo17a.getMapping(i, mapper, isghost);
                if (pablo17a.getIsNewC(i)){
                    for (int j=0; j<global2D.nchildren; j++){
                        if (isghost[j]){
                            oct_data_new[i] += ghost_data[mapper[j]]/global3D.nchildren;
                        }
                        else{
                            oct_data_new[i] += oct_data[mapper[j]]/global3D.nchildren;
                        }
                    }
                }
                else if (pablo17a.getIsNewR(i)){
                    oct_data_new[i] += oct_data[mapper[0]];
                }
                else{
                    oct_data_new[i] += oct_data[mapper[0]];
                }
            }
            oct_data = oct_data_new;
            oct_data_new.clear();

            pablo17a.updateConnectivity();
            pablo17a.writeTest("Pablo17a_iter"+to_string(static_cast<unsigned long long>(iter)), oct_data);

        }

#if NOMPI==0

        User_Data_LB<vector<double> > data_lb(oct_data,ghost_data);
        pablo17a.loadBalance(data_lb);
        pablo17b.loadBalance();

        pablo17a.updateConnectivity();
        pablo17a.writeTest("Pablo17a_iter"+to_string(static_cast<unsigned long long>(iter)), oct_data);

#endif


        vector<pair<pair<uint32_t,uint32_t>, pair<int, int> > > mapPablo;
        nocts = pablo17b.getNumOctants();
        mapPablo = pablo17b.mapPablos(pablo17a);
        vector<double> oct_data_b(nocts);
//      /**<Assign data to the new octree : process owner the first octant index in mapPablo.*/
//      for (int i=0; i<nocts; i++){
//          oct_data_b[i] = mapPablo[i].second.first;
//      }
        for (int i=0; i<nocts; i++){
            oct_data_b[i] = mapPablo[i].first.first;
        }

        pablo17b.updateConnectivity();
        pablo17b.writeTest("Pablo17b_iter"+to_string(static_cast<unsigned long long>(iter)), oct_data_b);

#if NOMPI==0
    }

    MPI::Finalize();
#endif
}