multiblocklattice2d.hh

open lattice boltzmann project www.openlb.org

    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::staticCollide (
        int x0_, int x1_, int y0_, int y1_,
        TensorFieldBase2D<T,2> const& u)
{
    OLB_ASSERT(false, "Method not yet implemented");
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::staticCollide(TensorFieldBase2D<T,2> const& u)
{
    OLB_ASSERT(false, "Method not yet implemented");
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::stream(int x0_, int x1_, int y0_, int y1_)
{
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getNonPeriodicEnvelope(), inters ) ) {
            inters = params.toLocal(inters);
            blockLattices[iBlock] -> stream(inters.x0, inters.x1, inters.y0, inters.y1);
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::stream(bool periodic) {
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    if (periodic) {
        for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
            int iBlock = relevantBlocks[rBlock];
            blockLattices[iBlock] -> stream();
        }
    }
    else {
        for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
            int iBlock = relevantBlocks[rBlock];
            BlockParameters2D const& params = getParameters(iBlock);
            BlockCoordinates2D npEnv = params.toLocal(params.getNonPeriodicEnvelope());
            blockLattices[iBlock] -> stream(npEnv.x0, npEnv.x1, npEnv.y0, npEnv.y1);
            blockLattices[iBlock] -> postProcess();
        }
    }
    postProcessMultiBlock();
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::collideAndStream(int x0_, int x1_, int y0_, int y1_) {
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getNonPeriodicEnvelope(), inters ) ) {
            inters = params.toLocal(inters);
            blockLattices[iBlock] -> stream(inters.x0, inters.x1, inters.y0, inters.y1);
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::collideAndStream(bool periodic) {
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    if (periodic) {
        for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
            int iBlock = relevantBlocks[rBlock];
            blockLattices[iBlock] -> collideAndStream();
        }
    }
    else {
        for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
            int iBlock = relevantBlocks[rBlock];
            BlockParameters2D const& params = getParameters(iBlock);
            BlockCoordinates2D npEnv = params.toLocal(params.getNonPeriodicEnvelope());
            blockLattices[iBlock] -> collideAndStream(npEnv.x0, npEnv.x1, npEnv.y0, npEnv.y1);
            blockLattices[iBlock] -> postProcess();
        }
    }
    postProcessMultiBlock();
}

template<typename T, template<typename U> class Lattice>
T MultiBlockLattice2D<T,Lattice>::computeAverageDensity(int x0_, int x1_, int y0_, int y1_) const {
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    T sumWeights = T(), sumDensities = T();
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getBulk(), inters ) ) {
            inters = params.toLocal(inters);
            T weight = (T) ( (inters.x1-inters.x0+1) * (inters.y1-inters.y0+1) ) /
                       (T) std::numeric_limits<int>::max();
            sumWeights += weight;
            T newDensity = blockLattices[iBlock] -> computeAverageDensity (
                               inters.x0, inters.x1, inters.y0, inters.y1 );
            sumDensities += newDensity * weight;
        }
    }
    if (sumWeights > 1.e-12) {
        sumDensities /= sumWeights;
    }
    multiBlockHandler->reduceAverage(sumDensities, sumWeights);
    return sumDensities;
}

template<typename T, template<typename U> class Lattice>
T MultiBlockLattice2D<T,Lattice>::computeAverageDensity() const {
    return computeAverageDensity(0, getNx()-1, 0, getNy()-1);
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::stripeOffDensityOffset(int x0_, int x1_, int y0_, int y1_, T offset) {
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getEnvelope(), inters ) ) {
            inters = params.toLocal(inters);
            blockLattices[iBlock] -> stripeOffDensityOffset (
                    inters.x0, inters.x1, inters.y0, inters.y1, offset );
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::stripeOffDensityOffset(T offset) {
    stripeOffDensityOffset(0, getNx()-1, 0,getNy()-1, offset);
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::forAll(int x0_, int x1_, int y0_, int y1_,
                                            WriteCellFunctional<T,Lattice> const& application)
{
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getEnvelope(), inters ) ) {
            inters = params.toLocal(inters);
            blockLattices[iBlock] -> forAll (
                    inters.x0, inters.x1, inters.y0, inters.y1, application );
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::forAll(WriteCellFunctional<T,Lattice> const& application)
{
    forAll(0, getNx()-1, 0,getNy()-1, application);
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::addPostProcessor(PostProcessorGenerator2D<T,Lattice> const& ppGen)
{
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters2D const& params = getParameters(iBlock);
        BlockCoordinates2D const& bulk = params.getBulk();
        PostProcessorGenerator2D<T,Lattice> *extractedPpGen = ppGen.clone();
        if (extractedPpGen->extract( bulk.x0, bulk.x1, bulk.y0, bulk.y1 ) ) {
            BlockCoordinates2D const& envelope = params.getEnvelope();
            extractedPpGen->shift(-envelope.x0, -envelope.y0);
            blockLattices[iBlock] -> addPostProcessor(*extractedPpGen);
        }
        delete extractedPpGen;
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::addLatticeCoupling (
                     LatticeCouplingGenerator2D<T,Lattice> const& lcGen,
                     std::vector<SpatiallyExtendedObject2D*> partners )
{
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        std::vector<SpatiallyExtendedObject2D*> extractedPartners(partners.size());
        for (unsigned iP=0; iP<extractedPartners.size(); ++iP) {
            extractedPartners[iP] = partners[iP]->getComponent(iBlock);
        }
        BlockParameters2D const& params = getParameters(iBlock);
        BlockCoordinates2D const& bulk = params.getBulk();
        LatticeCouplingGenerator2D<T,Lattice> *extractedLcGen = lcGen.clone();
        if (extractedLcGen->extract( bulk.x0, bulk.x1, bulk.y0, bulk.y1 ) ) {
            BlockCoordinates2D const& envelope = params.getEnvelope();
            extractedLcGen->shift(-envelope.x0, -envelope.y0);
            blockLattices[iBlock] -> addLatticeCoupling (*extractedLcGen, extractedPartners);
        }
        delete extractedLcGen;
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::resetPostProcessors() {
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        blockLattices[iBlock] -> resetPostProcessors();
    }
}

template<typename T, template<typename U> class Lattice>
LatticeStatistics<T>& MultiBlockLattice2D<T,Lattice>::getStatistics()
{
    return *statistics;
}

template<typename T, template<typename U> class Lattice>
LatticeStatistics<T> const&
        MultiBlockLattice2D<T,Lattice>::getStatistics() const
{
    return *statistics;
}

template<typename T, template<typename U> class Lattice>
DataAnalysisBase2D<T,Lattice> const& MultiBlockLattice2D<T,Lattice>::getDataAnalysis() const {
    dataAnalysis -> reset();
    return *dataAnalysis;
}

template<typename T, template<typename U> class Lattice>
DataSerializer<T> const& MultiBlockLattice2D<T,Lattice>::getSerializer(IndexOrdering::OrderingT ordering) const {
    delete serializer;
    serializer = new MultiSerializer2D<T>(serializerPolicy, ordering);
    return *serializer;
}

template<typename T, template<typename U> class Lattice>
DataUnSerializer<T>& MultiBlockLattice2D<T,Lattice>::getUnSerializer(IndexOrdering::OrderingT ordering) {
    delete unSerializer;
    unSerializer = new MultiUnSerializer2D<T>(unSerializerPolicy, ordering);
    return *unSerializer;
}

template<typename T, template<typename U> class Lattice>
DataSerializer<T> const& MultiBlockLattice2D<T,Lattice>::getSubSerializer (
            int x0_, int x1_, int y0_, int y1_,
            IndexOrdering::OrderingT ordering ) const
{
    delete serializer;
    serializer = new MultiSerializer2D<T> (
            serializerPolicy, x0_, x1_, y0_, y1_, ordering );
    return *serializer;
}

template<typename T, template<typename U> class Lattice>
DataUnSerializer<T>& MultiBlockLattice2D<T,Lattice>::getSubUnSerializer (
            int x0_, int x1_, int y0_, int y1_,
            IndexOrdering::OrderingT ordering )
{
    delete unSerializer;
    unSerializer = new MultiUnSerializer2D<T> (
            unSerializerPolicy, x0_, x1_, y0_, y1_, ordering );
    return *unSerializer;
}



template<typename T, template<typename U> class Lattice>
void MultiBlockLattice2D<T,Lattice>::postProcess(int x0_, int x1_, int y0_, int y1_) {
    BlockCoordinates2D domain(x0_, x1_, y0_, y1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();