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 {
484 } render_type;
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 }
rheolef::geo_basic
generic mesh with rerefence counting
Definition: domain_indirect.h:64
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 format format paraview
Definition: iorheo-members.h:125
usage
void usage()
Definition: geo.cc:438
rheolef::delete_suffix
string delete_suffix(const string &name, const string &suffix)
delete_suffix: see the rheostream page for the full documentation
Definition: rheostream.cc:222
mkgeo_contraction.ny
ny
Definition: mkgeo_contraction.sh:150
rheolef::point_basic< size_t >
check_macro
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
rheolef::append_dir_to_rheo_path
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
no_render
Definition: geo.cc:475
rheolef::is_float
bool is_float(const string &s)
is_float: see the rheostream page for the full documentation
Definition: rheostream.cc:476
plotmtv_render
Definition: geo.cc:478
atom_render
Definition: geo.cc:481
rheolef::get_basename
string get_basename(const string &name)
get_basename: see the rheostream page for the full documentation
Definition: rheostream.cc:254
rheolef::space_basic
the finite element space
Definition: space.h:352
vtk_render
Definition: geo.cc:479
rheolef::geo_element
see the geo_element page for the full documentation
Definition: geo_element.h:102
rheolef::pow
space_mult_list< T, M > pow(const space_basic< T, M > &X, size_t n)
Definition: space_mult.h:120
rheolef::odiststream::os
std::ostream & os()
Definition: diststream.h:236
paraview_render
Definition: geo.cc:477
rheolef.h
rheolef - reference manual
bamg
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format bamg
Definition: iorheo-members.h:71
vtk
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
Definition: iorheo-members.h:113
gmsh
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
Definition: iorheo-members.h:83
render_type
render_type
Definition: branch.cc:491
mkgeo_ball.clean
clean
Definition: mkgeo_ball.sh:335
rheolef::din
idiststream din
see the diststream page for the full documentation
Definition: diststream.h:427
extrema_element_measure
Float extrema_element_measure(const geo_basic< Float, sequential > &omega, bool is_min)
Definition: geo.cc:509
rheolef::field_basic
Definition: field_expr_utilities.h:38
render_type
render_type
Definition: geo.cc:474
rheolef::environment
see the environment page for the full documentation
Definition: environment.h:115
set_input_format
void set_input_format(idiststream &in, std::string input_format)
Definition: geo.cc:464
mkgeo_sector.m
m
Definition: mkgeo_sector.sh:118
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 format format format gnuplot
Definition: iorheo-members.h:121
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::details::max_abs
T max_abs(const T &x)
Definition: keller_details.h:47
main
int main(int argc, char **argv)
Definition: geo.cc:516
Float
see the Float page for the full documentation
rheolef::point
point_basic< Float > point
Definition: point.h:164
point
see the point page for the full documentation
mkgeo_grid.sys_coord
sys_coord
Definition: mkgeo_grid.sh:171
rheolef::form_basic
Definition: form.h:144
mkgeo_ball.verbose
verbose
Definition: mkgeo_ball.sh:133
rheolef::to_float
Float to_float(const string &s)
to_float: see the rheostream page for the full documentation
Definition: rheostream.cc:494
file_render
Definition: geo.cc:483
rheo
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color rheo
Definition: iorheo-members.h:65
rheolef::dout
odiststream dout(cout)
see the diststream page for the full documentation
Definition: diststream.h:430
rheolef::distributed
distributed
Definition: asr.cc:228
mkgeo_contraction.hmin
hmin
Definition: mkgeo_contraction.sh:129
mkgeo_contraction.nx
nx
Definition: mkgeo_contraction.sh:124
mkgeo_contraction.name
name
Definition: mkgeo_contraction.sh:133
x3d_render
Definition: geo.cc:482
gnuplot_render
Definition: geo.cc:476
gmsh_render
Definition: geo.cc:480
rheolef::dist
T dist(const point_basic< T > &x, const point_basic< T > &y)
Definition: point.h:299