26 #include <lifev/core/LifeV.hpp> 27 #include <lifev/core/util/LifeChrono.hpp> 28 #include <lifev/core/array/MapEpetra.hpp> 29 #include <lifev/core/mesh/MeshData.hpp> 30 #include <lifev/core/mesh/MeshPartitioner.hpp> 31 #include <lifev/core/fem/FESpace.hpp> 32 #include <lifev/navier_stokes/solver/OseenData.hpp> 33 #include <lifev/navier_stokes/solver/OseenSolver.hpp> 34 #include <lifev/core/filter/ExporterEnsight.hpp> 35 #include <lifev/core/filter/ExporterHDF5.hpp> 36 #include <lifev/core/filter/ExporterEmpty.hpp> 43 #include <Epetra_ConfigDefs.h> 45 #include "Epetra_MpiComm.h" 47 #include "Epetra_SerialComm.h" 51 using namespace LifeV;
96 int main (
int argc,
char** argv)
103 MPI_Init (&argc, &argv);
106 std::shared_ptr<Epetra_Comm> comm;
108 comm.reset (
new Epetra_MpiComm ( MPI_COMM_WORLD ) );
110 MPI_Comm_size (MPI_COMM_WORLD, &nproc);
112 comm.reset (
new Epetra_SerialComm() );
115 bool verbose (
false);
116 if (comm->MyPID() == 0)
120 <<
" +-----------------------------------------------+" << std::endl
121 <<
" | Cavity example for LifeV |" << std::endl
122 <<
" +-----------------------------------------------+" << std::endl
124 <<
" +-----------------------------------------------+" << std::endl
125 <<
" | Author: Gwenol Grandperrin |" << std::endl
126 <<
" | Date: October 25, 2010 |" << std::endl
127 <<
" +-----------------------------------------------+" << std::endl
130 std::cout <<
"[Initilization of MPI]" << std::endl;
132 std::cout <<
"Using MPI (" << nproc <<
" proc.)" << std::endl;
134 std::cout <<
"Using serial version" << std::endl;
147 std::cout << std::endl <<
"[Loading the data]" << std::endl;
150 const std::string dataFileName = command_line.follow (
"data", 2,
"-f",
"--file");
151 GetPot dataFile (dataFileName);
158 std::cout << std::endl <<
"[Loading the mesh]" << std::endl;
161 meshData.setup (dataFile,
"fluid/space_discretization");
164 std::cout <<
"Mesh file: " << meshData.meshDir() << meshData.meshFile() << std::endl;
166 std::shared_ptr< mesh_Type > fullMeshPtr (
new mesh_Type);
167 readMesh (*fullMeshPtr, meshData);
176 std::cout << std::endl <<
"[Creating the FE spaces]" << std::endl;
178 std::string uOrder = dataFile (
"fluid/space_discretization/vel_order",
"P2");
179 std::string pOrder = dataFile (
"fluid/space_discretization/press_order",
"P1");
180 if (verbose) std::cout <<
"FE for the velocity: " << uOrder << std::endl
181 <<
"FE for the pressure: " << pOrder << std::endl;
185 std::cout <<
"Building the velocity FE space... " << std::flush;
187 feSpacePtr_Type uFESpacePtr (
new feSpace_Type (meshPart.meshPartition(), uOrder, 3, comm) );
190 std::cout <<
"ok." << std::endl;
195 std::cout <<
"Building the pressure FE space... " << std::flush;
197 feSpacePtr_Type pFESpacePtr (
new feSpace_Type (meshPart.meshPartition(), pOrder, 1, comm) );
200 std::cout <<
"ok." << std::endl;
204 UInt totalVelDof = uFESpacePtr->map().map (Unique)->NumGlobalElements();
205 UInt totalPressDof = pFESpacePtr->map().map (Unique)->NumGlobalElements();
209 std::cout <<
"Total Velocity Dof: " << totalVelDof << std::endl;
213 std::cout <<
"Total Pressure Dof: " << totalPressDof << std::endl;
221 std::cout << std::endl <<
"[Boundary conditions]" << std::endl;
224 BCFunctionBase uZero (fZero);
225 BCFunctionBase uLid (lidBC);
227 std::vector<ID> xComp (1);
232 bcH.addBC (
"Top" , TOP , Essential, Full , uLid , 3 );
233 bcH.addBC (
"Left" , LEFT , Essential, Full , uZero, 3 );
234 bcH.addBC (
"Front" , FRONT , Essential, Component, uZero, xComp );
235 bcH.addBC (
"Right" , RIGHT , Essential, Full , uZero, 3 );
236 bcH.addBC (
"Back" , BACK , Essential, Component, uZero, xComp );
237 bcH.addBC (
"Bottom", BOTTOM, Essential, Full , uZero, 3 );
240 std::vector<bcName_Type> fluxVector = bcH.findAllBCWithType ( Flux );
241 UInt numLM =
static_cast<UInt> ( fluxVector.size() );
243 UInt offset = uFESpacePtr->map().map (Unique)->NumGlobalElements()
244 + pFESpacePtr->map().map (Unique)->NumGlobalElements();
246 for (
UInt i = 0; i < numLM; ++i )
248 bcH.setOffset ( fluxVector[i], offset + i );
256 std::cout << std::endl <<
"[Creating the problem]" << std::endl;
258 std::shared_ptr< OseenData > oseenData (
new OseenData);
259 oseenData->setup ( dataFile );
262 OseenSolver< mesh_Type > fluid (oseenData,
269 fluid.setUp (dataFile);
282 fluid.updateSystem (alpha, beta, rhs);
285 std::shared_ptr< ExporterHDF5<mesh_Type> > exporter;
287 std::string
const exporterType = dataFile (
"exporter/type",
"ensight");
289 exporter.reset (
new ExporterHDF5<mesh_Type> ( dataFile,
"cavity_stokes_example" ) );
290 exporter->setPostDir (
"./" );
291 exporter->setMeshProcId ( meshPart.meshPartition(), comm->MyPID() );
294 velAndPressure.reset (
new vector_Type (*fluid.solution(), exporter->mapType() ) );
296 exporter->addVariable ( ExporterData<mesh_Type>::VectorField,
"velocity", uFESpacePtr,
297 velAndPressure, UInt (0) );
299 exporter->addVariable ( ExporterData<mesh_Type>::ScalarField,
"pressure", pFESpacePtr,
300 velAndPressure, UInt (3 * uFESpacePtr->dof().numTotalDof() ) );
301 exporter->postProcess ( 0 );
306 std::cout <<
"Total simulation time: " << globalChrono.diff() <<
" s." << std::endl;
309 exporter->closeFile();
Real fZero(const Real &, const Real &, const Real &, const Real &, const ID &)
std::shared_ptr< vector_Type > vectorPtr_Type
void start()
Start the timer.
FESpace - Short description here please!
RegionMesh< LinearTetra > mesh_Type
GetPot(const int argc_, char **argv_, const char *FieldSeparator=0x0)
fluid_Type::vector_Type vector_Type
Real lidBC(const Real &, const Real &, const Real &, const Real &, const ID &i)
OseenSolver< mesh_Type > fluid_Type
Epetra_Import const & importer()
Getter for the Epetra_Import.
FESpace< mesh_Type, MapEpetra > feSpace_Type
std::shared_ptr< std::vector< Int > > M_isOnProc
MeshData - class for handling spatial discretization.
double Real
Generic real data.
void stop()
Stop the timer.
std::shared_ptr< feSpace_Type > feSpacePtr_Type
int main(int argc, char **argv)
uint32_type UInt
generic unsigned integer (used mainly for addressing)