libMesh
compare.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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 #include "libmesh/libmesh_config.h"
20 
21 // C++ includes
22 #include <iostream>
23 #include <vector>
24 #include <string>
25 #ifdef LIBMESH_HAVE_GETOPT_H
26 // GCC 2.95.3 (and maybe others) do not include
27 // getopt.h in unistd.h... However IBM xlC has no
28 // getopt.h! This works around that.
29 #include <getopt.h>
30 #endif
31 #include <stdio.h>
32 #include <fstream>
33 
34 // Local Includes
35 #include "libmesh/libmesh.h"
36 #include "libmesh/equation_systems.h"
37 #include "libmesh/mesh.h"
38 #include "libmesh/perfmon.h"
39 #include "libmesh/enum_xdr_mode.h"
40 
41 
42 using namespace libMesh;
43 
44 
48 void usage(const std::string & progName)
49 {
50  std::ostringstream helpList;
51  helpList << "usage:\n"
52  << " "
53  << progName
54  << " [options] ...\n"
55  << "\n"
56  << "options:\n"
57  << " -d <dim> <dim>-dimensional mesh\n"
58  << " -m <string> Mesh file name\n"
59  << " -l <string> Left Equation Systems file name\n"
60  << " -r <string> Right Equation Systems file name\n"
61  << " -t <float> threshold\n"
62  << " -a ASCII format (default)\n"
63  << " -b binary format\n"
64  << " -v Verbose\n"
65  << " -q really quiet\n"
66  << " -h Print help menu\n"
67  << "\n"
68  << "\n"
69  << " This program is used to compare equation systems to a user-specified\n"
70  << " threshold. Equation systems are imported in the libMesh format\n"
71  << " provided through the read and write methods in class EquationSystems.\n"
72  << " \n"
73  << " ./compare -d 3 -m grid.xda -l leftfile.dat -r rightfile.dat -b -t 1.e-8\n"
74  << "\n"
75  << " will read in the mesh grid.xda, the equation systems leftfile.dat and\n"
76  << " rightfile.dat in binary format and compare systems, and especially the\n"
77  << " floats stored in vectors. The comparison is said to be passed when the\n"
78  << " floating point values agree up to the given threshold. When no threshold\n"
79  << " is set the default libMesh tolerance is used. If neither -a or -b are set,\n"
80  << " ASCII format is assumed.\n"
81  << "\n";
82 
83  libmesh_error_msg(helpList.str());
84 }
85 
86 
87 
88 void process_cmd_line(int argc,
89  char ** argv,
90  std::vector<std::string> & names,
91  unsigned char & dim,
92  double & threshold,
93  XdrMODE & format,
94  bool & verbose,
95  bool & quiet)
96 {
97  char optionStr[] =
98  "d:m:l:r:t:abvq?h";
99 
100  int opt;
101 
102  bool format_set = false;
103  bool left_name_set = false;
104 
105  if (argc < 3)
106  usage(std::string(argv[0]));
107 
108 
109  while ((opt = getopt(argc, argv, optionStr)) != -1)
110  {
111  switch (opt)
112  {
113 
117  case 'm':
118  {
119  libmesh_error_msg_if(!names.empty(), "ERROR: Mesh file name must precede left file name!");
120 
121  names.push_back(optarg);
122  break;
123  }
124 
128  case 'd':
129  {
130  dim = cast_int<unsigned char>(atoi(optarg));
131  break;
132  }
133 
137  case 'l':
138  {
139  libmesh_error_msg_if(left_name_set, "ERROR: Mesh file name must precede right file name!");
140 
141  names.push_back(optarg);
142  left_name_set = true;
143  break;
144  }
145 
149  case 'r':
150  {
151  libmesh_error_msg_if(names.empty() || !left_name_set,
152  "ERROR: Mesh file name and left file name must precede right file name!");
153 
154  names.push_back(optarg);
155  break;
156  }
157 
161  case 't':
162  {
163  threshold = atof(optarg);
164  break;
165  }
166 
170  case 'a':
171  {
172  libmesh_error_msg_if(format_set, "ERROR: Equation system file format already set!");
173 
174  format = READ;
175  format_set = true;
176  break;
177  }
178 
182  case 'b':
183  {
184  libmesh_error_msg_if(format_set, "ERROR: Equation system file format already set!");
185 
186  format = DECODE;
187  format_set = true;
188  break;
189  }
190 
191 
195  case 'v':
196  {
197  verbose = true;
198  break;
199  }
200 
204  case 'q':
205  {
206  quiet = true;
207  break;
208  }
209 
210  case 'h':
211  case '?':
212  usage(argv[0]);
213 
214  default:
215  return;
216  }
217  }
218 
219 }
220 
221 
222 
223 
224 
231  EquationSystems & res,
232  double threshold,
233  bool verbose)
234 {
235 
236  if (verbose)
237  {
238  libMesh::out << "********* LEFT SYSTEM *********" << std::endl;
239  les.print_info ();
240  libMesh::out << "********* RIGHT SYSTEM *********" << std::endl;
241  res.print_info ();
242  libMesh::out << "********* COMPARISON PHASE *********" << std::endl
243  << std::endl;
244  }
245 
249  bool result = les.compare(res, threshold, verbose);
250  if (verbose)
251  {
252  libMesh::out << "********* FINISHED *********" << std::endl;
253  }
254  return result;
255 }
256 
257 
258 
259 
260 
261 
262 
263 
264 
265 
266 int main (int argc, char ** argv)
267 {
268  LibMeshInit init(argc, argv);
269 
270  // these should not be contained in the following braces
271  bool quiet = false;
272  bool are_equal;
273 
274  PerfMon perfmon(argv[0]);
275 
276  // default values
277  std::vector<std::string> names;
278  unsigned char dim = static_cast<unsigned char>(-1);
279  double threshold = double(TOLERANCE);
280  XdrMODE format = READ;
281  bool verbose = false;
282 
283  // get commands
284  process_cmd_line(argc, argv,
285  names,
286  dim,
287  threshold,
288  format,
289  verbose,
290  quiet);
291 
292  libmesh_error_msg_if(dim == static_cast<unsigned char>(-1),
293  "ERROR: you must specify the dimension on "
294  << "the command line!\n\n"
295  << argv[0]
296  << " -d 3 ... for example");
297 
298  if (quiet)
299  verbose = false;
300 
301  if (verbose)
302  {
303  libMesh::out << "Settings:" << std::endl
304  << " dimensionality = " << +dim << std::endl
305  << " mesh = " << names[0] << std::endl
306  << " left system = " << names[1] << std::endl
307  << " right system = " << names[2] << std::endl
308  << " threshold = " << threshold << std::endl
309  << " read format = " << format << std::endl
310  << std::endl;
311  }
312 
313 
317  Mesh left_mesh (init.comm(), dim);
318  Mesh right_mesh (init.comm(), dim);
319 
320 
321  if (!names.empty())
322  {
323  left_mesh.read (names[0]);
324  right_mesh.read (names[0]);
325 
326  if (verbose)
327  left_mesh.print_info();
328  }
329  else
330  {
331  libMesh::out << "No input specified." << std::endl;
332  return 1;
333  }
334 
335 
339  EquationSystems left_system (left_mesh);
340  EquationSystems right_system (right_mesh);
341 
342  libmesh_error_msg_if(names.size() != 3, "Bad input specified.");
343 
344  left_system.read (names[1], format);
345  right_system.read (names[2], format);
346 
347  are_equal = do_compare (left_system, right_system, threshold, verbose);
348 
349 
353  unsigned int our_result;
354 
355  if (are_equal)
356  {
357  if (!quiet)
358  libMesh::out << std::endl
359  << " Congrat's, up to the defined threshold, the two"
360  << std::endl
361  << " are identical."
362  << std::endl;
363  our_result=0;
364  }
365  else
366  {
367  if (!quiet)
368  libMesh::out << std::endl
369  << " Oops, differences occurred!"
370  << std::endl
371  << " Use -v to obtain more information where differences occurred."
372  << std::endl;
373  our_result=1;
374  }
375 
376  // return libMesh::close();
377  return our_result;
378 }
This is the EquationSystems class.
static constexpr Real TOLERANCE
unsigned int dim
void print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:90
The libMesh namespace provides an interface to certain functionality in the library.
XdrMODE
Defines an enum for read/write mode in Xdr format.
Definition: enum_xdr_mode.h:35
bool do_compare(EquationSystems &les, EquationSystems &res, double threshold, bool verbose)
everything that is identical for the systems, and should not go into EquationSystems::compare(), can go in this do_compare().
Definition: compare.C:230
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
void usage(const std::string &progName)
how to use this, and command line processor
Definition: compare.C:48
void read(std::string_view name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
Read & initialize the systems from disk using the XDR data format.
PAPI stands for Performance Application Programming Interface.
Definition: perfmon.h:57
int main(int argc, char **argv)
Definition: compare.C:266
OStreamProxy out
void process_cmd_line(int argc, char **argv, std::vector< std::string > &names, unsigned char &dim, double &threshold, XdrMODE &format, bool &verbose, bool &quiet)
Definition: compare.C:88
virtual bool compare(const EquationSystems &other_es, const Real threshold, const bool verbose) const
virtual void read(const std::string &name, void *mesh_data=nullptr, bool skip_renumber_nodes_and_elements=false, bool skip_find_neighbors=false) override
Reads the file specified by name.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50