GRINS-0.8.0
neumann_bc_function_base.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"
31 
32 // libMesh
33 #include "libmesh/libmesh_common.h"
34 #include "libmesh/point.h"
35 #include "libmesh/quadrature.h"
36 #include "libmesh/elem.h"
37 
38 namespace GRINS
39 {
40  template<typename FunctionType, typename FEShape>
42  AssemblyContext& context,
43  libMesh::Real sign,
44  bool is_axisymmetric )
45  {
46  libmesh_assert(!_vars.empty());
47 
48  unsigned int n_vars = this->_vars.size();
49 
50  // The number of local degrees of freedom in each variable.
51  // We're assuming that each var in our set is the same FEType
52  const unsigned int n_var_dofs = context.get_dof_indices(this->_vars[0]).size();
53 #ifndef NDEBUG
54  for( unsigned int v = 1; v < n_vars; v++ )
55  libmesh_assert_equal_to( n_var_dofs, context.get_dof_indices(this->_vars[v]).size() );
56 #endif
57 
58  libMesh::FEGenericBase<FEShape>* side_fe = NULL;
59  context.get_side_fe( this->_vars[0], side_fe );
60 
61  // Element Jacobian * quadrature weight for side integration.
62  const std::vector<libMesh::Real> &JxW_side = side_fe->get_JxW();
63 
64  // The var shape functions at side quadrature points.
65  const std::vector<std::vector<FEShape> >& var_phi_side =
66  side_fe->get_phi();
67 
68  // Physical location of the quadrature points
69  const std::vector<libMesh::Point>& var_qpoint =
70  side_fe->get_xyz();
71 
72  std::vector<libMesh::DenseSubVector<libMesh::Number>*> F_vars(_vars.size(),NULL);
73 
74  for( unsigned int v = 0; v < n_vars; v++ )
75  F_vars[v] = &(context.get_elem_residual(this->_vars[v])); // residual
76 
77  unsigned int n_qpoints = context.get_side_qrule().n_points();
78 
79  FEShape value = 0.0;
80 
81  for (unsigned int qp=0; qp != n_qpoints; qp++)
82  {
83  libMesh::Real jac = JxW_side[qp];
84 
85  if( is_axisymmetric )
86  {
87  const libMesh::Number r = var_qpoint[qp](0);
88  jac *= r;
89  }
90 
91  for( unsigned int v = 0; v < n_vars; v++ )
92  {
93  value = this->eval_func( context, var_qpoint[qp], context.get_time(), this->_vars[v], *_func );
94 
95  for (unsigned int i=0; i != n_var_dofs; i++)
96  (*F_vars[v])(i) += sign*(value*var_phi_side[i][qp])*jac;
97  }
98  }
99 
101  libmesh_error_msg("ERROR: Analytical Jacobians of FEMFunctionBase objects currently not supported!");
102 
103  return compute_jacobian;
104  }
105 
106  // Instantiate
107  template class NeumannBCFunctionBase<libMesh::FunctionBase<libMesh::Number> ,libMesh::Number>;
109  template class NeumannBCFunctionBase<libMesh::FunctionBase<libMesh::Gradient> ,libMesh::Gradient>;
111 
112 } // end namespace GRINS
GRINS namespace.
virtual bool eval_flux(bool compute_jacobian, AssemblyContext &context, libMesh::Real sign, bool is_axisymmetric)

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