libMesh
Public Types | Public Member Functions | Protected Attributes | List of all members
libMesh::DTKEvaluator Class Reference

Implements the evaluate() function to compute FE solution values at points requested by DTK. More...

#include <dtk_evaluator.h>

Inheritance diagram for libMesh::DTKEvaluator:
[legend]

Public Types

typedef DataTransferKit::MeshContainer< intMeshContainerType
 
typedef DataTransferKit::FieldContainer< NumberFieldContainerType
 

Public Member Functions

 DTKEvaluator (System &in_sys, std::string var_name)
 
virtual FieldContainerType evaluate (const Teuchos::ArrayRCP< int > &elements, const Teuchos::ArrayRCP< double > &coords) override
 

Protected Attributes

Systemsys
 
NumericVector< Number > & current_local_solution
 
EquationSystemses
 
MeshBasemesh
 
unsigned int dim
 
DofMapdof_map
 
unsigned int var_num
 
const FETypefe_type
 

Detailed Description

Implements the evaluate() function to compute FE solution values at points requested by DTK.

Author
Derek Gaston
Date
2013

Definition at line 60 of file dtk_evaluator.h.

Member Typedef Documentation

◆ FieldContainerType

typedef DataTransferKit::FieldContainer<Number> libMesh::DTKEvaluator::FieldContainerType

Definition at line 64 of file dtk_evaluator.h.

◆ MeshContainerType

typedef DataTransferKit::MeshContainer<int> libMesh::DTKEvaluator::MeshContainerType

Definition at line 63 of file dtk_evaluator.h.

Constructor & Destructor Documentation

◆ DTKEvaluator()

libMesh::DTKEvaluator::DTKEvaluator ( System in_sys,
std::string  var_name 
)

Definition at line 34 of file dtk_evaluator.C.

35  :
36  sys(in_sys),
38  es(in_sys.get_equation_systems()),
39  mesh(sys.get_mesh()),
42  var_num(sys.variable_number(var_name)),
44 {}
NumericVector< Number > & current_local_solution
Definition: dtk_evaluator.h:73
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2144
const FEType & fe_type
Definition: dtk_evaluator.h:79
const MeshBase & get_mesh() const
Definition: system.h:2277
unsigned int variable_number(std::string_view var) const
Definition: system.C:1557
unsigned int var_num
Definition: dtk_evaluator.h:78
unsigned int mesh_dimension() const
Definition: mesh_base.C:324
std::unique_ptr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
Definition: system.h:1585
const DofMap & get_dof_map() const
Definition: system.h:2293
EquationSystems & es
Definition: dtk_evaluator.h:74

Member Function Documentation

◆ evaluate()

DTKEvaluator::FieldContainerType libMesh::DTKEvaluator::evaluate ( const Teuchos::ArrayRCP< int > &  elements,
const Teuchos::ArrayRCP< double > &  coords 
)
overridevirtual

Definition at line 47 of file dtk_evaluator.C.

References libMesh::FEInterface::compute_data(), current_local_solution, dim, libMesh::DofMap::dof_indices(), dof_map, libMesh::MeshBase::elem_ptr(), es, fe_type, libMesh::index_range(), libMesh::FEMap::inverse_map(), mesh, libMesh::FEComputeData::shape, value, and var_num.

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(FEMap::inverse_map(dim, 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 (auto j : index_range(dof_indices))
75  value += current_local_solution(dof_indices[j]) * data.shape[j];
76 
77  values[i] = value;
78  }
79 
80  return evaluations;
81 }
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:1992
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1626
NumericVector< Number > & current_local_solution
Definition: dtk_evaluator.h:73
const FEType & fe_type
Definition: dtk_evaluator.h:79
unsigned int var_num
Definition: dtk_evaluator.h:78
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...
virtual const Elem * elem_ptr(const dof_id_type i) const =0
static const bool value
Definition: xdr_io.C:54
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111
EquationSystems & es
Definition: dtk_evaluator.h:74

Member Data Documentation

◆ current_local_solution

NumericVector<Number>& libMesh::DTKEvaluator::current_local_solution
protected

Definition at line 73 of file dtk_evaluator.h.

Referenced by evaluate().

◆ dim

unsigned int libMesh::DTKEvaluator::dim
protected

Definition at line 76 of file dtk_evaluator.h.

Referenced by evaluate().

◆ dof_map

DofMap& libMesh::DTKEvaluator::dof_map
protected

Definition at line 77 of file dtk_evaluator.h.

Referenced by evaluate().

◆ es

EquationSystems& libMesh::DTKEvaluator::es
protected

Definition at line 74 of file dtk_evaluator.h.

Referenced by evaluate().

◆ fe_type

const FEType& libMesh::DTKEvaluator::fe_type
protected

Definition at line 79 of file dtk_evaluator.h.

Referenced by evaluate().

◆ mesh

MeshBase& libMesh::DTKEvaluator::mesh
protected

Definition at line 75 of file dtk_evaluator.h.

Referenced by evaluate().

◆ sys

System& libMesh::DTKEvaluator::sys
protected

Definition at line 72 of file dtk_evaluator.h.

◆ var_num

unsigned int libMesh::DTKEvaluator::var_num
protected

Definition at line 78 of file dtk_evaluator.h.

Referenced by evaluate().


The documentation for this class was generated from the following files: