Rheolef  7.1
an efficient C++ finite element environment
geo.cc
Go to the documentation of this file.
1 //
2 // This file is part of Rheolef.
3 //
4 // Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
5 //
6 // Rheolef is free software; you can redistribute it and/or modify
7 // it under the terms of the GNU General Public License as published by
8 // the Free Software Foundation; either version 2 of the License, or
9 // (at your option) any later version.
10 //
11 // Rheolef is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 // GNU General Public License for more details.
15 //
16 // You should have received a copy of the GNU General Public License
17 // along with Rheolef; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 //
20 // =========================================================================
21 // author: Pierre.Saramito@imag.fr
22 // date: 16 september 2011
23 
24 namespace rheolef { } // namespace rheolef
434 
435 # include "rheolef.h"
436 # include "rheolef/iofem.h"
437 using namespace rheolef;
438 void usage() {
439  std::cerr
440  << "geo: usage:" << std::endl
441  << "geo "
442  << "{-|file[.geo[.gz]]}"
443  << "[-Idir|-I dir] "
444  << "[-check] "
445  << "[-upgrade] "
446  << "[-subdivide int] "
447  << "[-if {geo,bamg,vtk}] "
448  << "[-geo|-gmsh|-gnuplot|-paraview] "
449  << "[-name string] "
450  << "[-[no]full|-[no]fill|-[no]lattice] "
451  << "[-[no]stereo|-[no]full|-[no]cut|-[no]shrink|-[no]showlabel] "
452  << "[-image-format string] "
453  << "[-resolution int [int]] "
454  << "[-add-boundary] "
455  << "[-zr|-rz] "
456  << "[-size] [-n-vertex] "
457  << "[-hmin] [-hmax] "
458  << "[-[min|max]-element-measure] "
459  << "[-round [float]] "
460  << "[-[no]clean] [-[no]execute] [-[no]verbose] "
461  << std::endl;
462  exit (1);
463 }
464 void set_input_format (idiststream& in, std::string input_format)
465 {
466  if (input_format == "bamg") in.is() >> bamg;
467  else if (input_format == "vtk") in.is() >> vtk;
468  else if (input_format == "geo") in.is() >> rheo;
469  else {
470  std::cerr << "geo: invalid input format \""<<input_format<<"\"" << std::endl;
471  usage();
472  }
473 }
474 typedef enum {
485 
486 #ifdef TO_CLEAN
487 void h_extrema (const geo_basic<Float,sequential>& omega, Float& hmin, Float& hmax)
488 {
489  hmin = std::numeric_limits<Float>::max();
490  hmax = 0;
491  for (size_t iedg = 0, nedg = omega.geo_element_ownership(1).size(); iedg < nedg; ++iedg) {
492  const geo_element& E = omega.get_geo_element (1, iedg);
493  const point& x0 = omega.dis_node(E[0]);
494  const point& x1 = omega.dis_node(E[1]);
495  Float hloc = dist(x0,x1);
496  hmin = min(hmin, hloc);
497  hmax = max(hmax, hloc);
498  }
499 #ifdef TODO
500  // M=sequential here (just if this fct goes into the geo<M,T> class later)
501 #ifdef _RHEOLEF_HAVE_MPI
502  hmin = mpi::all_reduce (comm(), hmin, mpi::minimum<Float>());
503  hmax = mpi::all_reduce (comm(), hmax, mpi::maximum<Float>());
504 #endif // _RHEOLEF_HAVE_MPI
505 #endif // TODO
506 }
507 #endif // TO_CLEAN
508 
510 {
511  space_basic<Float,sequential> Xh (omega, "P0");
512  form_basic<Float,sequential> m (Xh, Xh, "mass");
513  field_basic<Float,sequential> one (Xh, 1);
514  return is_min ? (m*one).min_abs() : (m*one).max_abs();
515 }
516 int main(int argc, char**argv) {
517  environment distributed(argc, argv);
518  check_macro (communicator().size() == 1, "geo: command may be used as mono-process only");
519  // ----------------------------
520  // default options
521  // ----------------------------
522  bool do_upgrade = false ;
523  bool do_check = false ;
524  bool do_add_bdry = false ;
525  bool do_stereo = false ;
526  bool show_n_element = false ;
527  bool show_n_vertex = false ;
528  bool show_sys_coord = false ;
529  bool show_hmin = false ;
530  bool show_hmax = false ;
531  bool show_xmin = false ;
532  bool show_xmax = false ;
533  bool show_min_element_measure = false ;
534  bool show_max_element_measure = false ;
535  bool do_round = false ;
536  Float round_prec = pow(10., -std::numeric_limits<Float>::digits10/2);
537  std::string name = "output";
538  dout.os() << verbose; bool bverbose = true;
539  render_type render = no_render;
541  dout.os() << lattice;
542  dout.os() << nofill;
543  dout.os() << grid;
544  dout.os() << showlabel;
545  // this normal is not so bad for the dirichlet.cc demo on the cube:
546  dout.os() << setnormal(point(-0.015940197423022637, -0.9771157601293953, -0.21211011624358989));
547  dout.os() << setorigin(point(std::numeric_limits<Float>::max()));
548  std::string filename;
549  std::string input_format = "geo";
550  std::string sys_coord = "";
551  // ----------------------------
552  // scan the command line
553  // ----------------------------
554  for (int i = 1; i < argc; i++) {
555 
556  // general options:
557  if (strcmp (argv[i], "-clean") == 0) dout.os() << clean;
558  else if (strcmp (argv[i], "-noclean") == 0) dout.os() << noclean;
559  else if (strcmp (argv[i], "-execute") == 0) dout.os() << execute;
560  else if (strcmp (argv[i], "-noexecute") == 0) dout.os() << noexecute;
561  else if (strcmp (argv[i], "-verbose") == 0) { bverbose = true; dout.os() << verbose; }
562  else if (strcmp (argv[i], "-noverbose") == 0) { bverbose = false; dout.os() << noverbose; }
563  else if (strcmp (argv[i], "-add-boundary") == 0) { do_add_bdry = true; }
564  else if (strcmp (argv[i], "-rz") == 0) { sys_coord = "rz"; }
565  else if (strcmp (argv[i], "-zr") == 0) { sys_coord = "zr"; }
566  else if (strcmp (argv[i], "-size") == 0) { show_n_element = true; }
567  else if (strcmp (argv[i], "-n-element") == 0) { show_n_element = true; }
568  else if (strcmp (argv[i], "-n-vertex") == 0) { show_n_vertex = true; }
569  else if (strcmp (argv[i], "-sys-coord") == 0) { show_sys_coord = true; }
570  else if (strcmp (argv[i], "-hmin") == 0) { show_hmin = true; }
571  else if (strcmp (argv[i], "-hmax") == 0) { show_hmax = true; }
572  else if (strcmp (argv[i], "-xmin") == 0) { show_xmin = true; }
573  else if (strcmp (argv[i], "-xmax") == 0) { show_xmax = true; }
574  else if (strcmp (argv[i], "-min-element-measure") == 0) { show_min_element_measure = true; }
575  else if (strcmp (argv[i], "-max-element-measure") == 0) { show_max_element_measure = true; }
576  else if (strcmp (argv[i], "-I") == 0) {
577  if (i == argc-1) { std::cerr << "geo -I: option argument missing" << std::endl; usage(); }
578  append_dir_to_rheo_path (argv[++i]);
579  }
580  else if (argv [i][0] == '-' && argv [i][1] == 'I') { append_dir_to_rheo_path (argv[i]+2); }
581  // output file option:
582  else if (strcmp (argv[i], "-geo") == 0) { dout.os() << rheo; render = file_render; }
583  else if (strcmp (argv[i], "-gmsh") == 0) { dout.os() << gmsh; render = gmsh_render; }
584 
585  // render spec:
586  else if (strcmp (argv[i], "-gnuplot") == 0) { dout.os() << gnuplot; render = gnuplot_render; }
587  else if (strcmp (argv[i], "-paraview") == 0) { dout.os() << paraview; render = paraview_render; }
588 
589  // render option:
590  else if (strcmp (argv[i], "-full") == 0) { dout.os() << full; }
591  else if (strcmp (argv[i], "-nofull") == 0) { dout.os() << nofull; }
592  else if (strcmp (argv[i], "-fill") == 0) { dout.os() << fill; }
593  else if (strcmp (argv[i], "-nofill") == 0) { dout.os() << nofill; }
594  else if (strcmp (argv[i], "-stereo") == 0) { dout.os() << stereo; do_stereo = true; }
595  else if (strcmp (argv[i], "-nostereo") == 0) { dout.os() << nostereo; }
596  else if (strcmp (argv[i], "-cut") == 0) { dout.os() << cut; }
597  else if (strcmp (argv[i], "-nocut") == 0) { dout.os() << nocut; }
598  else if (strcmp (argv[i], "-shrink") == 0) { dout.os() << shrink; }
599  else if (strcmp (argv[i], "-noshrink") == 0) { dout.os() << noshrink; }
600  else if (strcmp (argv[i], "-lattice") == 0) { dout.os() << lattice; }
601  else if (strcmp (argv[i], "-nolattice") == 0) { dout.os() << nolattice; }
602  else if (strcmp (argv[i], "-showlabel") == 0) { dout.os() << showlabel; }
603  else if (strcmp (argv[i], "-noshowlabel") == 0){ dout.os() << noshowlabel; }
604  else if (strcmp (argv[i], "-subdivide") == 0) {
605  if (i == argc-1) { std::cerr << "geo -subdivide: option argument missing" << std::endl; usage(); }
606  size_t nsub = atoi(argv[++i]);
607  dout.os() << setsubdivide (nsub);
608  }
609  else if (strcmp (argv[i], "-image-format") == 0) {
610  if (i == argc-1) { std::cerr << "geo -image-format: option argument missing" << std::endl; usage(); }
611  std::string format = argv[++i];
612  if (format == "jpeg") format = "jpg";
613  if (format == "postscript") format = "ps";
614  dout.os() << setimage_format(format);
615  }
616  else if (strcmp (argv[i], "-resolution") == 0) {
617  if (i == argc-1 || !isdigit(argv[i+1][0])) { std::cerr << "geo -resolution: option argument missing" << std::endl; usage(); }
618  size_t nx = atoi(argv[++i]);
619  size_t ny = (i < argc-1 && isdigit(argv[i+1][0])) ? atoi(argv[++i]) : nx;
620  dout.os() << setresolution(point_basic<size_t>(nx,ny));
621  }
622  else if (strcmp (argv[i], "-round") == 0) {
623 
624  do_round = true;
625  if (i+1 < argc && is_float(argv[i+1])) {
626  round_prec = to_float (argv[++i]);
627  }
628  }
629  // input options:
630  else if (strcmp (argv[i], "-upgrade") == 0) { do_upgrade = true; dout.os() << rheo; render = file_render; }
631  else if (strcmp (argv[i], "-check") == 0) { do_check = true; }
632  else if (strcmp (argv[i], "-name") == 0) {
633  if (i == argc-1) { std::cerr << "geo -name: option argument missing" << std::endl; usage(); }
634  name = argv[++i];
635  }
636  else if (strcmp (argv[i], "-if") == 0 ||
637  strcmp (argv[i], "-input-format") == 0) {
638  if (i == argc-1) { std::cerr << "geo "<<argv[i]<<": option argument missing" << std::endl; usage(); }
639  input_format = argv[++i];
640  }
641  else if (strcmp (argv[i], "-") == 0) {
642  filename = "-"; // geo on stdin
643  }
644  else if (argv[i][0] != '-') {
645  filename = argv[i]; // input on file
646  }
647  else { usage(); }
648  }
649  // ----------------------------
650  // geo treatment
651  // ----------------------------
652  if (filename == "") {
653  std::cerr << "geo: no input file specified" << std::endl;
654  usage();
655  } else if (filename == "-") {
656  // geo on stdin
657  if (do_upgrade) std::cin >> upgrade;
658  std::string thename;
659  if (name != "output") thename = name;
660  std::cin >> setbasename(thename);
661  set_input_format (din, input_format);
662  din >> omega;
663  dout.os() << setbasename(name) << reader_on_stdin;
664  } else {
665  // input geo on file
666  std::string suffix = input_format;
667  if (name == "output") {
668  name = get_basename(delete_suffix (delete_suffix (filename, "gz"), suffix));
669  }
670  idiststream ips;
671  ips.open (filename, suffix);
672  check_macro(ips.good(), "\"" << filename << "[."<<suffix<<"[.gz]]\" not found.");
673  if (do_upgrade) ips.is() >> upgrade;
674  std::string root_name = delete_suffix (delete_suffix(filename, "gz"), suffix);
675  name = get_basename (root_name);
676  ips.is() >> setbasename(name);
677  set_input_format (ips, input_format);
678  ips >> omega;
679  dout.os() << setbasename(name);
680  }
681  if (sys_coord != "") { omega.set_coordinate_system(sys_coord); }
682  if (do_add_bdry) { omega.boundary(); }
683 
684  if (render == file_render) {
685  size_t nsub = iorheo::getsubdivide (std::cout);
686  if (nsub >= 2) {
687  geo_basic<Float,sequential> new_omega;
688  new_omega.build_by_subdividing (omega, nsub);
689  omega = new_omega;
690  }
691  }
692  if (do_check) {
693  check_macro (omega.check(bverbose), "geo check failed");
694  }
695  bool show = show_min_element_measure
696  || show_max_element_measure
697  || show_hmin || show_hmax
698  || show_xmin || show_xmax
699  || show_n_element || show_n_vertex || show_sys_coord;
700  if (show) {
701  if (show_hmin) { dout << "hmin = " << omega.hmin() << std::endl; }
702  if (show_hmax) { dout << "hmax = " << omega.hmax() << std::endl; }
703  if (show_xmin) { dout << "xmin = " << omega.xmin() << std::endl; }
704  if (show_xmax) { dout << "xmax = " << omega.xmax() << std::endl; }
705  if (show_n_element) { dout << "size = " << omega.dis_size() << std::endl; }
706  if (show_n_vertex) { dout << "n_vertex = " << omega.dis_size(0) << std::endl; }
707  if (show_sys_coord) { dout << "sys_coord = " << omega.coordinate_system_name() << std::endl; }
708  if (show_min_element_measure) {
709  dout << "min_element_measure = " << extrema_element_measure (omega, true) << std::endl;
710  }
711  if (show_max_element_measure) {
712  dout << "max_element_measure = " << extrema_element_measure (omega, false) << std::endl;
713  }
714  return 0;
715  }
716  if (do_round) {
717  dout.os() << setrounding_precision(round_prec);
718  }
719  if (render == no_render) {
720  // try to choose the best render from dimension
721 #if (_RHEOLEF_PARAVIEW_VERSION_MAJOR >= 5) && (_RHEOLEF_PARAVIEW_VERSION_MINOR >= 5)
722  // paraview version >= 5.5 has high order elements
723  if (do_stereo || omega.dimension() >= 2) {
724 #else
725  if (do_stereo || omega.dimension() == 3) {
726 #endif
727  dout.os() << paraview;
728  } else {
729  dout.os() << gnuplot;
730  }
731  }
732  dout << omega;
733 }
render_type
Definition: branch.cc:491
void usage()
Definition: geo.cc:438
int main(int argc, char **argv)
Definition: geo.cc:516
void set_input_format(idiststream &in, std::string input_format)
Definition: geo.cc:464
render_type
Definition: geo.cc:474
@ file_render
Definition: geo.cc:483
@ gnuplot_render
Definition: geo.cc:476
@ no_render
Definition: geo.cc:475
@ paraview_render
Definition: geo.cc:477
@ gmsh_render
Definition: geo.cc:480
@ vtk_render
Definition: geo.cc:479
@ x3d_render
Definition: geo.cc:482
@ plotmtv_render
Definition: geo.cc:478
@ atom_render
Definition: geo.cc:481
Float extrema_element_measure(const geo_basic< Float, sequential > &omega, bool is_min)
Definition: geo.cc:509
see the Float page for the full documentation
see the point page for the full documentation
see the environment page for the full documentation
Definition: environment.h:104
generic mesh with rerefence counting
Definition: geo.h:1089
see the geo_element page for the full documentation
Definition: geo_element.h:102
std::ostream & os()
Definition: diststream.h:236
the finite element space
Definition: space.h:352
point_basic< Float > point
Definition: point.h:164
idiststream din
see the diststream page for the full documentation
Definition: diststream.h:427
odiststream dout(cout)
see the diststream page for the full documentation
Definition: diststream.h:430
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format bamg
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format format format format format format format format format format format format format format format format format format format format paraview
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format format format format format format format format format format format format format format format format format format gnuplot
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format format format format format format format format format format format format format format format vtk
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color rheo
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format format format format format gmsh
string sys_coord
Definition: mkgeo_grid.sh:171
T max_abs(const T &x)
This file is part of Rheolef.
T dist(const point_basic< T > &x, const point_basic< T > &y)
Definition: point.h:299
string delete_suffix(const string &name, const string &suffix)
delete_suffix: see the rheostream page for the full documentation
Definition: rheostream.cc:222
string get_basename(const string &name)
get_basename: see the rheostream page for the full documentation
Definition: rheostream.cc:254
bool is_float(const string &s)
is_float: see the rheostream page for the full documentation
Definition: rheostream.cc:476
void append_dir_to_rheo_path(const string &dir)
append_dir_to_rheo_path: see the rheostream page for the full documentation
Definition: rheostream.cc:330
space_mult_list< T, M > pow(const space_basic< T, M > &X, size_t n)
Definition: space_mult.h:120
Float to_float(const string &s)
to_float: see the rheostream page for the full documentation
Definition: rheostream.cc:494
rheolef - reference manual