www.mooseframework.org
DTKInterpolationEvaluator.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 #include "libmesh/libmesh_config.h"
15 
16 #ifdef LIBMESH_TRILINOS_HAVE_DTK
17 
18 // DTK includes
20 #include "DTKInterpolationHelper.h"
21 
22 #include "libmesh/dof_map.h"
23 #include "libmesh/fe_interface.h"
24 #include "libmesh/fe_compute_data.h"
25 #include "libmesh/numeric_vector.h"
26 #include "libmesh/equation_systems.h"
27 #include "libmesh/mesh.h"
28 #include "libmesh/system.h"
29 
30 namespace libMesh
31 {
32 
34  std::string var_name,
35  const Point & offset)
36  : sys(in_sys),
37  _offset(offset),
38  current_local_solution(*sys.current_local_solution),
39  es(in_sys.get_equation_systems()),
40  mesh(sys.get_mesh()),
41  dim(mesh.mesh_dimension()),
42  dof_map(sys.get_dof_map()),
43  var_num(sys.variable_number(var_name)),
44  fe_type(dof_map.variable_type(var_num))
45 {
46 }
47 
49 DTKInterpolationEvaluator::evaluate(const Teuchos::ArrayRCP<GlobalOrdinal> & elements,
50  const Teuchos::ArrayRCP<double> & coords)
51 {
52  GlobalOrdinal num_values = elements.size();
53 
54  Teuchos::ArrayRCP<Number> values(num_values);
55  DataTransferKit::FieldContainer<Number> evaluations(values, 1);
56 
57  for (GlobalOrdinal i = 0; i < num_values; i++)
58  {
59  Elem * elem = mesh.elem(elements[i]);
60 
61  Point p;
62 
63  for (unsigned int j = 0; j < dim; j++)
64  p(j) = coords[(j * num_values) + i] - _offset(j);
65 
66  const Point mapped_point(FEInterface::inverse_map(dim, dof_map.variable_type(0), elem, p));
67 
68  FEComputeData data(es, mapped_point);
69  FEInterface::compute_data(dim, fe_type, elem, data);
70 
71  std::vector<dof_id_type> dof_indices;
72  dof_map.dof_indices(elem, dof_indices, var_num);
73 
74  Number value = 0;
75 
76  for (unsigned int j = 0; j < dof_indices.size(); j++)
77  value += current_local_solution(dof_indices[j]) * data.shape[j];
78 
79  values[i] = value;
80  }
81 
82  return evaluations;
83 }
84 
85 } // namespace libMesh
86 
87 #endif // LIBMESH_TRILINOS_HAVE_DTK
NumericVector< Number > & current_local_solution
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
DofMap & dof_map
DataTransferKit::FieldContainer< Number > FieldContainerType
DataTransferKit::MeshTraits< MeshContainerType >::global_ordinal_type GlobalOrdinal
FieldContainerType evaluate(const Teuchos::ArrayRCP< GlobalOrdinal > &elements, const Teuchos::ArrayRCP< double > &coords)
DTKInterpolationEvaluator(System &in_sys, std::string var_name, const Point &offset)