multiblocklattice3d.hh

来自「open lattice boltzmann project www.open」· HH 代码 · 共 839 行 · 第 1/3 页

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

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

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

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<T,Lattice>::stream(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_)
{
    BlockCoordinates3D domain(x0_, x1_, y0_, y1_, z0_, z1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters3D 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, inters.z0, inters.z1);
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<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];
            BlockParameters3D const& params = getParameters(iBlock);
            BlockCoordinates3D npEnv = params.toLocal(params.getNonPeriodicEnvelope());
            blockLattices[iBlock] -> stream(npEnv.x0, npEnv.x1, npEnv.y0, npEnv.y1, npEnv.z0, npEnv.z1);
            blockLattices[iBlock] -> postProcess();
        }
    }
    postProcessMultiBlock();
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<T,Lattice>::collideAndStream(int x0_, int x1_, int y0_, int y1_, int z0_, int z1_) {
    BlockCoordinates3D domain(x0_, x1_, y0_, y1_, z0_, z1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters3D 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, inters.z0, inters.z1);
        }
    }
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<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];
            BlockParameters3D const& params = getParameters(iBlock);
            BlockCoordinates3D npEnv = params.toLocal(params.getNonPeriodicEnvelope());
            blockLattices[iBlock] ->
                collideAndStream(npEnv.x0, npEnv.x1, npEnv.y0, npEnv.y1, npEnv.z0, npEnv.z1);
            blockLattices[iBlock] -> postProcess();
        }
    }
    postProcessMultiBlock();
}

template<typename T, template<typename U> class Lattice>
T MultiBlockLattice3D<T,Lattice>::computeAverageDensity(int x0_, int x1_, int y0_, int y1_,
                                                        int z0_, int z1_) const
{
    BlockCoordinates3D domain(x0_, x1_, y0_, y1_, z0_, z1_), inters;
    T sumWeights = T(), sumDensities = T();
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters3D const& params = getParameters(iBlock);
        if (util::intersect(domain, params.getEnvelope(), inters ) ) {
            inters = params.toLocal(inters);
            T weight = (T) ( (inters.x1-inters.x0+1) * (inters.y1-inters.y0+1) * (inters.z1-inters.z0+1) ) /
                       (T) std::numeric_limits<int>::max();
            sumWeights += weight;
            T newDensity = blockLattices[iBlock] -> computeAverageDensity (
                               inters.x0, inters.x1, inters.y0, inters.y1, inters.z0, inters.z1 );
            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 MultiBlockLattice3D<T,Lattice>::computeAverageDensity() const {
    return computeAverageDensity(0, getNx()-1, 0, getNy()-1, 0, getNz()-1);
}

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<T,Lattice>::stripeOffDensityOffset (
        int x0_, int x1_, int y0_, int y1_, int z0_, int z1_, T offset )
{
    BlockCoordinates3D domain(x0_, x1_, y0_, y1_, z0_, z1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters3D 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, inters.z0, inters.z1, offset );
        }
    }
}

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

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<T,Lattice>::forAll (
        int x0_, int x1_, int y0_, int y1_, int z0_, int z1_,
        WriteCellFunctional<T,Lattice> const& application )
{
    BlockCoordinates3D domain(x0_, x1_, y0_, y1_, z0_, z1_), inters;
    std::vector<int> const& relevantBlocks = getRelevantBlocks();
    for (int rBlock=0; rBlock < getNumRelevantBlocks(); ++rBlock) {
        int iBlock = relevantBlocks[rBlock];
        BlockParameters3D 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, inters.z0, inters.z1, application );
        }
    }
}

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

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

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

template<typename T, template<typename U> class Lattice>
void MultiBlockLattice3D<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>& MultiBlockLattice3D<T,Lattice>::getStatistics()
{
    return *statistics;
}

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

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

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

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

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

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