GRINS-0.6.0
test_axi_ns_poiseuille_flow.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // GRINS - General Reacting Incompressible Navier-Stokes
5 //
6 // Copyright (C) 2014-2015 Paul T. Bauman, Roy H. Stogner
7 // Copyright (C) 2010-2013 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 
26 #include "grins_config.h"
27 
28 #include <iostream>
29 
30 // GRINS
31 #include "grins/simulation.h"
33 #include "grins/multiphysics_sys.h"
35 
36 //libMesh
37 #include "libmesh/exact_solution.h"
38 
39 // GRVY
40 #ifdef GRINS_HAVE_GRVY
41 #include "grvy.h"
42 #endif
43 
44 libMesh::Number
45 exact_solution( const libMesh::Point& p,
46  const libMesh::Parameters&, // parameters, not needed
47  const std::string&, // sys_name, not needed
48  const std::string&); // unk_name, not needed);
49 
50 libMesh::Gradient
51 exact_derivative( const libMesh::Point& p,
52  const libMesh::Parameters&, // parameters, not needed
53  const std::string&, // sys_name, not needed
54  const std::string&); // unk_name, not needed);
55 
57 {
58 public:
59 
62  { return; };
63 
65 
66  std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > build_dirichlet( );
67 };
68 
69 int main(int argc, char* argv[])
70 {
71 #ifdef GRINS_USE_GRVY_TIMERS
72  GRVY::GRVY_Timer_Class grvy_timer;
73  grvy_timer.Init("GRINS Timer");
74 #endif
75 
76  // Check command line count.
77  if( argc < 2 )
78  {
79  // TODO: Need more consistent error handling.
80  std::cerr << "Error: Must specify libMesh input file." << std::endl;
81  exit(1); // TODO: something more sophisticated for parallel runs?
82  }
83 
84  // libMesh input file should be first argument
85  std::string libMesh_input_filename = argv[1];
86 
87  // Create our GetPot object.
88  GetPot libMesh_inputfile( libMesh_input_filename );
89 
90 #ifdef GRINS_USE_GRVY_TIMERS
91  grvy_timer.BeginTimer("Initialize Solver");
92 #endif
93 
94  // Initialize libMesh library.
95  libMesh::LibMeshInit libmesh_init(argc, argv);
96 
97  GRINS::SimulationBuilder sim_builder;
98 
99  std::tr1::shared_ptr<AxiParabolicBCFactory> bc_factory( new AxiParabolicBCFactory );
100 
101  sim_builder.attach_bc_factory(bc_factory);
102 
103  GRINS::Simulation grins( libMesh_inputfile,
104  sim_builder,
105  libmesh_init.comm() );
106 
107 #ifdef GRINS_USE_GRVY_TIMERS
108  grvy_timer.EndTimer("Initialize Solver");
109 
110  // Attach GRVY timer to solver
111  grins.attach_grvy_timer( &grvy_timer );
112 #endif
113 
114  // Get equation systems to create ExactSolution object
115  std::tr1::shared_ptr<libMesh::EquationSystems> es = grins.get_equation_system();
116 
117  // Do solve here
118  grins.run();
119 
120 
121  // Create Exact solution object and attach exact solution quantities
122  libMesh::ExactSolution exact_sol(*es);
123 
124  exact_sol.attach_exact_value(&exact_solution);
125  exact_sol.attach_exact_deriv(&exact_derivative);
126 
127  // Compute error and get it in various norms
128  exact_sol.compute_error("GRINS", "z_vel");
129 
130  double l2error = exact_sol.l2_error("GRINS", "z_vel");
131  double h1error = exact_sol.h1_error("GRINS", "z_vel");
132 
133  int return_flag = 0;
134 
135  double tol = 1.0e-9;
136 
137  if( l2error > tol || h1error > tol )
138  {
139  return_flag = 1;
140 
141  std::cout << "Tolerance exceeded for velocity in Poiseuille test." << std::endl
142  << "l2 error = " << l2error << std::endl
143  << "h1 error = " << h1error << std::endl;
144  }
145 
146  /*
147  // Compute error and get it in various norms
148  exact_sol.compute_error("GRINS", "p");
149 
150  l2error = exact_sol.l2_error("GRINS", "p");
151  h1error = exact_sol.h1_error("GRINS", "p");
152 
153  if( l2error > 1.0e-11 || h1error > 1.0e-11 )
154  {
155  return_flag = 1;
156 
157  std::cout << "Tolerance exceeded for pressure in Poiseuille test." << std::endl
158  << "l2 error = " << l2error << std::endl
159  << "h1 error = " << h1error << std::endl;
160  }
161  */
162 
163  return return_flag;
164 }
165 
166 std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > AxiParabolicBCFactory::build_dirichlet( )
167 {
168  GRINS::DBCContainer cont;
169  cont.add_var_name( "z_vel" );
170  cont.add_bc_id( 0 );
171  cont.add_bc_id( 2 );
172 
173  std::tr1::shared_ptr<libMesh::FunctionBase<libMesh::Number> >
174  vel_func( new GRINS::ParabolicProfile( -100.0/4.0, 0.0, 0.0, 0.0, 0.0, 100.0/4.0 ) );
175 
176  cont.set_func( vel_func );
177 
178  std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > mymap;
179 
180  mymap.insert( std::pair<GRINS::PhysicsName, GRINS::DBCContainer >(GRINS::incompressible_navier_stokes, cont) );
181 
182  return mymap;
183 }
184 
185 libMesh::Number
186 exact_solution( const libMesh::Point& p,
187  const libMesh::Parameters&, // parameters, not needed
188  const std::string&, // sys_name, not needed
189  const std::string& var ) // unk_name, not needed);
190 {
191  const double r = p(0);
192 
193  const double r0 = 1.0;
194  const double mu = 1.0;
195  const double dpdx = -1.0;
196 
197  libMesh::Number f = 0;
198  // Hardcoded to velocity in input file.
199  if( var == "z_vel" ) f = -dpdx*100.0/(4.0*mu)*(r0*r0 - r*r);
200  else libmesh_assert(false);
201 
202  return f;
203 }
204 
205 libMesh::Gradient
206 exact_derivative( const libMesh::Point& p,
207  const libMesh::Parameters&, // parameters, not needed
208  const std::string&, // sys_name, not needed
209  const std::string& var)
210 {
211  const double r = p(0);
212 
213  libMesh::Gradient g;
214 
215  // Hardcoded to velocity in input file.
216  if( var == "z_vel" )
217  {
218  g(0) = -100.0/2.0*r;
219  g(1) = 0.0;
220  }
221 
222  return g;
223 }
const PhysicsName incompressible_navier_stokes
Object for constructing boundary condition function objects.
Definition: bc_factory.h:50
std::multimap< GRINS::PhysicsName, GRINS::DBCContainer > build_dirichlet()
Builds all required libMesh::DirichletBoundary objects and adds them to DofMap.
Simple helper class to setup general Dirichlet boundary conditions.
Definition: dbc_container.h:49
Parabolic profile.
void add_var_name(const GRINS::VariableName &var)
Add variables that are constrained by the Dirichlet bc.
Definition: dbc_container.C:46
GRINS namespace.
libMesh::Gradient exact_derivative(const libMesh::Point &p, const libMesh::Parameters &, const std::string &, const std::string &)
int main(int argc, char *argv[])
std::tr1::shared_ptr< libMesh::EquationSystems > get_equation_system()
Definition: simulation.C:422
libMesh::Number exact_solution(const libMesh::Point &p, const libMesh::Parameters &, const std::string &, const std::string &)
void add_bc_id(const GRINS::BoundaryID bc_id)
Add boundary id's for which this Dirichlet bc is to be applied.
Definition: dbc_container.C:52
void set_func(std::tr1::shared_ptr< libMesh::FunctionBase< libMesh::Number > > func)
Add the Dirichlet bc functor.
Definition: dbc_container.C:58
void attach_bc_factory(std::tr1::shared_ptr< BoundaryConditionsFactory > bc_factory)

Generated on Mon Jun 22 2015 21:32:20 for GRINS-0.6.0 by  doxygen 1.8.9.1