File: mshrable.cpp

package info (click to toggle)
mshr 2018.1.0%2Bdfsg1-7
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 1,136 kB
  • sloc: cpp: 9,808; python: 680; makefile: 242; sh: 62; ansic: 11
file content (172 lines) | stat: -rw-r--r-- 5,845 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
// Copyright (C) 2014 Benjamin Kehlet
//
// This file is part of mshr.
//
// mshr is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
// 
// mshr is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
// 
// You should have received a copy of the GNU General Public License
// along with mshr.  If not, see <http://www.gnu.org/licenses/>.
//

#include <mshr.h>
#include <dolfin.h>

#include <boost/program_options/options_description.hpp>
#include <boost/program_options/variables_map.hpp>
#include <boost/program_options/parsers.hpp>
#include <boost/filesystem.hpp>

#include <string>
#include <iostream>

namespace po = boost::program_options;

// This program reads in a surface from file and generates mesh


namespace
{
void print_mesh_statistics(const dolfin::Mesh& m)
{
  std::cout << "Dolfin mesh of topological dimension " << m.topology().dim() << std::endl;
  std::cout << "  " << m.num_vertices() << " vertices" << std::endl;
  std::cout << "  " << m.num_cells() << " cells" << std::endl;

  const std::pair<double, double> volume_min_max = mshr::DolfinMeshUtils::cell_volume_min_max(m);
  std::cout << "Min cell volume: " << volume_min_max.first << std::endl;
  std::cout << "Max cell volume: " << volume_min_max.second << std::endl;
  const std::pair<double, double> radii_ratio = dolfin::MeshQuality::radius_ratio_min_max(m);
  std::cout << "Minimum cell radii ratio: " << radii_ratio.first  << std::endl;
  std::cout << "Maximum cell radii ratio: " << radii_ratio.second << std::endl;
}
//-----------------------------------------------------------------------------
// Define options and parse command line
void handle_commandline(int argc, char** argv, po::variables_map &vm)
{
  // Command line options
  po::options_description visible("Generate mesh from surface file");
  visible.add_options()
    ("outfile,o",    po::value<std::string>(), "Filename of generated Dolfin mesh")
    ("resolution,r", po::value<double>()->default_value(15.0), "Resolution of result mesh")
    ("stats,s", "Write some statistics of the mesh to stdout")
    ("polyout", po::value<std::string>(), "Write the polyhedron to a .off file (and do not create a mesh)")
    ("polystats", "Write statistics of polyhedron (and do not create a mesh")
    ("backend,b", po::value<std::string>()->default_value("cgal"), "Use 3D mesh generation backend [tetgen|cgal]")
    ("degenerate_tolerance", po::value<double>()->default_value(1e-12), "Tolerance for considering a facet as degenerate. Set to 0 to not remove degenerate facets")
    ("check-mesh", "Check consistency of output mesh (most for debugging/testing")
    ("verbose,v", "Output more information about what is going on")
    ("help,h",   "write help message");

  // Options not shown to the user
  po::options_description hidden("Hidden options");
  hidden.add_options()
    ("input-file", po::value<std::string>(), "Input file")
    ;

  po::options_description cmdline("Generate mesh from surface file");
  cmdline.add(visible).add(hidden);

  //Command line positional arguments
  po::positional_options_description p;
  p.add("input-file", 1);

  // parse command line
  po::store(po::command_line_parser(argc, argv)
            .options(cmdline)
            .positional(p).run(), vm);

  // Print help message if requested 
  // (before notify to avoid error messages if only --help is given)
  if (vm.count("help"))
  {
    std::cout << visible << std::endl;
    exit(EXIT_SUCCESS);
  }

  po::notify(vm);

}
} //end anonymous namespace
//-----------------------------------------------------------------------------
int main(int argc, char** argv)
{
  // Ensure Dolfin initializes MPI correctly.
  #ifdef HAS_MPI
  dolfin::SubSystemsManager::init_mpi();
  #endif

  po::variables_map vm;
  handle_commandline(argc, argv, vm);

  if (vm.count("verbose"))
    dolfin::set_log_level(dolfin::TRACE);

  // Read the infile
  if (!boost::filesystem::exists(vm["input-file"].as<std::string>()))
  {
    std::cerr << "File " << vm["input-file"].as<std::string>() << "does not exist" << std::endl;
    exit(1);
  }


  std::shared_ptr<mshr::Surface3D> surf(new mshr::Surface3D(vm["input-file"].as<std::string>()));
  surf->degenerate_tolerance = vm["degenerate_tolerance"].as<double>();

  // Operations that disable mesh generation
  if (vm.count("polyout") || vm.count("polystats"))
  {
    mshr::CSGCGALDomain3D domain(surf);

    if (vm.count("polyout"))
    {
      std::string extension = boost::filesystem::extension(vm["polyout"].as<std::string>());

      if (extension != ".off")
      {
        std::cerr << "Unknown file type: " << extension << std::endl;
        exit(1);
      }

      domain.save_off(vm["polyout"].as<std::string>());
    }

    if (vm.count("polystats"))
      std::cout << domain.str(true) << std::endl;

    exit(EXIT_SUCCESS);
  }

  // Generate the mesh
  std::shared_ptr<dolfin::Mesh> m = mshr::generate_mesh(surf,
                                                        vm["resolution"].as<double>(),
                                                        vm["backend"].as<std::string>());

  // Output mesh if requested
  if (vm.count("outfile"))
  {
    dolfin::File f(vm["outfile"].as<std::string>());
    f << *m;
  }
    
  if (vm.count("stats"))
    print_mesh_statistics(*m);

  if (vm.count("check-mesh"))
  {
    if (!mshr::DolfinMeshUtils::check_mesh(*m))
    {
      std::cout << "  Error: Mesh check failed" << std::endl;
      return EXIT_FAILURE;
    }
  }

  return EXIT_SUCCESS;
}