apps/tutorials/lbm/01_BasicLBM.cpp

@@ -28,14 +28,15 @@
 #include "gui/all.h"
 #include "lbm/all.h"
 #include "timeloop/all.h"
-
+#include "GenLM.h"
 
 
 using namespace walberla;
 
-typedef lbm::D2Q9< lbm::collision_model::SRT >  LatticeModel_T;
-typedef LatticeModel_T::Stencil                 Stencil_T;
-typedef LatticeModel_T::CommunicationStencil    CommunicationStencil_T;
+//typedef lbm::D2Q9< lbm::collision_model::SRT >  LatticeModel_T;
+typedef lbm::GenLM                           LatticeModel_T;
+typedef LatticeModel_T::Stencil              Stencil_T;
+typedef LatticeModel_T::CommunicationStencil CommunicationStencil_T;
 
 typedef lbm::PdfField< LatticeModel_T >  PdfField_T;
 
@@ -49,7 +50,6 @@ int main( int argc, char ** argv )
    walberla::Environment walberlaEnv( argc, argv );
 
    auto blocks = blockforest::createUniformBlockGridFromConfig( walberlaEnv.config() );
-
    // read parameters
    auto parameters = walberlaEnv.config()->getOneBlock( "Parameters" );
 
@@ -59,8 +59,11 @@ int main( int argc, char ** argv )
 
    const double remainingTimeLoggerFrequency = parameters.getParameter< double >( "remainingTimeLoggerFrequency", 3.0 ); // in seconds
 
+
    // create fields
-   LatticeModel_T latticeModel = LatticeModel_T( lbm::collision_model::SRT( omega ) );
+   BlockDataID forceFieldId = field::addToStorage<GhostLayerField<real_t, 3> >(blocks, "force", 0.0);
+   LatticeModel_T latticeModel = LatticeModel_T(forceFieldId, omega);
+
    BlockDataID pdfFieldId = lbm::addPdfFieldToStorage( blocks, "pdf field", latticeModel, initialVelocity, real_t(1) );
    BlockDataID flagFieldId = field::addFlagFieldToStorage< FlagField_T >( blocks, "flag field" );
 
@@ -77,6 +80,7 @@ int main( int argc, char ** argv )
                                                                              boundariesConfig.getParameter< real_t > ( "pressure0", real_c( 1.0 ) ),
                                                                              boundariesConfig.getParameter< real_t > ( "pressure1", real_c( 1.0 ) ) );
 
+
    geometry::initBoundaryHandling<BHFactory::BoundaryHandling>( *blocks, boundaryHandlingId, boundariesConfig );
    geometry::setNonBoundaryCellsToDomain<BHFactory::BoundaryHandling> ( *blocks, boundaryHandlingId );
 
@@ -87,20 +91,26 @@ int main( int argc, char ** argv )
    blockforest::communication::UniformBufferedScheme< CommunicationStencil_T > communication( blocks );
    communication.addPackInfo( make_shared< lbm::PdfFieldPackInfo< LatticeModel_T > >( pdfFieldId ) );
 
+
+   LatticeModel_T::Sweep streamCollide(pdfFieldId);
+
    // add LBM sweep and communication to time loop
    timeloop.add() << BeforeFunction( communication, "communication" )
                   << Sweep( BHFactory::BoundaryHandling::getBlockSweep( boundaryHandlingId ), "boundary handling" );
-   timeloop.add() << Sweep( makeSharedSweep( lbm::makeCellwiseSweep< LatticeModel_T, FlagField_T >( pdfFieldId, flagFieldId, fluidFlagUID ) ), "LB stream & collide" );
+   timeloop.add() << Sweep( streamCollide );
 
    // LBM stability check
+   /*
    timeloop.addFuncAfterTimeStep( makeSharedFunctor( field::makeStabilityChecker< PdfField_T, FlagField_T >( walberlaEnv.config(), blocks, pdfFieldId,
                                                                                                              flagFieldId, fluidFlagUID ) ),
                                   "LBM stability check" );
+   */
 
    // log remaining time
    timeloop.addFuncAfterTimeStep( timing::RemainingTimeLogger( timeloop.getNrOfTimeSteps(), remainingTimeLoggerFrequency ), "remaining time logger" );
 
    // add VTK output to time loop
+   //TODO VTK leads to crash?
    lbm::VTKOutput< LatticeModel_T, FlagField_T >::addToTimeloop( timeloop, blocks, walberlaEnv.config(), pdfFieldId, flagFieldId, fluidFlagUID );
 
    // create adaptors, so that the GUI also displays density and velocity