GRINS-0.8.0
constant_source_term.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 // libMesh
32 #include "libmesh/getpot.h"
33 #include "libmesh/fem_system.h"
34 #include "libmesh/quadrature.h"
35 
36 namespace GRINS
37 {
38  ConstantSourceTerm::ConstantSourceTerm( const std::string& physics_name, const GetPot& input )
39  : SourceTermBase(physics_name,input),
40  _value(0.0)
41  {
42  if( !input.have_variable("Physics/"+physics_name+"/Function/value") )
43  {
44  libMesh::err << "Error: Must specify value for ConstantSourceTerm." << std::endl
45  << " Please specify Physics/"+physics_name+"/Function/value" << std::endl;
46  libmesh_error();
47  }
48 
49  this->set_parameter
50  (_value, input,
51  "Physics/"+physics_name+"/Function/value", _value);
52  }
53 
55  ( bool /*compute_jacobian*/, AssemblyContext & context )
56  {
57  for( std::vector<VariableIndex>::const_iterator v_it = _vars.begin();
58  v_it != _vars.end(); ++v_it )
59  {
60  VariableIndex var = *v_it;
61 
62  // The number of local degrees of freedom in each variable.
63  const unsigned int n_dofs = context.get_dof_indices(var).size();
64 
65  // Element Jacobian * quadrature weights for interior integration.
66  const std::vector<libMesh::Real> &JxW =
67  context.get_element_fe(var)->get_JxW();
68 
69  // The temperature shape functions at interior quadrature points.
70  const std::vector<std::vector<libMesh::Real> >& phi =
71  context.get_element_fe(var)->get_phi();
72 
73  // Get residuals
74  libMesh::DenseSubVector<libMesh::Number> &F_var = context.get_elem_residual(var);
75 
76  // Now we will build the element Jacobian and residual.
77  // Constructing the residual requires the solution and its
78  // gradient from the previous timestep. This must be
79  // calculated at each quadrature point by summing the
80  // solution degree-of-freedom values by the appropriate
81  // weight functions.
82  unsigned int n_qpoints = context.get_element_qrule().n_points();
83 
84  for (unsigned int qp=0; qp != n_qpoints; qp++)
85  {
86  for (unsigned int i=0; i != n_dofs; i++)
87  {
88  F_var(i) += (this->_value)*phi[i][qp]*JxW[qp];
89  }
90  }
91 
92  } // Variable loop
93 
94  return;
95  }
96 
97 } // end namespace GRINS
virtual void set_parameter(libMesh::Number &param_variable, const GetPot &input, const std::string &param_name, libMesh::Number param_default)
Each subclass can simultaneously read a parameter value from.
unsigned int VariableIndex
More descriptive name of the type used for variable indices.
Definition: var_typedefs.h:42
GRINS namespace.
Base class for generic source function term, f(x,y,z,t).
virtual void element_time_derivative(bool compute_jacobian, AssemblyContext &context)
Time dependent part(s) of physics for element interiors.

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