PABLO_example_00001.cpp

/*---------------------------------------------------------------------------*\
 *
 *  bitpit
 *
 *  Copyright (C) 2015-2021 OPTIMAD engineering Srl
 *
 *  -------------------------------------------------------------------------
 *  License
 *  This file is part of bitpit.
 *
 *  bitpit is free software: you can redistribute it and/or modify it
 *  under the terms of the GNU Lesser General Public License v3 (LGPL)
 *  as published by the Free Software Foundation.
 *
 *  bitpit is distributed in the hope that it will be useful, but WITHOUT
 *  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 *  FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
 *  License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with bitpit. If not, see <http://www.gnu.org/licenses/>.
 *
\*---------------------------------------------------------------------------*/
 
#if BITPIT_ENABLE_MPI==1
#include <mpi.h>
#endif
 
#include "bitpit_PABLO.hpp"
 
using namespace std;
using namespace bitpit;
 
// =================================================================================== //
// =================================================================================== //
 
void run()
{
    PabloUniform pablo1(2);
 
    pablo1.computeConnectivity();
    pablo1.write("pablo00001_iter0");
 
    pablo1.adaptGlobalRefine();
    pablo1.updateConnectivity();
    pablo1.write("pablo00001_iter1");
 
    double xc, yc;
    xc = yc = 0.5;
 
    pablo1.setBalanceCodimension(1);
 
    uint32_t nocts = pablo1.getNumOctants();
    for (unsigned int i=0; i<nocts; i++){
        array<double,3> center = pablo1.getCenter(i);
        double x = center[0];
        if (x>xc)
            pablo1.setBalance(i,false);
    }
 
    double radius = 0.4;
 
    int nref1 = 6;
    for (int iter=0; iter<nref1; iter++){
        nocts = pablo1.getNumOctants();
        for (unsigned int i=0; i<nocts; i++){
            Octant *oct = pablo1.getOctant(i);
            array<double,3> center = pablo1.getCenter(oct);
            double x = center[0];
            double y = center[1];
 
            if ((pow((x-xc),2.0)+pow((y-yc),2.0) < pow(radius,2.0)) &&
                    (y<yc)){
                pablo1.setMarker(oct, 1);
            }
        }
        pablo1.adapt();
        pablo1.updateConnectivity();
        pablo1.write("pablo00001_iter"+to_string(static_cast<unsigned long long>(iter+2)));
    }
 
    int nref2 = 5;
    int iter = 0;
    bool done = true;
    while(iter<=nref2){
        done = true;
        while(done)
        {
            nocts = pablo1.getNumOctants();
            for (unsigned int i=0; i<nocts; i++){
                array<double,3> center = pablo1.getCenter(i);
                double x = center[0];
                double y = center[1];
                if ((pow((x-xc),2.0)+pow((y-yc),2.0) < pow(radius,2.0)) &&
                        (y>yc) && iter<=nref2 && pablo1.getLevel(i)<=iter+1){
 
                    pablo1.setMarker(i, 1);
                }
            }
            done = pablo1.adapt();
            pablo1.updateConnectivity();
            pablo1.write("pablo00001_iter"+to_string(static_cast<unsigned long long>(iter+nref1+2)));
        }
        iter++;
    }
    pablo1.adaptGlobalRefine();
    pablo1.updateConnectivity();
    pablo1.write("pablo00001_iter"+to_string(static_cast<unsigned long long>(iter+nref1+3)));
}
 
int main(int argc, char *argv[])
{
#if BITPIT_ENABLE_MPI==1
    MPI_Init(&argc,&argv);
#else
    BITPIT_UNUSED(argc);
    BITPIT_UNUSED(argv);
#endif
 
    int nProcs;
    int rank;
#if BITPIT_ENABLE_MPI==1
    MPI_Comm_size(MPI_COMM_WORLD, &nProcs);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
#else
    nProcs = 1;
    rank   = 0;
#endif
 
    // Initialize the logger
    log::manager().initialize(log::MODE_SEPARATE, false, nProcs, rank);
    log::cout() << log::fileVerbosity(log::LEVEL_INFO);
    log::cout() << log::disableConsole();
 
    // Run the example
    try {
        run();
    } catch (const std::exception &exception) {
        log::cout() << exception.what();
        exit(1);
    }
 
#if BITPIT_ENABLE_MPI==1
    MPI_Finalize();
#endif
}