levelSetKernel.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/>.
 *
\*---------------------------------------------------------------------------*/
 
# include "levelSetKernel.hpp"
 
# include "bitpit_CG.hpp"
 
namespace bitpit {


LevelSetKernel::LevelSetKernel() {
    m_mesh           = NULL ;
    m_expectedFillIn = LevelSetFillIn::SPARSE ;
 
#if BITPIT_ENABLE_MPI
    m_communicator = MPI_COMM_NULL;
# endif
 
}

LevelSetKernel::LevelSetKernel( VolumeKernel *patch, LevelSetFillIn expectedFillIn ): LevelSetKernel() {
    m_mesh           = patch ;
    m_expectedFillIn = expectedFillIn ;
 
#if BITPIT_ENABLE_MPI
    // Initialize the communicator
    if (m_mesh->isPartitioned()) {
        initializeCommunicator();
    }
#endif
}

LevelSetKernel::~LevelSetKernel(){
    m_mesh = NULL ;
 
# if BITPIT_ENABLE_MPI
    freeCommunicator();
# endif
 
}

VolumeKernel * LevelSetKernel::getMesh() const{
    return m_mesh ;
}

LevelSetFillIn LevelSetKernel::getExpectedFillIn() const
{
    return m_expectedFillIn;
}

double LevelSetKernel::getDistanceTolerance() const
{
    return m_mesh->getTol();
}

bool LevelSetKernel::update( const std::vector<adaption::Info> &adaptionData ) {
 
    // Early return if the kernel doesn't need to be updated
#if BITPIT_ENABLE_MPI
    VolumeKernel *mesh = getMesh();
    bool isMeshPartitioned = mesh->isPartitioned();
#endif
 
    // Check if the levelset needs to be updated
    bool isDirty = false;
    for( const adaption::Info &adaptionInfo : adaptionData){
        if( adaptionInfo.entity != adaption::Entity::ENTITY_CELL){
            continue;
        }
 
        isDirty = true;
        break;
    }
 
#if BITPIT_ENABLE_MPI
    if (isMeshPartitioned) {
        MPI_Allreduce(MPI_IN_PLACE, &isDirty, 1, MPI_C_BOOL, MPI_LOR, getCommunicator());
    }
#endif
 
    // Early return if no changes were detected
    if (!isDirty) {
        return false;
    }
 
#if BITPIT_ENABLE_MPI
    // Set the communicator
    //
    // The mesh may have been partitioned after being set.
    if (m_mesh->isPartitioned() && !isCommunicatorSet()) {
        initializeCommunicator();
    }
#endif
 
    // The kernel has been updated
    return true;
}

bool LevelSetKernel::intersectCellPlane( long id, const std::array<double,3> &root, const std::array<double,3> &normal, double tolerance ) {
 
    std::array<double,3> minPoint;
    std::array<double,3> maxPoint;
    m_mesh->evalCellBoundingBox(id, &minPoint, &maxPoint);
 
    int dim = m_mesh->getDimension();
    return CGElem::intersectPlaneBox( root, normal, minPoint, maxPoint, dim, tolerance);
}

bool LevelSetKernel::isPointInCell(long id, const std::array<double,3> &pointCoords) const {
    return m_mesh->isPointInside(id,pointCoords);
}

double LevelSetKernel::isCellInsideBoundingBox( long id, const std::array<double, 3> &boxMin, const std::array<double, 3> &boxMax ) const {
 
    double distanceTolerance = getDistanceTolerance();
 
    std::array<double,3> cellBoxMin;
    std::array<double,3> cellBoxMax;
    m_mesh->evalCellBoundingBox(id, &cellBoxMin, &cellBoxMax);
 
    return CGElem::intersectBoxBox(cellBoxMin, cellBoxMax, boxMin, boxMax, 3, distanceTolerance);
}
 
# if BITPIT_ENABLE_MPI

void LevelSetKernel::initializeCommunicator()
{
    // The communicator has to be valid
    MPI_Comm communicator = m_mesh->getCommunicator();
    if (communicator == MPI_COMM_NULL) {
        throw std::runtime_error ("Levelset communicator is not valid");
    }
 
    // Communication can be set just once
    if (isCommunicatorSet()) {
        int comparison;
        MPI_Comm_compare(communicator, m_communicator, &comparison);
        if (comparison == MPI_IDENT || comparison == MPI_CONGRUENT) {
            return;
        } else {
            throw std::runtime_error ("Levelset communicator cannot be modified.");
        }
    }
 
    // Create a duplicate of the patch communicator
    MPI_Comm_dup(communicator, &m_communicator);
}

void LevelSetKernel::freeCommunicator() {
 
    if (!isCommunicatorSet()) {
        return;
    }
 
    int finalizedCalled;
    MPI_Finalized(&finalizedCalled);
    if (finalizedCalled) {
        return;
    }
 
    MPI_Comm_free(&m_communicator);
}

MPI_Comm LevelSetKernel::getCommunicator() const {
 
    return m_communicator;
}

bool LevelSetKernel::isCommunicatorSet() const {
 
    return (getCommunicator() != MPI_COMM_NULL);
}

std::unique_ptr<DataCommunicator> LevelSetKernel::createDataCommunicator( ) const {
 
    std::unique_ptr<DataCommunicator> dataCommunicator;
    if (isCommunicatorSet()) {
        dataCommunicator = std::unique_ptr<DataCommunicator>(new DataCommunicator(getCommunicator())) ;
    }
 
    return dataCommunicator;
}
 
#endif


LevelSetCachedKernel::LevelSetCachedKernel( VolumeKernel *patch, LevelSetFillIn fillIn)
    : LevelSetKernel(patch, fillIn),
      m_cellCacheCollection(std::unique_ptr<CellCacheCollection>(new CellCacheCollection(&(patch->getCells()))))
{
}

void LevelSetCachedKernel::clearCaches(bool release)
{
    for (CellCacheCollection::Item &cacheItem : *m_cellCacheCollection) {
        if (!cacheItem.hasCache()) {
            continue;
        }
 
        if (release) {
            cacheItem.destroyCache();
        } else {
            cacheItem.getCache()->clear();
        }
    }
}

LevelSetCachedKernel::CellCacheCollection & LevelSetCachedKernel::getCellCacheCollection()
{
    return *m_cellCacheCollection;
}

const LevelSetCachedKernel::CellCacheCollection & LevelSetCachedKernel::getCellCacheCollection() const
{
    return *m_cellCacheCollection;
}

bool LevelSetCachedKernel::update( const std::vector<adaption::Info> &adaptionData ) {
 
    // Update base class
    bool updated = LevelSetKernel::update( adaptionData );
    if (!updated) {
        return false;
    }
 
    // Update cell cache collection
    for (const adaption::Info &adaptionInfo : adaptionData) {
        if (adaptionInfo.entity != adaption::Entity::ENTITY_CELL) {
            continue;
        }
 
        for (CellCacheCollection::Item &cacheItem : *m_cellCacheCollection) {
            if (!cacheItem.hasCache()) {
                continue;
            }
 
            CellCacheCollection::Cache *cache = cacheItem.getCache();
            cache->erase(adaptionInfo.previous);
        }
    }
 
    for (CellCacheCollection::Item &cacheItem : *m_cellCacheCollection) {
        if (!cacheItem.hasCache()) {
            continue;
        }
 
        CellCacheCollection::Cache *cache = cacheItem.getCache();
        if (m_expectedFillIn == LevelSetFillIn::DENSE) {
            cache->reserve(getMesh()->getCells().size());
        } else {
            cache->shrink_to_fit();
        }
    }
 
    // The kernel has been updated
    return true;
 
}
 
}