BoundaryFromDomainBorder.impl.h 7.27 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108
//======================================================================================================================
//
//  This file is part of waLBerla. waLBerla is free software: you can 
//  redistribute it and/or modify it under the terms of the GNU General Public
//  License as published by the Free Software Foundation, either version 3 of 
//  the License, or (at your option) any later version.
//  
//  waLBerla 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 General Public License 
//  for more details.
//  
//  You should have received a copy of the GNU General Public License along
//  with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
//
//! \file BoundaryFromDomainBorder.impl.h
//! \ingroup geometry
//! \author Martin Bauer <martin.bauer@fau.de>
//
//======================================================================================================================


#include <boost/algorithm/string.hpp>

namespace walberla {
namespace geometry {
namespace initializer {


//*******************************************************************************************************************
/*! Constructor for BoundaryFromDomainBorder
*
* \param blocks                 the structured block storage
* \param handlerBlockDataID     the id of the BoundaryHandler, which type is the template parameter Handling
* \param globalDomain           This parameter can be used if the simulation domain should be smaller than the
*                               number of total cells the StructuredBlockStorage has. Instead of using
*                               blockStorage.getDomainCellBB() here a custom cell interval can be given.
*                               If the empty CellInterval is given (default), the complete domain
*                               is used.
*/
//*******************************************************************************************************************
template<typename Handling>
BoundaryFromDomainBorder<Handling>::BoundaryFromDomainBorder( StructuredBlockStorage & blocks,
                                                              BlockDataID handlerBlockDataID,
                                                              CellInterval globalDomain )
   : blocks_ ( blocks), handlerBlockDataID_( handlerBlockDataID )
{
   if ( globalDomain.empty() )
      globalDomain_ = blocks.getDomainCellBB();
   else
      globalDomain_ = globalDomain;
}


template<typename Handling>
void BoundaryFromDomainBorder<Handling>::init( BlockStorage & /*blockstorage*/, const Config::BlockHandle & blockHandle )
{
   init( blockHandle );
}

template<typename Handling>
void BoundaryFromDomainBorder<Handling>::init( const Config::BlockHandle & blockHandle )
{
   BoundarySetter<Handling> boundarySetter;
   boundarySetter.setConfigBlock( blockHandle );

   std::string directionStr = blockHandle.getParameter<std::string>( "direction" );
   cell_idx_t  wallDistance            = blockHandle.getParameter<cell_idx_t>( "walldistance", 0 );
   cell_idx_t  ghostLayersToInitialize = blockHandle.getParameter<cell_idx_t>( "ghostLayersToInitialize", std::numeric_limits<cell_idx_t>::max() );

   std::set<std::string> directionStrings;
   boost::split( directionStrings, directionStr,boost::is_any_of(","));

   bool atLeastOneBoundarySet = false;
   using stencil::D3Q7;
   for( auto dirIt = D3Q7::beginNoCenter(); dirIt != D3Q7::end(); ++dirIt )
   {
      bool isAll = boost::algorithm::iequals( directionStr, "all" );
      bool isInDirectionStrings = directionStrings.find( stencil::dirToString[*dirIt] ) != directionStrings.end();

      if( isAll || isInDirectionStrings )
      {
         init( *dirIt, wallDistance, ghostLayersToInitialize, boundarySetter );
         atLeastOneBoundarySet = true;
      }
   }

   if ( ! atLeastOneBoundarySet )
      WALBERLA_ABORT( "Invalid Direction " << directionStr << ". Allowed values: all, N,S,W,E,T,B ");
}


template<typename Handling>
void BoundaryFromDomainBorder<Handling>::init( stencil::Direction direction,
                                               cell_idx_t wallDistance,
                                               cell_idx_t ghostLayersToInitialize,
                                               BoundarySetter<Handling> & boundarySetter )
{
   for( auto blockIt = blocks_.begin(); blockIt != blocks_.end(); ++blockIt )
   {
      boundarySetter.configure( *blockIt, handlerBlockDataID_ );

      const cell_idx_t gls = std::min( ghostLayersToInitialize, cell_idx_c( boundarySetter.getFlagField()->nrOfGhostLayers() ) );

      CellInterval ci;

      const cell_idx_t wd = wallDistance;

109
      CellInterval dBB = blocks_.getDomainCellBB( blocks_.getLevel(*blockIt) );
110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172

      for( uint_t dim = 0; dim< 3; ++dim )
         switch ( stencil::c[dim][direction] )
         {
            case -1: ci.min()[dim] =          wd;        ci.max()[dim] =  wd ;                   break;
            case  0: ci.min()[dim] =       - gls;        ci.max()[dim] =  dBB.max()[dim] + gls;  break;
            case  1: ci.min()[dim] = dBB.max()[dim]-wd;  ci.max()[dim] =  dBB.max()[dim] - wd;   break;
         }


      auto blockCellBB = blocks_.getBlockCellBB( *blockIt );
      blockCellBB.expand( cell_idx_c( boundarySetter.getFlagField()->nrOfGhostLayers() )  );
      ci.intersect( blockCellBB );
      blocks_.transformGlobalToBlockLocalCellInterval( ci, *blockIt );
      boundarySetter.set( ci );
   }
}

template<typename Handling>
void BoundaryFromDomainBorder<Handling>::init( FlagUID flagUID, stencil::Direction direction,
                                               cell_idx_t wallDistance, cell_idx_t ghostLayersToInitialize )
{
   BoundarySetter<Handling> boundarySetter;
   boundarySetter.setFlagUID( flagUID );
   init( direction, wallDistance, ghostLayersToInitialize, boundarySetter );
}

template<typename Handling>
void BoundaryFromDomainBorder<Handling>::init( BoundaryUID boundaryUID, stencil::Direction direction,
                                               const shared_ptr<BoundaryConfiguration> & conf,
                                               cell_idx_t wallDistance, cell_idx_t ghostLayersToInitialize )
{
   BoundarySetter<Handling> boundarySetter;
   boundarySetter.setBoundaryConfig( boundaryUID, conf );
   init( direction, wallDistance, ghostLayersToInitialize, boundarySetter );
}


template<typename Handling>
void BoundaryFromDomainBorder<Handling>::initAllBorders( FlagUID flagUID, cell_idx_t wallDistance,
                                                         cell_idx_t ghostLayersToInitialize )
{
   using stencil::D3Q7;
   for( auto dirIt = D3Q7::beginNoCenter(); dirIt != D3Q7::end(); ++dirIt )
      init( flagUID, *dirIt, wallDistance, ghostLayersToInitialize );

}

template<typename Handling>
void BoundaryFromDomainBorder<Handling>::initAllBorders( BoundaryUID boundaryUID, const shared_ptr<BoundaryConfiguration> & conf,
                                                         cell_idx_t wallDistance, cell_idx_t ghostLayersToInitialize )
{
   using stencil::D3Q7;
   for( auto dirIt = D3Q7::beginNoCenter(); dirIt != D3Q7::end(); ++dirIt )
      init( boundaryUID, *dirIt, conf, wallDistance, ghostLayersToInitialize );
}



} // namespace initializer
} // namespace geometry
} // namespace walberla