libMesh
meshtool.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 #include "libmesh/libmesh_config.h"
19 
20 // C++ includes
21 #include <algorithm>
22 #include <fstream>
23 #ifdef LIBMESH_HAVE_GETOPT_H
24 // GCC 2.95.3 (and maybe others) do not include
25 // getopt.h in unistd.h... However IBM xlC has no
26 // getopt.h! This works around that.
27 #include <getopt.h>
28 #endif
29 #include <iostream>
30 #include <math.h>
31 #include <stdio.h>
32 #include <stdlib.h>
33 #include <string.h>
34 #include <string>
35 #include <unistd.h>
36 #include <vector>
37 
38 // Local Includes
39 #include "libmesh/boundary_info.h"
40 #include "libmesh/boundary_mesh.h"
41 #include "libmesh/dof_map.h"
42 #include "libmesh/elem.h"
43 #include "libmesh/elem_quality.h"
44 #include "libmesh/gmv_io.h"
45 #include "libmesh/inf_elem_builder.h"
46 #include "libmesh/libmesh.h"
47 #include "libmesh/mesh.h"
48 #include "libmesh/mesh_modification.h"
49 #include "libmesh/mesh_refinement.h"
50 #include "libmesh/perfmon.h"
51 #include "libmesh/statistics.h"
52 #include "libmesh/string_to_enum.h"
53 
54 
55 using namespace libMesh;
56 
57 
58 /*
59  * convenient enum for the mode in which the boundary mesh
60  * should be written
61  */
63 
64 
65 void usage(const std::string & progName)
66 {
67  std::ostringstream helpList;
68  helpList << "usage:\n"
69  << " "
70  << progName
71  << " [options] ...\n"
72  << "\n"
73  << "options:\n"
74  << " -d <dim> <dim>-dimensional mesh\n"
75  << " -i <string> Input file name\n"
76  << " -o <string> Output file name\n"
77  << " -s <string> Solution file name\n"
78  << "\n -b Write the boundary conditions\n"
79  << " -D <factor> Randomly move interior nodes by D*hmin\n"
80  << " -h Print help menu\n"
81  << " -p <count> Partition into <count> subdomains\n"
82 #ifdef LIBMESH_ENABLE_AMR
83  << " -r <count> Globally refine <count> times\n"
84 #endif
85  << " -t (-d 2 only) Convert to triangles first\n"
86  << " (allows to write .unv file of the\n"
87  << " boundary with the correct node ids)\n"
88  << " -v Verbose\n"
89  << " -q <metric> Evaluates the named element quality metric\n"
90  << " -1 Converts a mesh of higher order elements\n"
91  << " to their first-order counterparts:\n"
92  << " Quad8 -> Quad4, Tet10 -> Tet4 etc\n"
93  << " -2 Converts a mesh of linear elements\n"
94  << " to their second-order counterparts:\n"
95  << " Quad4 -> Quad8, Tet4 -> Tet10 etc\n"
96  << " -3 Same, but to the highest possible:\n"
97  << " Quad4 -> Quad9, Hex8 -> Hex27 etc\n"
98 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
99  << "\n -a Add infinite elements\n"
100  << " -x <coord> Specify infinite element origin\n"
101  << " -y <coord> coordinates. If none given, origin\n"
102  << " -z <coord> is determined automatically.\n"
103  << " -X When building infinite elements \n"
104  << " -Y treat mesh as x/y/z-symmetric.\n"
105  << " -Z When -X is given, -x <coord> also\n"
106  << " has to be given. Similar for y,z.\n"
107 #endif
108  // << "\n -l Build the L connectivity matrix \n"
109  // << " -L Build the script L connectivity matrix \n"
110  << "\n"
111  << "\n"
112  << " This program is used to convert and partition from/to a variety of\n"
113  << " formats. File types are inferred from file extensions. For example,\n"
114  << " the command:\n"
115  << "\n"
116  << " ./meshtool -d 2 -i in.e -o out.plt\n"
117  << "\n"
118  << " will read a 2D mesh in the ExodusII format (from Cubit, for example)\n"
119  << " from the file in.e. It will then write the mesh in the Tecplot\n"
120  << " binary format to out.plt.\n"
121  << "\n"
122  << " and\n"
123  << "\n"
124  << " ./meshtool -d 3 -i bench12.mesh.0000 -o out.gmv -s bench12.soln.0137\n"
125  << "\n"
126  << " will read a 3D MGF mesh from the file bench12.mesh.0000, read a\n"
127  << " solution from bench12.soln.0137, and write the output in GMV format\n"
128  << " to out.gmv\n"
129 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
130  << "\n"
131  << " and\n"
132  << "\n"
133  << " ./meshtool -d 3 -i dry.unv -o packed.gmv -a -x 30.5 -y -10.5 -X\n"
134  << "\n"
135  << " will read a 3D Universal file, determine the z-coordinate of the origin\n"
136  << " automatically, e.g. z_origin = 3., build infinite elements with the\n"
137  << " origin (30.5, -10.5, 3.) on top of volume elements, while preserving\n"
138  << " a symmetry plane through (30.5, -10.5, 3.) perpendicular to x.\n"
139  << " It is imperative that the origin lies _inside_ the given volume mesh.\n"
140  << " If not, infinite elements are not correctly built!\n"
141 #endif
142  << "\n"
143  << " Currently this program supports the following formats:\n"
144  << "\n"
145  << "INPUT:\n"
146  << " .e -- Sandia's ExodusII binary grid format\n"
147  << " .ucd -- AVS unstructured ASCII grid format\n"
148  << " .unv -- SDRC I-Deas Universal File ASCII format\n"
149  << " .xda -- libMesh human-readable ASCII format\n"
150  << " .xdr -- libMesh binary format\n"
151  << "\n"
152  << "OUTPUT:\n"
153  << " .dat -- Tecplot ASCII format\n"
154  << " .e -- Sandia's ExodusII format\n"
155  << " .exd -- Sandia's ExodusII format\n"
156  << " .fro -- ACDL's .fro format\n"
157  << " .gmv -- LANL's General Mesh Viewer format\n"
158  << " .mesh -- MEdit mesh format\n"
159  << " .msh -- GMSH ASCII file\n"
160  << " .plt -- Tecplot binary format\n"
161  << " .poly -- TetGen ASCII file\n"
162  << " .pvtu -- Paraview VTK format\n"
163  << " .ucd -- AVS's ASCII UCD format\n"
164  << " .unv -- I-deas Universal format\n"
165  << " .xda -- libMesh ASCII format\n"
166  << " .xdr -- libMesh binary format\n"
167  << " .gz -- any above format gzipped\n"
168  << " .bz2 -- any above format bzip2'ed\n"
169  << "\n";
170 
171  libmesh_error_msg(helpList.str());
172 }
173 
174 
175 
176 void process_cmd_line(int argc,
177  char ** argv,
178  std::vector<std::string> & names,
179  unsigned int & n_subdomains,
180  unsigned int & n_rsteps,
181  unsigned char & dim,
182  double & dist_fact,
183  bool & verbose,
184  BoundaryMeshWriteMode & write_bndry,
185  bool & convert_first_order,
186  unsigned int & convert_second_order,
187  bool & triangulate,
188  bool & do_quality,
189  ElemQuality & quality_type,
190  bool & addinfelems,
191 
192 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
196 #endif
197 
198  bool & x_sym,
199  bool & y_sym,
200  bool & z_sym
201  )
202 {
203 
204 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
205 
206  /*
207  * initialize these to some values,
208  * so that the compiler does not complain
209  */
210  triangulate = false;
211  do_quality = false;
212  quality_type = DIAGONAL;
213  addinfelems = false;
214  x_sym = y_sym = z_sym = false;
215 
216  char optionStr[] =
217  "i:o:s:d:D:r:p:tbB123vlLm?h";
218 
219 #else
220 
221  char optionStr[] =
222  "i:o:q:s:d:D:r:p:tbB123a::x:y:z:XYZvlLm?h";
223 
224 #endif
225 
226  bool b_mesh_B_given = false;
227 
228  int opt;
229 
230  if (argc < 3)
231  usage(std::string(argv[0]));
232 
233 
234 
235  while ((opt = getopt(argc, argv, optionStr)) != -1)
236  {
237  switch (opt)
238  {
239 
243  case 'i':
244  {
245  if (names.empty())
246  names.push_back(optarg);
247  else
248  libmesh_error_msg("ERROR: Input name must precede output name!");
249  break;
250  }
251 
255  case 'o':
256  {
257  if (!names.empty())
258  names.push_back(optarg);
259  else
260  libmesh_error_msg("ERROR: Input name must precede output name!");
261  break;
262  }
263 
267  case 's':
268  {
269  if (names.size() == 2)
270  names.push_back(optarg);
271  else
272  libmesh_error_msg("ERROR: Input and output names must precede solution name!");
273  break;
274  }
275 
276 
280  case 'd':
281  {
282  dim = cast_int<unsigned char>(atoi(optarg));
283  break;
284  }
285 
289  case 'D':
290  {
291  dist_fact = atof(optarg);
292  break;
293  }
294 
298  case 'r':
299  {
300  n_rsteps = atoi(optarg);
301  break;
302  }
303 
307  case 'p':
308  {
309  n_subdomains = atoi(optarg);
310  break;
311  }
312 
316  case 't':
317  {
318  triangulate = true;
319  break;
320  }
321 
325  case 'q':
326  {
327  do_quality = true;
328  quality_type = Utility::string_to_enum<ElemQuality>(optarg);
329  break;
330  }
331 
335  case 'v':
336  {
337  verbose = true;
338  break;
339  }
340 
344  case 'b':
345  {
346  if (b_mesh_B_given)
347  libmesh_error_msg("ERROR: Do not use -b and -B concurrently!");
348 
349  write_bndry = BM_MESH_ONLY;
350  break;
351  }
352 
357  case '1':
358  {
359  convert_first_order = true;
360  break;
361  }
362 
363 
368  case '2':
369  {
370  convert_second_order = 2;
371  break;
372  }
373 
378  case '3':
379  {
380  convert_second_order = 22;
381  break;
382  }
383 
384 
385 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
386 
390  case 'a':
391  {
392  addinfelems = true;
393  break;
394  }
395 
399  case 'x':
400  {
401  origin_x.first = true;
402  origin_x.second = atof(optarg);
403  break;
404  }
405 
406  case 'y':
407  {
408  origin_y.first = true;
409  origin_y.second = atof(optarg);
410  break;
411  }
412 
413  case 'z':
414  {
415  origin_z.first = true;
416  origin_z.second = atof(optarg);
417  break;
418  }
419 
423  case 'X':
424  {
425  x_sym = true;
426  break;
427  }
428  case 'Y':
429  {
430  y_sym = true;
431  break;
432  }
433  case 'Z':
434  {
435  z_sym = true;
436  break;
437  }
438 
439 #endif //ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
440 
441 
442 
443  case 'h':
444  case '?':
445  usage(argv[0]);
446 
447  default:
448  return;
449  };
450  };
451 
452 }
453 
454 
455 
456 // A helper function for creating a submesh of active elements. Why do we need this?
457 // Well, the create submesh function is expecting const_element_iterators,
458 // and this is just one way to get them...
460 {
463  mesh.create_submesh(new_mesh, it, it_end);
464 }
465 
466 
467 
468 
469 
470 
471 int main (int argc, char ** argv)
472 {
473  LibMeshInit init(argc, argv);
474 
475  PerfMon perfmon(argv[0]);
476 
477  unsigned int n_subdomains = 1;
478  unsigned int n_rsteps = 0;
479  unsigned char dim = static_cast<unsigned char>(-1); // invalid dimension
480  double dist_fact = 0.;
481  bool verbose = false;
482  BoundaryMeshWriteMode write_bndry = BM_DISABLED;
483  bool convert_first_order = false;
484  unsigned int convert_second_order = 0;
485  bool addinfelems = false;
486  bool triangulate = false;
487  bool do_quality = false;
488  ElemQuality quality_type = DIAGONAL;
489 
490 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
491  InfElemBuilder::InfElemOriginValue origin_x(false, 0.);
492  InfElemBuilder::InfElemOriginValue origin_y(false, 0.);
493  InfElemBuilder::InfElemOriginValue origin_z(false, 0.);
494 #endif
495 
496  bool x_sym=false;
497  bool y_sym=false;
498  bool z_sym=false;
499 
500 
501  std::vector<std::string> names;
502  std::vector<std::string> var_names;
503  std::vector<Number> soln;
504 
505  process_cmd_line(argc, argv, names,
506  n_subdomains, n_rsteps, dim,
507  dist_fact, verbose, write_bndry,
508  convert_first_order,
509  convert_second_order,
510 
511  triangulate,
512 
513  do_quality,
514  quality_type,
515 
516  addinfelems,
517 
518 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
519  origin_x, origin_y, origin_z,
520 #endif
521 
522  x_sym, y_sym, z_sym);
523 
524  UniquePtr<Mesh> mesh_ptr;
525  if (dim == static_cast<unsigned char>(-1))
526  {
527  mesh_ptr.reset(new Mesh(init.comm()));
528  }
529  else
530  {
531  mesh_ptr.reset(new Mesh(init.comm(),dim));
532  }
533 
534  Mesh & mesh = *mesh_ptr;
535 
539  if (!names.empty())
540  {
541  mesh.read(names[0]);
542 
543  if (verbose)
544  {
545  mesh.print_info();
547  }
548 
549  }
550 
551  else
552  {
553  libMesh::out << "No input specified." << std::endl;
554  return 1;
555  }
556 
557 
558 
559 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
560 
561  if (addinfelems)
562  {
563  if (names.size() == 3)
564  libmesh_error_msg("ERROR: Invalid combination: Building infinite elements\n"
565  << "not compatible with solution import.");
566 
567  if (write_bndry != BM_DISABLED)
568  libmesh_error_msg("ERROR: Invalid combination: Building infinite elements\n"
569  << "not compatible with writing boundary conditions.");
570 
571  /*
572  * Sanity checks: -X/Y/Z can only be used, when the
573  * corresponding coordinate is also given (using -x/y/z)
574  */
575  if ((x_sym && !origin_x.first) || // claim x-symmetry, but x-coordinate of origin not given!
576  (y_sym && !origin_y.first) || // the same for y
577  (z_sym && !origin_z.first)) // the same for z
578  libmesh_error_msg("ERROR: When x-symmetry is requested using -X, then\n"
579  << "the option -x <coord> also has to be given.\n"
580  << "This holds obviously for y and z, too.");
581 
582  // build infinite elements
583  InfElemBuilder(mesh).build_inf_elem(origin_x, origin_y, origin_z,
584  x_sym, y_sym, z_sym,
585  verbose);
586 
587  if (verbose)
588  {
589  mesh.print_info();
591  }
592 
593  }
594 
595  // sanity check
596  else if ((origin_x.first || origin_y.first || origin_z.first) ||
597  (x_sym || y_sym || z_sym))
598  libmesh_error_msg("ERROR: -x/-y/-z/-X/-Y/-Z is only to be used when\n"
599  << "the option -a is also specified!");
600 
601 #endif
602 
603 
607  if (names.size() == 3)
608  {
609  // TODO: Read XDR/A mesh file, construct an EquationSystems
610  // object, read XDR/A solution file by calling
611  // es.read(file, READ_HEADER|READ_DATA|READ_ADDITIONAL_DATA);
612  // then store a localized copy of the solution vector into 'soln'.
613  libmesh_error_msg("Importing an XDA solution file with -s is not supported.");
614  }
615 
616 
617 
621  // if (dim == 2 && triangulate)
622  if (triangulate)
623  {
624  if (verbose)
625  libMesh::out << "...Converting to all simplices...\n";
626 
628  }
629 
633  if (do_quality)
634  {
636  sv.reserve(mesh.n_elem());
637 
638  libMesh::out << "Quality type is: " << Quality::name(quality_type) << std::endl;
639 
640  // What are the quality bounds for this element?
641  std::pair<Real, Real> bounds = mesh.elem_ref(0).qual_bounds(quality_type);
642  libMesh::out << "Quality bounds for this element type are: ("
643  << bounds.first
644  << ", "
645  << bounds.second
646  << ") "
647  << std::endl;
648 
649  for (const auto & elem : mesh.active_element_ptr_range())
650  sv.push_back(elem->quality(quality_type));
651 
652  const unsigned int n_bins = 10;
653  libMesh::out << "Avg. shape quality: " << sv.mean() << std::endl;
654 
655  // Find element indices below the specified cutoff.
656  // These might be considered "bad" elements which need refinement.
657  std::vector<dof_id_type> bad_elts = sv.cut_below(0.8);
658  libMesh::out << "Found " << bad_elts.size()
659  << " of " << mesh.n_elem()
660  << " elements below the cutoff." << std::endl;
661 
662  // Compute the histogram for this distribution
663  std::vector<dof_id_type> histogram;
664  sv.histogram(histogram, n_bins);
665 
666  /*
667  for (unsigned int i=0; i<n_bins; i++)
668  histogram[i] = histogram[i] / mesh.n_elem();
669  */
670 
671  const bool do_matlab = true;
672 
673  if (do_matlab)
674  {
675  std::ofstream out ("histo.m");
676 
677  out << "% This is a sample histogram plot for Matlab." << std::endl;
678  out << "bin_members = [" << std::endl;
679  for (unsigned int i=0; i<n_bins; i++)
680  out << static_cast<Real>(histogram[i]) / static_cast<Real>(mesh.n_elem())
681  << std::endl;
682  out << "];" << std::endl;
683 
684  std::vector<Real> bin_coords(n_bins);
685  const Real max = *(std::max_element(sv.begin(), sv.end()));
686  const Real min = *(std::min_element(sv.begin(), sv.end()));
687  const Real delta = (max - min) / static_cast<Real>(n_bins);
688  for (unsigned int i=0; i<n_bins; i++)
689  bin_coords[i] = min + (i * delta) + delta / 2.0 ;
690 
691  out << "bin_coords = [" << std::endl;
692  for (unsigned int i=0; i<n_bins; i++)
693  out << bin_coords[i] << std::endl;
694  out << "];" << std::endl;
695 
696  out << "bar(bin_coords, bin_members, 1);" << std::endl;
697  out << "hold on" << std::endl;
698  out << "plot (bin_coords, 0, 'kx');" << std::endl;
699  out << "xlabel('Quality (0=Worst, 1=Best)');" << std::endl;
700  out << "ylabel('Percentage of elements in each bin');" << std::endl;
701  out << "axis([" << min << "," << max << ",0, max(bin_members)]);" << std::endl;
702 
703  out << "title('" << Quality::name(quality_type) << "');" << std::endl;
704 
705  }
706  }
707 
708 
713  if (convert_first_order)
714  {
715  if (verbose)
716  libMesh::out << "Converting elements to first order counterparts\n";
717 
718  mesh.all_first_order();
719 
720  if (verbose)
721  {
722  mesh.print_info();
724  }
725  }
726 
731  if (convert_second_order > 0)
732  {
733  bool second_order_mode = true;
734  std:: string message = "Converting elements to second order counterparts";
735  if (convert_second_order == 2)
736  {
737  second_order_mode = false;
738  message += ", lower version: Quad4 -> Quad8, not Quad9";
739  }
740 
741  else if (convert_second_order == 22)
742  {
743  second_order_mode = true;
744  message += ", highest version: Quad4 -> Quad9";
745  }
746 
747  else
748  libmesh_error_msg("Invalid value, convert_second_order = " << convert_second_order);
749 
750  if (verbose)
751  libMesh::out << message << std::endl;
752 
753  mesh.all_second_order(second_order_mode);
754 
755  if (verbose)
756  {
757  mesh.print_info();
759  }
760  }
761 
762 
763 #ifdef LIBMESH_ENABLE_AMR
764 
768  if (n_rsteps > 0)
769  {
770  if (verbose)
771  libMesh::out << "Refining the mesh "
772  << n_rsteps << " times"
773  << std::endl;
774 
775  MeshRefinement mesh_refinement (mesh);
776  mesh_refinement.uniformly_refine(n_rsteps);
777 
778  if (verbose)
779  {
780  mesh.print_info();
782  }
783  }
784 
785 
789  if (dist_fact > 0.)
790  {
791  libMesh::out << "Distorting the mesh by a factor of "
792  << dist_fact
793  << std::endl;
794 
795  MeshTools::Modification::distort(mesh,dist_fact);
796  };
797 
798 
799  /*
800  char filechar[81];
801  sprintf(filechar,"%s-%04d.plt", "out", 0);
802  std::string oname(filechar);
803 
804  mesh.write(oname);
805 
806  for (unsigned int step=0; step<100; step++)
807  {
808  // const Real x = .5 + .25*cos((((Real) step)/100.)*3.1415927);
809  // const Real y = .5 + .25*sin((((Real) step)/100.)*3.1415927);
810  const Real x = 2.5*cos((((Real) step)/100.)*3.1415927);
811  const Real y = 2.5*sin((((Real) step)/100.)*3.1415927);
812 
813  const Point p(x,y);
814 
815  for (unsigned int e=0; e<mesh.n_elem(); e++)
816  if (mesh.elem_ref(e).active())
817  mesh.elem_ref(e).set_refinement_flag() = -1;
818 
819 
820 
821  for (unsigned int e=0; e<mesh.n_elem(); e++)
822  if (mesh.elem_ref(e).active())
823  {
824  const Point diff = mesh.elem_ref(e).centroid(mesh) - p;
825 
826  if (diff.size() < .5)
827  {
828  if (mesh.elem_ref(e).level() < 4)
829  mesh.elem_ref(e).set_refinement_flag() = 1;
830  else if (mesh.elem_ref(e).level() == 4)
831  mesh.elem_ref(e).set_refinement_flag() = 0;
832  }
833  }
834 
835 
836  mesh.mesh_refinement.refine_and_coarsen_elements();
837 
838  char filechar[81];
839  sprintf(filechar,"%s-%04d.plt", "out", step+1);
840  std::string oname(filechar);
841 
842  mesh.write(oname);
843  }
844  */
845 
846 #endif
847 
848 
849 
850  // /**
851  // * Possibly partition the mesh
852  // */
853  if (n_subdomains > 1)
854  mesh.partition(n_subdomains);
855 
856 
860  {
861  if (names.size() >= 2)
862  {
863  /*
864  * When the mesh got refined, it is likely that
865  * the user does _not_ want to write also
866  * the coarse elements, but only the active ones.
867  * Use Mesh::create_submesh() to create a mesh
868  * of only active elements, and then write _this_
869  * new mesh.
870  */
871  if (n_rsteps > 0)
872  {
873  if (verbose)
874  libMesh::out << " Mesh got refined, will write only _active_ elements." << std::endl;
875 
876  Mesh new_mesh (init.comm(), mesh.mesh_dimension());
877 
878  construct_mesh_of_active_elements(new_mesh, mesh);
879 
880  // now write the new_mesh
881  if (names.size() == 2)
882  new_mesh.write(names[1]);
883  else if (names.size() == 3)
884  new_mesh.write(names[1], soln, var_names);
885  else
886  libmesh_error_msg("Invalid names.size() = " << names.size());
887  }
888  else
889  {
890  if (names.size() == 2)
891  mesh.write(names[1]);
892  else if (names.size() == 3)
893  mesh.write(names[1], soln, var_names);
894  else
895  libmesh_error_msg("Invalid names.size() = " << names.size());
896  }
897 
898 
899 
903  if (write_bndry != BM_DISABLED)
904  {
905  BoundaryMesh boundary_mesh
906  (mesh.comm(), cast_int<unsigned char>(mesh.mesh_dimension()-1));
907 
908  std::string boundary_name = "bndry_";
909  boundary_name += names[1];
910 
911  if (write_bndry == BM_MESH_ONLY)
912  mesh.get_boundary_info().sync(boundary_mesh);
913 
914  else
915  libmesh_error_msg("Invalid value write_bndry = " << write_bndry);
916 
917  if (names.size() == 2)
918  boundary_mesh.write(boundary_name);
919  else if (names.size() == 3)
920  boundary_mesh.write(boundary_name, soln, var_names);
921  }
922  }
923  };
924 
925  return 0;
926 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition: mesh_base.h:117
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.
void sync(UnstructuredMesh &boundary_mesh)
Generates boundary_mesh data structures corresponding to the mesh data structures.
virtual dof_id_type n_elem() const libmesh_override
virtual void all_first_order() libmesh_override
Converts a mesh with higher-order elements into a mesh with linear elements.
void print_summary(std::ostream &out=libMesh::out) const
Prints a summary of the boundary information.
The definition of the const_element_iterator struct.
Definition: mesh_base.h:1494
unsigned int dim
void distort(MeshBase &mesh, const Real factor, const bool perturb_boundary=false)
Randomly perturb the nodal locations.
virtual void all_second_order(const bool full_ordered=true) libmesh_override
Converts a (conforming, non-refined) mesh with linear elements into a mesh with second-order elements...
MeshBase & mesh
virtual Real mean() const
Definition: statistics.C:74
The StatisticsVector class is derived from the std::vector<> and therefore has all of its useful feat...
Definition: statistics.h:67
int main(int argc, char **argv)
Definition: meshtool.C:471
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition: libmesh.h:62
The libMesh namespace provides an interface to certain functionality in the library.
long double max(long double a, double b)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual void partition(const unsigned int n_parts)
Call the default partitioner (currently metis_partition()).
Definition: mesh_base.C:462
void create_submesh(UnstructuredMesh &new_mesh, const_element_iterator &it, const const_element_iterator &it_end) const
Constructs a mesh called "new_mesh" from the current mesh by iterating over the elements between it a...
void all_tri(MeshBase &mesh)
Converts the 2D quadrilateral elements of a Mesh into triangular elements.
void construct_mesh_of_active_elements(Mesh &new_mesh, const Mesh &mesh)
Definition: meshtool.C:459
The BoundaryMesh is a Mesh in its own right, but it contains a description of the boundary of some ot...
Definition: boundary_mesh.h:39
std::pair< bool, double > InfElemOriginValue
Useful typedef.
This is the MeshRefinement class.
void usage(const std::string &progName)
Definition: meshtool.C:65
virtual std::vector< dof_id_type > cut_below(Real cut) const
Definition: statistics.C:325
virtual void write(const std::string &name) libmesh_override
Write the file specified by name.
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
virtual void histogram(std::vector< dof_id_type > &bin_members, unsigned int n_bins=10)
Definition: statistics.C:178
BoundaryMeshWriteMode
Definition: meshtool.C:62
const Parallel::Communicator & comm() const
Definition: libmesh.h:81
PAPI stands for Performance Application Programming Interface.
Definition: perfmon.h:52
This class is used to build infinite elements on top of an existing mesh.
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:490
virtual element_iterator active_elements_begin() libmesh_override
Active, local, and negation forms of the element iterators described above.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Point build_inf_elem(const bool be_verbose=false)
Build infinite elements atop a volume-based mesh, determine origin automatically. ...
const Parallel::Communicator & comm() const
OStreamProxy out
void process_cmd_line(int argc, char **argv, std::vector< std::string > &names, unsigned int &n_subdomains, unsigned int &n_rsteps, unsigned char &dim, double &dist_fact, bool &verbose, BoundaryMeshWriteMode &write_bndry, bool &convert_first_order, unsigned int &convert_second_order, bool &triangulate, bool &do_quality, ElemQuality &quality_type, bool &addinfelems,#ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS InfElemBuilder::InfElemOriginValue &origin_x, InfElemBuilder::InfElemOriginValue &origin_y, InfElemBuilder::InfElemOriginValue &origin_z,#endif bool &x_sym, bool &y_sym, bool &z_sym)
Definition: meshtool.C:176
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
virtual SimpleRange< element_iterator > active_element_ptr_range() libmesh_override
ElemQuality
Defines an enum for element quality metrics.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition: mesh.h:50
long double min(long double a, double b)
void print_info(std::ostream &os=libMesh::out) const
Prints relevant information about the mesh.
Definition: mesh_base.C:448
virtual std::pair< Real, Real > qual_bounds(const ElemQuality) const
Definition: elem.h:875
virtual element_iterator active_elements_end() libmesh_override
void uniformly_refine(unsigned int n=1)
Uniformly refines the mesh n times.