GRINS-0.8.0
postprocessed_quantities.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // GRINS - General Reacting Incompressible Navier-Stokes
5 //
6 // Copyright (C) 2014-2017 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 // This class
27 
28 // GRINS
29 #include "grins/assembly_context.h"
30 
31 namespace GRINS
32 {
33  template<class NumericType>
35  : libMesh::FEMFunctionBase<NumericType>()
36  {
37  if( input.have_variable("vis-options/output_vars") )
38  {
39 
40  std::cerr << "================================================================================" << std::endl
41  << "WARNING: Detected input variable 'vis-options/output_vars." << std::endl
42  << "WARNING: THIS IS OUTDATED. output_vars is now input within" << std::endl
43  << "WARNING: each Physics section." << std::endl
44  << " Erroring out to make you aware." << std::endl
45  << "================================================================================" << std::endl;
46  libmesh_error();
47  }
48  return;
49  }
50 
51  template<class NumericType>
53  {
54  return;
55  }
56 
57  template<class NumericType>
59  {
60  // Check if this quantity has already been registered
61  if( _quantity_name_index_map.find(name) != _quantity_name_index_map.end() )
62  {
63  std::cerr << "Error: trying to add existing quantity: " << name << std::endl;
64  libmesh_error();
65  }
66 
67  /* We add 1 so that the locally cached indices in each of the Physics
68  can safely initialize the index to zero. In particular this ensures
69  we count starting from 1. */
70  unsigned int new_index = _quantity_name_index_map.size()+1;
71 
72  _quantity_name_index_map.insert( std::make_pair( name, new_index ) );
73 
74  return new_index;
75  }
76 
77  template<class NumericType>
79  libMesh::EquationSystems& equation_systems )
80  {
81  // Only need to initialize if the user requested any output quantities.
82  if( !_quantity_name_index_map.empty() )
83  {
84  // Need to cache the MultiphysicsSystem
85  _multiphysics_sys = &system;
86 
87  libMesh::System& output_system = equation_systems.add_system<libMesh::System>("interior_output");
88 
89  for( std::map<std::string, unsigned int>::const_iterator it = _quantity_name_index_map.begin();
90  it != _quantity_name_index_map.end();
91  ++it )
92  {
93  unsigned int var = output_system.add_variable( it->first, libMesh::FIRST );
94  _quantity_index_var_map.insert( std::make_pair(var,it->second) );
95  }
96  }
97 
98  return;
99  }
100 
101  template<class NumericType>
102  void PostProcessedQuantities<NumericType>::update_quantities( libMesh::EquationSystems& equation_systems )
103  {
104  // Only do the projection if the user actually added any quantities to compute.
105  if( !_quantity_name_index_map.empty() )
106  {
107  libMesh::System& output_system = equation_systems.get_system<libMesh::System>("interior_output");
108  output_system.project_solution(this);
109  }
110 
111  return;
112  }
113 
114 
115  template<class NumericType>
116  NumericType PostProcessedQuantities<NumericType>::component( const libMesh::FEMContext& context,
117  unsigned int component,
118  const libMesh::Point& p,
119  libMesh::Real /*time*/ )
120  {
121  // Check if the Elem is the same between the incoming context and the cached one.
122  // If not, reinit the cached MultiphysicsSystem context
123  if(context.has_elem() && _multiphysics_context->has_elem())
124  {
125  if( &(context.get_elem()) != &(_multiphysics_context->get_elem()) )
126  {
127  _multiphysics_context->pre_fe_reinit(*_multiphysics_sys,&context.get_elem());
128  _multiphysics_context->elem_fe_reinit();
129  }
130  }
131  else if( !context.has_elem() && _multiphysics_context->has_elem() )
132  {
133  // Incoming context has NULL elem ==> SCALAR variables
134  _multiphysics_context->pre_fe_reinit(*_multiphysics_sys,NULL);
135  _multiphysics_context->elem_fe_reinit();
136  }
137  else if( context.has_elem() && !_multiphysics_context->has_elem() )
138  {
139  _multiphysics_context->pre_fe_reinit(*_multiphysics_sys,&context.get_elem());
140  _multiphysics_context->elem_fe_reinit();
141  }
142  //else
143  /* If has_elem() is false for both contexts, we're still dealing with SCALAR variables
144  and therefore don't need to reinit. */
145 
146  libMesh::Real value = 0.0;
147 
148  // Quantity we want had better be there.
149  libmesh_assert(_quantity_index_var_map.find(component) != _quantity_index_var_map.end());
150  unsigned int quantity_index = _quantity_index_var_map.find(component)->second;
151 
152  _multiphysics_sys->compute_postprocessed_quantity( quantity_index,
153  *(this->_multiphysics_context),
154  p, value );
155 
156  return value;
157  }
158 
159  template<class NumericType>
160  void PostProcessedQuantities<NumericType>::init_context( const libMesh::FEMContext& context )
161  {
162  // Make sure we prepare shape functions for our output variables.
164  for( typename std::map<VariableIndex,unsigned int>::const_iterator it = _quantity_index_var_map.begin();
165  it != _quantity_index_var_map.end(); it++ )
166  {
167  libMesh::FEBase* elem_fe = NULL;
168  context.get_element_fe( it->first, elem_fe );
169  elem_fe->get_phi();
170  elem_fe->get_dphi();
171  elem_fe->get_JxW();
172  elem_fe->get_xyz();
173 
174  libMesh::FEBase* side_fe = NULL;
175  context.get_side_fe( it->first, side_fe );
176  side_fe->get_phi();
177  side_fe->get_dphi();
178  side_fe->get_JxW();
179  side_fe->get_xyz();
180  }
181 
182  // Create the context we'll be using to compute MultiphysicsSystem quantities
183  _multiphysics_context.reset( new AssemblyContext( *_multiphysics_sys ) );
184  _multiphysics_sys->init_context(*_multiphysics_context);
185  return;
186  }
187 
188  template<class NumericType>
189  void PostProcessedQuantities<NumericType>::operator()( const libMesh::FEMContext& context, const libMesh::Point& p,
190  const libMesh::Real time,
191  libMesh::DenseVector<NumericType>& output )
192  {
193  for( unsigned int i = 0; i != output.size(); i++ )
194  {
195  output(i) = this->component(context,i,p,time);
196  }
197  return;
198  }
199 
200  template<class NumericType>
201  NumericType PostProcessedQuantities<NumericType>::operator()( const libMesh::FEMContext&,
202  const libMesh::Point&,
203  const libMesh::Real )
204  {
205  libmesh_error();
206  return 0.0; //dummy
207  }
208 
209  // Instantiate
211 
212 } // namespace GRINS
unsigned int register_quantity(std::string name)
Register quantity to be postprocessed.
virtual void initialize(MultiphysicsSystem &system, libMesh::EquationSystems &equation_systems)
virtual NumericType operator()(const libMesh::FEMContext &context, const libMesh::Point &p, const libMesh::Real time=0.)
GRINS namespace.
Interface with libMesh for solving Multiphysics problems.
virtual NumericType component(const libMesh::FEMContext &context, unsigned int i, const libMesh::Point &p, libMesh::Real time=0.)
virtual void init_context(const libMesh::FEMContext &context)
virtual void update_quantities(libMesh::EquationSystems &equation_systems)

Generated on Tue Dec 19 2017 12:47:28 for GRINS-0.8.0 by  doxygen 1.8.9.1