libMesh
meshfree_solution_transfer.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 
19 
20 #include "libmesh/meshfree_solution_transfer.h"
21 
22 #include "libmesh/mesh.h"
23 #include "libmesh/system.h"
24 #include "libmesh/numeric_vector.h"
25 #include "libmesh/threads.h"
26 #include "libmesh/meshfree_interpolation.h"
27 #include "libmesh/function_base.h"
28 #include "libmesh/node.h"
29 
30 // C++ includes
31 #include <cstddef>
32 
33 namespace libMesh
34 {
35 
36 // Forward Declarations
37 template <typename T>
38 class DenseVector;
39 
40 // Helper function for doing the projection
41 class MeshlessInterpolationFunction : public FunctionBase<Number>
42 {
43 public:
45  Threads::spin_mutex & mutex) :
46  _mfi(mfi),
47  _mutex(mutex)
48  {}
49 
50  void init () {}
51  void clear () {}
52 
54  {
56  }
57 
59  const Real /*time*/)
60  {
61  _pts.clear();
62  _pts.push_back(p);
63  _vals.resize(1);
64 
65  Threads::spin_mutex::scoped_lock lock(_mutex);
66 
68 
69  return _vals.front();
70  }
71 
72 
73  void operator() (const Point & p,
74  const Real time,
75  DenseVector<Number> & output)
76  {
77  output.resize(1);
78  output(0) = (*this)(p,time);
79  return;
80  }
81 
82 private:
84  mutable std::vector<Point> _pts;
85  mutable std::vector<Number> _vals;
87 };
88 
89 void
91  const Variable & to_var)
92 {
93  libmesh_experimental();
94 
95  System * from_sys = from_var.system();
96  System * to_sys = to_var.system();
97 
98  EquationSystems & from_es = from_sys->get_equation_systems();
99 
100  MeshBase & from_mesh = from_es.get_mesh();
101 
103  (from_mesh.comm(), 4, 2);
104 
105  std::vector<Point> & src_pts (idi.get_source_points());
106  std::vector<Number> & src_vals (idi.get_source_vals());
107 
108  std::vector<std::string> field_vars;
109  field_vars.push_back(from_var.name());
110  idi.set_field_variables(field_vars);
111 
112  // We now will loop over every node in the source mesh
113  // and add it to a source point list, along with the solution
114  for (const auto & node : from_mesh.local_node_ptr_range())
115  {
116  src_pts.push_back(*node);
117  src_vals.push_back((*from_sys->solution)(node->dof_number(from_sys->number(),from_var.number(),0)));
118  }
119 
120  // We have only set local values - prepare for use by gathering remote data
121  idi.prepare_for_use();
122 
123  // Create a MeshlessInterpolationFunction that uses our
124  // InverseDistanceInterpolation object. Since each
125  // MeshlessInterpolationFunction shares the same
126  // InverseDistanceInterpolation object in a threaded environment we
127  // must also provide a locking mechanism.
128  Threads::spin_mutex mutex;
129  MeshlessInterpolationFunction mif(idi, mutex);
130 
131  // project the solution
132  to_sys->project_solution(&mif);
133 }
134 
135 } // namespace libMesh
This is the EquationSystems class.
Inverse distance interpolation.
const std::string & name() const
Definition: variable.h:100
virtual SimpleRange< node_iterator > local_node_ptr_range()=0
const std::vector< std::string > & field_variables() const
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
void init()
The actual initialization process.
void project_solution(FunctionBase< Number > *f, FunctionBase< Gradient > *g=libmesh_nullptr) const
Projects arbitrary functions onto the current solution.
The libMesh namespace provides an interface to certain functionality in the library.
This is the MeshBase class.
Definition: mesh_base.h:68
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
Number operator()(const Point &p, const Real time=0.)
unsigned int number() const
Definition: variable.h:106
virtual void transfer(const Variable &from_var, const Variable &to_var) libmesh_override
Transfer the values of a variable to another.
This class defines the notion of a variable in the system.
Definition: variable.h:49
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
Base class to support various mesh-free interpolation methods.
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
MeshlessInterpolationFunction(const MeshfreeInterpolation &mfi, Threads::spin_mutex &mutex)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned int number() const
Definition: system.h:2006
const Parallel::Communicator & comm() const
System * system() const
Definition: variable.h:92
void set_field_variables(const std::vector< std::string > &names)
Defines the field variable(s) we are responsible for, and importantly their assumed ordering...
const EquationSystems & get_equation_systems() const
Definition: system.h:712
const MeshBase & get_mesh() const
virtual void interpolate_field_data(const std::vector< std::string > &field_names, const std::vector< Point > &tgt_pts, std::vector< Number > &tgt_vals) const =0
Interpolate source data at target points.
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
virtual UniquePtr< FunctionBase< Number > > clone() const