26 #include "grins_config.h"
38 #include "libmesh/exact_solution.h"
39 #include "libmesh/zero_function.h"
42 #ifdef GRINS_HAVE_GRVY
48 const libMesh::Parameters&,
54 const libMesh::Parameters&,
68 std::multimap< GRINS::PhysicsName, GRINS::DBCContainer >
build_dirichlet( );
71 int main(
int argc,
char* argv[])
73 #ifdef GRINS_USE_GRVY_TIMERS
74 GRVY::GRVY_Timer_Class grvy_timer;
75 grvy_timer.Init(
"GRINS Timer");
82 std::cerr <<
"Error: Must specify libMesh input file." << std::endl;
87 std::string libMesh_input_filename = argv[1];
90 GetPot libMesh_inputfile( libMesh_input_filename );
92 #ifdef GRINS_USE_GRVY_TIMERS
93 grvy_timer.BeginTimer(
"Initialize Solver");
97 libMesh::LibMeshInit libmesh_init(argc, argv);
100 libmesh_example_requires
101 (libMesh::default_solver_package() != libMesh::LASPACK_SOLVERS,
112 libmesh_init.comm() );
114 #ifdef GRINS_USE_GRVY_TIMERS
115 grvy_timer.EndTimer(
"Initialize Solver");
118 grins.attach_grvy_timer( &grvy_timer );
125 std::tr1::shared_ptr<libMesh::EquationSystems> es = grins.get_equation_system();
128 libMesh::ExactSolution exact_sol(*es);
134 exact_sol.compute_error(
"GRINS",
"z_vel");
136 double l2error = exact_sol.l2_error(
"GRINS",
"z_vel");
137 double h1error = exact_sol.h1_error(
"GRINS",
"z_vel");
141 if( l2error > 2.0e-9 || h1error > 4.0e-7 )
145 std::cout <<
"Tolerance exceeded for velocity in Poiseuille test." << std::endl
146 <<
"l2 error = " << l2error << std::endl
147 <<
"h1 error = " << h1error << std::endl;
187 std::tr1::shared_ptr<libMesh::FunctionBase<libMesh::Number> >
188 vel_func2(
new libMesh::ZeroFunction<libMesh::Number> );
192 std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > mymap;
203 const libMesh::Parameters& ,
205 const std::string& var )
207 const double r = p(0);
209 const double r0 = 1.0;
210 const double r1 = 2.0;
211 const double u0 = 2.0;
213 libMesh::Number f = 0;
215 if( var ==
"z_vel" ) f = u0*std::log( r1/r )/std::log( r1/r0 );
216 else libmesh_assert(
false);
223 const libMesh::Parameters& ,
225 const std::string& var)
227 const double r = p(0);
229 const double r0 = 1.0;
230 const double r1 = 2.0;
231 const double u0 = 2.0;
238 g(0) = -u0/std::log(r1/r0)*r/r1*r1/(r*r);
const PhysicsName incompressible_navier_stokes
Object for constructing boundary condition function objects.
libMesh::Number exact_solution(const libMesh::Point &p, const libMesh::Parameters &, const std::string &, const std::string &)
Simple helper class to setup general Dirichlet boundary conditions.
libMesh::Gradient exact_derivative(const libMesh::Point &p, const libMesh::Parameters &, const std::string &, const std::string &)
std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > build_dirichlet()
Builds all required libMesh::DirichletBoundary objects and adds them to DofMap.
BoundaryConditionsFactory()
void add_var_name(const GRINS::VariableName &var)
Add variables that are constrained by the Dirichlet bc.
Profile for flow between axially moving concentric cylinders.
int main(int argc, char *argv[])
void add_bc_id(const GRINS::BoundaryID bc_id)
Add boundary id's for which this Dirichlet bc is to be applied.
void set_func(std::tr1::shared_ptr< libMesh::FunctionBase< libMesh::Number > > func)
Add the Dirichlet bc functor.
void attach_bc_factory(std::tr1::shared_ptr< BoundaryConditionsFactory > bc_factory)