libMesh
dtk_evaluator.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 #include "libmesh/libmesh_config.h"
19 
20 #ifdef LIBMESH_TRILINOS_HAVE_DTK
21 
22 #include "libmesh/dtk_evaluator.h"
23 
24 #include "libmesh/dof_map.h"
25 #include "libmesh/fe_interface.h"
26 #include "libmesh/fe_compute_data.h"
27 #include "libmesh/numeric_vector.h"
28 #include "libmesh/equation_systems.h"
29 #include "libmesh/system.h"
30 
31 namespace libMesh
32 {
33 
35  std::string var_name):
36  sys(in_sys),
37  current_local_solution(*sys.current_local_solution),
38  es(in_sys.get_equation_systems()),
39  mesh(sys.get_mesh()),
40  dim(mesh.mesh_dimension()),
41  dof_map(sys.get_dof_map()),
42  var_num(sys.variable_number(var_name)),
43  fe_type(dof_map.variable_type(var_num))
44 {}
45 
47 DTKEvaluator::evaluate(const Teuchos::ArrayRCP<int> & elements,
48  const Teuchos::ArrayRCP<double> & coords)
49 {
50  unsigned int num_values = elements.size();
51 
52  Teuchos::ArrayRCP<Number> values(num_values);
53  DataTransferKit::FieldContainer<Number> evaluations(values, 1);
54 
55  for (unsigned int i=0; i<num_values; i++)
56  {
57  Elem * elem = mesh.elem_ptr(elements[i]);
58 
59  Point p;
60 
61  for (unsigned int j=0; j<dim; j++)
62  p(j) = coords[(j*num_values)+i];
63 
64  const Point mapped_point(FEInterface::inverse_map(dim, dof_map.variable_type(0), elem, p));
65 
66  FEComputeData data (es, mapped_point);
67  FEInterface::compute_data (dim, fe_type, elem, data);
68 
69  std::vector<dof_id_type> dof_indices;
70  dof_map.dof_indices(elem, dof_indices, var_num);
71 
72  Number value = 0;
73 
74  for (std::size_t j=0; j<dof_indices.size(); j++)
75  value += current_local_solution(dof_indices[j]) * data.shape[j];
76 
77  values[i] = value;
78  }
79 
80  return evaluations;
81 }
82 
83 } // namespace libMesh
84 
85 #endif
DTKEvaluator(System &in_sys, std::string var_name)
Definition: dtk_evaluator.C:34
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
unsigned int dim
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
NumericVector< Number > & current_local_solution
Definition: dtk_evaluator.h:73
ImplicitSystem & sys
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
DataTransferKit::FieldContainer< Number > FieldContainerType
Definition: dtk_evaluator.h:64
const FEType & fe_type
Definition: dtk_evaluator.h:79
The libMesh namespace provides an interface to certain functionality in the library.
unsigned int var_num
Definition: dtk_evaluator.h:78
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
Definition: fe_interface.C:827
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:569
std::vector< Number > shape
Storage for the computed shape function values.
virtual FieldContainerType evaluate(const Teuchos::ArrayRCP< int > &elements, const Teuchos::ArrayRCP< double > &coords) libmesh_override
Definition: dtk_evaluator.C:47
static const bool value
Definition: xdr_io.C:108
IterBase * data
Ideally this private member data should have protected access.
virtual const Elem * elem_ptr(const dof_id_type i) const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
sys get_dof_map()
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:1917
EquationSystems & es
Definition: dtk_evaluator.h:74