19 #include "libmesh/libmesh_config.h" 25 #ifdef LIBMESH_HAVE_GETOPT_H 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" 48 void usage(
const std::string & progName)
50 std::ostringstream helpList;
51 helpList <<
"usage:\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" 65 <<
" -q really quiet\n" 66 <<
" -h Print help menu\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" 73 <<
" ./compare -d 3 -m grid.xda -l leftfile.dat -r rightfile.dat -b -t 1.e-8\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" 83 libmesh_error_msg(helpList.str());
90 std::vector<std::string> & names,
102 bool format_set =
false;
103 bool left_name_set =
false;
106 usage(std::string(argv[0]));
109 while ((opt = getopt(argc, argv, optionStr)) != -1)
119 libmesh_error_msg_if(!names.empty(),
"ERROR: Mesh file name must precede left file name!");
121 names.push_back(optarg);
130 dim = cast_int<unsigned char>(atoi(optarg));
139 libmesh_error_msg_if(left_name_set,
"ERROR: Mesh file name must precede right file name!");
141 names.push_back(optarg);
142 left_name_set =
true;
151 libmesh_error_msg_if(names.empty() || !left_name_set,
152 "ERROR: Mesh file name and left file name must precede right file name!");
154 names.push_back(optarg);
163 threshold = atof(optarg);
172 libmesh_error_msg_if(format_set,
"ERROR: Equation system file format already set!");
184 libmesh_error_msg_if(format_set,
"ERROR: Equation system file format already set!");
238 libMesh::out <<
"********* LEFT SYSTEM *********" << std::endl;
240 libMesh::out <<
"********* RIGHT SYSTEM *********" << std::endl;
242 libMesh::out <<
"********* COMPARISON PHASE *********" << std::endl
249 bool result = les.
compare(res, threshold, verbose);
252 libMesh::out <<
"********* FINISHED *********" << std::endl;
266 int main (
int argc,
char ** argv)
277 std::vector<std::string> names;
278 unsigned char dim =
static_cast<unsigned char>(-1);
281 bool verbose =
false;
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" 296 <<
" -d 3 ... for example");
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
323 left_mesh.
read (names[0]);
324 right_mesh.read (names[0]);
327 left_mesh.print_info();
342 libmesh_error_msg_if(names.size() != 3,
"Bad input specified.");
344 left_system.
read (names[1], format);
345 right_system.
read (names[2], format);
347 are_equal =
do_compare (left_system, right_system, threshold, verbose);
353 unsigned int our_result;
359 <<
" Congrat's, up to the defined threshold, the two" 369 <<
" Oops, differences occurred!" 371 <<
" Use -v to obtain more information where differences occurred." This is the EquationSystems class.
static constexpr Real TOLERANCE
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.
The libMesh namespace provides an interface to certain functionality in the library.
XdrMODE
Defines an enum for read/write mode in Xdr format.
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().
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
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.
int main(int argc, char **argv)
void process_cmd_line(int argc, char **argv, std::vector< std::string > &names, unsigned char &dim, double &threshold, XdrMODE &format, bool &verbose, bool &quiet)
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.