libMesh
Functions
miscellaneous_ex8.C File Reference

Go to the source code of this file.

Functions

void create_random_point_cloud (const unsigned int Npts, std::vector< Point > &pts, const Real max_range=10)
 
Real exact_solution_u (const Point &p)
 
Real exact_solution_v (const Point &p)
 
Number exact_value (const Point &p, const Parameters &, const std::string &, const std::string &)
 
void init_sys (EquationSystems &es, const std::string &system_name)
 
int main (int argc, char **argv)
 

Function Documentation

void create_random_point_cloud ( const unsigned int  Npts,
std::vector< Point > &  pts,
const Real  max_range = 10 
)

Definition at line 47 of file miscellaneous_ex8.C.

References libMesh::out, and libMesh::Real.

Referenced by main().

50 {
51  libMesh::out << "Generating "<< Npts << " point cloud...";
52  pts.resize(Npts);
53 
54  for (size_t i=0;i<Npts;i++)
55  {
56  pts[i](0) = max_range * (std::rand() % 1000) / Real(1000);
57  pts[i](1) = max_range * (std::rand() % 1000) / Real(1000);
58  pts[i](2) = max_range * (std::rand() % 1000) / Real(1000);
59  }
60  libMesh::out << "done\n";
61 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OStreamProxy out
Real exact_solution_u ( const Point p)

Definition at line 65 of file miscellaneous_ex8.C.

References libMesh::Real, and libMesh::x.

Referenced by main().

66 {
67  const Real
68  x = p(0),
69  y = p(1),
70  z = p(2);
71 
72  return (x*x*x +
73  y*y*y*y +
74  z*z*z*z*z);
75 }
PetscErrorCode Vec x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real exact_solution_v ( const Point p)

Definition at line 79 of file miscellaneous_ex8.C.

References libMesh::Real, and libMesh::x.

Referenced by exact_value(), and main().

80 {
81  const Real
82  x = p(0),
83  y = p(1),
84  z = p(2);
85 
86  return (x*x +
87  y*y +
88  z*z*z);
89 }
PetscErrorCode Vec x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Number exact_value ( const Point p,
const Parameters ,
const std::string &  ,
const std::string &   
)

Definition at line 91 of file miscellaneous_ex8.C.

References exact_solution_v().

Referenced by init_sys().

95 {
96  return exact_solution_v(p);
97 }
Real exact_solution_v(const Point &p)
void init_sys ( EquationSystems es,
const std::string &  system_name 
)

Definition at line 103 of file miscellaneous_ex8.C.

References exact_value(), libMesh::EquationSystems::get_system(), libmesh_nullptr, libMesh::EquationSystems::parameters, and libMesh::System::project_solution().

Referenced by main().

105 {
106  // Get a reference to the Convection-Diffusion system object.
107  System & system =
108  es.get_system<System>(system_name);
109 
111 }
const class libmesh_nullptr_t libmesh_nullptr
void project_solution(FunctionBase< Number > *f, FunctionBase< Gradient > *g=libmesh_nullptr) const
Projects arbitrary functions onto the current solution.
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
Parameters parameters
Data structure holding arbitrary parameters.
const T_sys & get_system(const std::string &name) const
Number exact_value(const Point &p, const Parameters &, const std::string &, const std::string &)
int main ( int  argc,
char **  argv 
)

Definition at line 116 of file miscellaneous_ex8.C.

References libMesh::EquationSystems::add_system(), libMesh::System::add_variable(), libMesh::System::attach_init_function(), libMesh::LibMeshInit::comm(), create_random_point_cloud(), libMesh::System::current_solution(), exact_solution_u(), exact_solution_v(), libMesh::FIRST, libMesh::TriangleWrapper::init(), libMesh::EquationSystems::init(), init_sys(), libMesh::out, libMesh::UnstructuredMesh::read(), libMesh::MeshfreeInterpolation::set_field_variables(), and libMesh::MeshOutput< MT >::write_equation_systems().

117 {
118  // Skip this example if we do not meet certain requirements
119  libmesh_example_requires(3 <= LIBMESH_DIM, "3D support");
120 #ifndef LIBMESH_HAVE_EIGEN
121  libmesh_example_requires(false, "--enable-eigen");
122 #endif
123 #ifndef LIBMESH_HAVE_ZLIB_H
124  libmesh_example_requires(false, "--enable-zlib");
125 #endif
126 
127  // Initialize libMesh.
128  LibMeshInit init (argc, argv);
129  {
130  // Demonstration case 1
131  {
132  std::vector<Point> tgt_pts;
133  std::vector<Number> tgt_data_idi, tgt_data_rbi;
134  std::vector<std::string> field_vars;
135 
136  field_vars.push_back("u");
137  field_vars.push_back("v");
138 
140  /* n_interp_pts = */ 8,
141  /* power = */ 2);
142 
143  RadialBasisInterpolation<3> rbi (init.comm());
144 
145  idi.set_field_variables (field_vars);
146  rbi.set_field_variables (field_vars);
147 
149  idi.get_source_points());
150 
151 
152  // Explicitly set the data values we will interpolate from
153  {
154  const std::vector<Point> & src_pts (idi.get_source_points());
155  std::vector<Number> & src_vals (idi.get_source_vals());
156 
157  src_vals.clear(); src_vals.reserve(2*src_pts.size());
158 
159  for (std::vector<Point>::const_iterator pt_it=src_pts.begin();
160  pt_it != src_pts.end(); ++pt_it)
161  {
162  src_vals.push_back (exact_solution_u (*pt_it));
163  src_vals.push_back (exact_solution_v (*pt_it));
164  }
165  }
166 
167  // give rbi the same info as idi
168  rbi.get_source_points() = idi.get_source_points();
169  rbi.get_source_vals() = idi.get_source_vals();
170 
171  idi.prepare_for_use();
172  rbi.prepare_for_use();
173 
174  libMesh::out << idi;
175 
176  // Interpolate to some other random points, and evaluate the result
177  {
179  tgt_pts);
180 
181  //tgt_pts = rbi.get_source_points();
182 
183  idi.interpolate_field_data (field_vars,
184  tgt_pts,
185  tgt_data_idi);
186 
187  rbi.interpolate_field_data (field_vars,
188  tgt_pts,
189  tgt_data_rbi);
190 
191  std::vector<Number>::const_iterator
192  v_idi = tgt_data_idi.begin(),
193  v_rbi = tgt_data_rbi.begin();
194 
195  for (std::vector<Point>::const_iterator p_it=tgt_pts.begin();
196  p_it!=tgt_pts.end(); ++p_it)
197  {
198  libMesh::out << "\nAt target point " << *p_it
199  << "\n u_interp_idi=" << *v_idi
200  << ", u_interp_rbi=" << *v_rbi
201  << ", u_exact=" << exact_solution_u(*p_it);
202  ++v_idi;
203  ++v_rbi;
204  libMesh::out << "\n v_interp_idi=" << *v_idi
205  << ", v_interp_rbi=" << *v_rbi
206  << ", v_exact=" << exact_solution_v(*p_it)
207  << std::endl;
208  ++v_idi;
209  ++v_rbi;
210  }
211  }
212  }
213 
214 
215  // Demonstration case 2
216  {
217  Mesh mesh_a(init.comm()), mesh_b(init.comm());
218 
219  mesh_a.read("struct.ucd.gz");
220  mesh_b.read("unstruct.ucd.gz");
221 
222  // Create equation systems objects.
224  es_a(mesh_a), es_b(mesh_b);
225 
226  System & sys_a = es_a.add_system<System>("src_system");
227  System & sys_b = es_b.add_system<System>("dest_system");
228 
229  sys_a.add_variable ("Cp", FIRST);
230  sys_b.add_variable ("Cp", FIRST);
231 
233  es_a.init();
234 
235  // Write out the initial conditions.
236  TecplotIO(mesh_a).write_equation_systems ("src.dat",
237  es_a);
238 
240  /* n_interp_pts = */ 4,
241  /* power = */ 2);
242  RadialBasisInterpolation<3> rbi (init.comm());
243 
244  std::vector<Point> & src_pts (idi.get_source_points());
245  std::vector<Number> & src_vals (idi.get_source_vals());
246  std::vector<std::string> field_vars;
247  field_vars.push_back("Cp");
248  idi.set_field_variables(field_vars);
249 
250  // We now will loop over every node in the source mesh
251  // and add it to a source point list, along with the solution
252  for (const auto & node : mesh_a.local_node_ptr_range())
253  {
254  src_pts.push_back(*node);
255  src_vals.push_back(sys_a.current_solution(node->dof_number(0, 0, 0)));
256  }
257 
258  rbi.set_field_variables(field_vars);
259  rbi.get_source_points() = idi.get_source_points();
260  rbi.get_source_vals() = idi.get_source_vals();
261 
262  // We have only set local values - prepare for use by gathering remote data
263  idi.prepare_for_use();
264  rbi.prepare_for_use();
265 
266  // Create a MeshlessInterpolationFunction that uses our InverseDistanceInterpolation
267  // object. Since each MeshlessInterpolationFunction shares the same InverseDistanceInterpolation
268  // object in a threaded environment we must also provide a locking mechanism.
269  {
270  Threads::spin_mutex mutex;
271  MeshlessInterpolationFunction mif(idi, mutex);
272 
273  // project the solution onto system b
274  es_b.init();
275  sys_b.project_solution (&mif);
276 
277  // Write the result
278  TecplotIO(mesh_b).write_equation_systems ("dest_idi.dat",
279  es_b);
280  }
281 
282  // Create a MeshlessInterpolationFunction that uses our RadialBasisInterpolation
283  // object. Since each MeshlessInterpolationFunction shares the same RadialBasisInterpolation
284  // object in a threaded environment we must also provide a locking mechanism.
285  {
286  Threads::spin_mutex mutex;
287  MeshlessInterpolationFunction mif(rbi, mutex);
288 
289  // project the solution onto system b
290  sys_b.project_solution (&mif);
291 
292  // Write the result
293  TecplotIO(mesh_b).write_equation_systems ("dest_rbi.dat",
294  es_b);
295  }
296  }
297  }
298  return 0;
299 }
This is the EquationSystems class.
Inverse distance interpolation.
virtual void read(const std::string &name, void *mesh_data=libmesh_nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false) libmesh_override
Reads the file specified by name.
Radial Basis Function interpolation.
unsigned int add_variable(const std::string &var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=libmesh_nullptr)
Adds the variable var to the list of variables for this system.
Definition: system.C:1101
Real exact_solution_v(const Point &p)
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
void init_sys(EquationSystems &es, const std::string &system_name)
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
void attach_init_function(void fptr(EquationSystems &es, const std::string &name))
Register a user function to use in initializing the system.
Definition: system.C:1760
virtual void write_equation_systems(const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
This method implements writing a mesh with data to a specified file where the data is taken from the ...
Definition: mesh_output.C:31
Real exact_solution_u(const Point &p)
Number current_solution(const dof_id_type global_dof_number) const
Definition: system.C:192
OStreamProxy out
void set_field_variables(const std::vector< std::string > &names)
Defines the field variable(s) we are responsible for, and importantly their assumed ordering...
This class implements writing meshes in the Tecplot format.
Definition: tecplot_io.h:43
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
void create_random_point_cloud(const unsigned int Npts, std::vector< Point > &pts, const Real max_range=10)