~ubuntu-branches/ubuntu/raring/rheolef/raring-proposed

« back to all changes in this revision

Viewing changes to nfem/pbin/geo.cc

  • Committer: Package Import Robot
  • Author(s): Pierre Saramito
  • Date: 2012-04-06 09:12:21 UTC
  • mfrom: (1.1.5)
  • Revision ID: package-import@ubuntu.com-20120406091221-m58me99p1nxqui49
Tags: 6.0-1
* New upstream release 6.0 (major changes):
  - massively distributed and parallel support
  - full FEM characteristic method (Lagrange-Gakerkin method) support
  - enhanced users documentation 
  - source code supports g++-4.7 (closes: #667356)
* debian/control: dependencies for MPI distributed solvers added
* debian/rules: build commands simplified
* debian/librheolef-dev.install: man1/* to man9/* added
* debian/changelog: package description rewritted (closes: #661689)

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
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
/*Prog:geo
 
22
NAME: @code{geo} - plot a finite element mesh 
 
23
@pindex geo
 
24
@cindex plotting
 
25
@cindex plotting mesh
 
26
@cindex mesh graphic representation
 
27
@clindex geo
 
28
SYNOPSIS:
 
29
  @example
 
30
        geo @var{options} @var{mesh}[.geo[.gz]]
 
31
  @end example
 
32
DESCRIPTION:       
 
33
  @noindent
 
34
  Plot or upgrade a finite element mesh.
 
35
EXAMPLES:
 
36
  Plot a mesh:
 
37
  @example
 
38
        geo square.geo
 
39
        geo box.geo
 
40
        geo box.geo -full
 
41
  @end example
 
42
  Plot a mesh into a file:
 
43
  @example
 
44
        geo square.geo -image-format png
 
45
  @end example
 
46
  Convert from a old geo file format to the new one:
 
47
  @example
 
48
        geo -upgrade - < square-old.geo > square.geo
 
49
  @end example
 
50
  See below for the geo file format scpecification.
 
51
  The old file format does not contains edges and faces connectivity in 3d geometries,
 
52
  or edges connectivity in 2d geometries.
 
53
  The converter add it automatically into the upgraded file format.
 
54
  Conversely, the old file format is usefull when combined with a
 
55
  translator from another file format that do not provides edges and faces connectivity.
 
56
 
 
57
INPUT FILE SPECIFICATION:
 
58
  @table @code
 
59
  @itemx @var{filename}
 
60
        specifies the name of the file containing
 
61
        the input mesh.
 
62
        The ".geo" suffix extension is assumed.
 
63
  @itemx -
 
64
        read mesh on standard input instead on a file.
 
65
  @itemx -name
 
66
        when mesh comes from standard input, the mesh name
 
67
        is not known and is set to "output" by default.
 
68
        This option allows to change this default.
 
69
        Useful when dealing with output formats (graphic, format conversion)
 
70
        that creates auxilliary files, based on this name.
 
71
@cindex RHEOPATH environment variable
 
72
  @itemx -I@var{dir} 
 
73
  @itemx -I @var{dir} 
 
74
        Add @var{dir} to the rheolef file search path.
 
75
        This mechanism initializes a search path given by the environment variable @samp{RHEOPATH}.
 
76
        If the environment variable @samp{RHEOPATH} is not set, the default value is the current directory.
 
77
  @itemx -check
 
78
        Check that element orientation is positive.
 
79
  @end table
 
80
 
 
81
INPUT FORMAT OPTIONS:
 
82
@toindex @code{bamg}
 
83
@fiindex @file{.bamg} mesh file
 
84
@toindex @code{vtk}
 
85
@fiindex @file{.vtk} mesh file
 
86
  @table @code
 
87
  @itemx -if @var{format}
 
88
  @itemx -input-format @var{format}
 
89
        load mesh in @var{format} file format,
 
90
        instead of plotting it.
 
91
        Supported output formats are: @code{geo}, @code{bamg}, @code{vtk}.
 
92
        When loading from a file, the corresponding suffix extension is assumed.
 
93
  @end table
 
94
 
 
95
RENDER SPECIFICATION:
 
96
@cindex graphic render
 
97
@toindex @code{gnuplot}
 
98
  @table @code
 
99
  @itemx -gnuplot
 
100
        use gnuplot tool.
 
101
        This is the default.
 
102
@toindex @code{mayavi}
 
103
  @itemx -mayavi
 
104
        use @code{mayavi} tool.
 
105
  @end table
 
106
 
 
107
RENDER OPTIONS:
 
108
  @table @code
 
109
  @itemx -lattice
 
110
  @itemx -nolattice
 
111
        When using a high order geometry, the lattice inside any element appears.
 
112
        Default is on;
 
113
  @itemx -subdivide int
 
114
        When using a high order geometry, the number of points per edge used to draw
 
115
        a curved element. Default value is the mesh order.
 
116
  @itemx -full
 
117
  @itemx -nofull
 
118
        All internal edges appears, for 3d meshes.
 
119
        Default is off;
 
120
  @itemx -fill
 
121
  @itemx -nofill
 
122
        Fill mesh faces using light effects, when available.
 
123
@cindex @code{stereo 3D rendering}
 
124
  @itemx -stereo
 
125
  @itemx -[no]stereo
 
126
        Rendering mode suitable for red-blue anaglyph 3D stereoscopic glasses.
 
127
        Option only available with @code{mayavi}.
 
128
  @itemx -cut
 
129
  @itemx -nocut
 
130
        cut by plane and clip (with @code{mayavi} only).
 
131
  @end table
 
132
 
 
133
OUTPUT FILE FORMAT OPTIONS:
 
134
  @table @code
 
135
  @itemx -geo
 
136
        output mesh on standard output stream in geo text file format,
 
137
        instead of plotting it.
 
138
  @itemx -upgrade
 
139
        Convert from a old geo file format to the new one.
 
140
@toindex @code{gmsh}
 
141
@fiindex @file{.gmsh} mesh file
 
142
  @itemx -gmsh
 
143
        output mesh on standard output stream in gmsh text file format,
 
144
        instead of plotting it.
 
145
@cindex graphic render
 
146
@cindex image file format
 
147
@cindex graphic render
 
148
@fiindex @file{.png} image 
 
149
@fiindex @file{.gif} image
 
150
@fiindex @file{.jpg} image
 
151
@fiindex @file{.ps} image
 
152
@fiindex @file{.pdf} image
 
153
  @itemx -image-format @code{string}
 
154
        The argument is any valid image format, such as
 
155
        bitmap @code{png}, @code{jpg}, @code{gif} or vectorial @code{ps} or @code{pdf}
 
156
        image file formats, and that could be handled by the corresponding render.
 
157
        The output file is e.g. @file{@emph{basename}.png} when @emph{basename}
 
158
        is the name of the mesh, or can be set with the @code{-name} option.
 
159
  @end table
 
160
 
 
161
OTHERS OPTIONS:
 
162
  @table @code
 
163
    @itemx -add-boundary
 
164
        check for a domain named "boundary"; If this domain
 
165
        does not exists, extract the boundary of the geometry
 
166
        and append it to the domain list.
 
167
@toindex @code{bamg}
 
168
        This command is usefull for mesh converted from 
 
169
        generators, as @code{bamg}, that cannot have more than
 
170
        one domain specification per boundary edge.
 
171
 
 
172
    @itemx -verbose
 
173
        print messages related to graphic files created and
 
174
        command system calls (this is the default).
 
175
 
 
176
    @itemx -noverbose
 
177
        does not print previous messages.
 
178
 
 
179
    @itemx -clean
 
180
        clear temporary graphic files (this is the default).
 
181
    
 
182
    @itemx -noclean
 
183
        does not clear temporary graphic files.
 
184
 
 
185
    @itemx -execute
 
186
        execute graphic command (this is the default).
 
187
    
 
188
    @itemx -noexecute
 
189
        does not execute graphic command. Generates only 
 
190
        graphic files. This is usefull in conjuction with the
 
191
        "-noclean" command.
 
192
 
 
193
  @itemx -check
 
194
  @itemx -dump
 
195
        used for debug purpose.
 
196
  @end table
 
197
 
 
198
GEO FILE FORMAT:       
 
199
        TODO
 
200
 
 
201
INQUIRE OPTIONS:
 
202
        TODO
 
203
 
 
204
AUTHOR: Pierre.Saramito@imag.fr
 
205
DATE: 16 september 2011
 
206
End:
 
207
*/
 
208
 
 
209
# include "rheolef/geo.h"
 
210
# include "rheolef/iorheo.h"
 
211
# include "rheolef/iofem.h"
 
212
# include "rheolef/rheostream.h"
 
213
using namespace rheolef;
 
214
using namespace std;
 
215
void usage() {
 
216
      cerr << "geo: usage:" << endl
 
217
           << "geo "
 
218
           << "{-|file[.geo[.gz]]}"
 
219
           << "[-Idir|-I dir] "
 
220
           << "[-check] "
 
221
           << "[-upgrade] "
 
222
           << "[-if {geo,bamg,vtk}] "
 
223
           << "[-geo|-gmsh|-gnuplot|-mayavi] "
 
224
           << "[-name string] "
 
225
           << "[-[no]full|-[no]fill|-[no]lattice] "
 
226
           << "[-subdivide int] "
 
227
           << "[-[no]stereo|-[no]full|-[no]cut] "
 
228
           << "[-image-format string] "
 
229
           << "[-add-boundary] "
 
230
           << "[-[no]clean] [-[no]execute] [-[no]verbose] "
 
231
           << endl;
 
232
      exit (1);
 
233
}
 
234
void set_input_format (idiststream& in, std::string input_format) 
 
235
{
 
236
       if (input_format == "bamg") in.is() >> bamg;
 
237
  else if (input_format == "vtk")  in.is() >> vtk;
 
238
  else if (input_format == "geo")  in.is() >> rheo;
 
239
  else {
 
240
      cerr << "geo: invalid input format \""<<input_format<<"\"" << endl;
 
241
      usage();
 
242
  }
 
243
}
 
244
typedef enum {
 
245
        no_render,
 
246
        gnuplot_render,
 
247
        mayavi_render,
 
248
        plotmtv_render,
 
249
        vtk_render,
 
250
        gmsh_render,
 
251
        atom_render,
 
252
        x3d_render,
 
253
        file_render
 
254
} render_type;
 
255
 
 
256
template <class T>
 
257
void add_boundary (const geo_basic<T,sequential>& omega) 
 
258
{
 
259
  using namespace std;
 
260
  typedef typename geo_basic<T,sequential>::size_type size_type;
 
261
  typedef pair<size_type,size_type>                   pair_t;
 
262
  typedef heap_allocator<size_type>                   alloc_t;
 
263
  typedef geo_element_auto<alloc_t>                   element_t;
 
264
  typedef array<element_t, sequential, alloc_t>       e_array_t;
 
265
 
 
266
  if (omega.have_domain_indirect ("boundary")) return;
 
267
  size_type map_dim = omega.map_dimension();
 
268
  if (map_dim == 0) return;
 
269
 
 
270
  vector<index_set> ball[4];
 
271
  // compute bdry_tmp:
 
272
  // 1) vertice-to-element map:
 
273
  ball[map_dim].resize (omega.n_vertex());
 
274
  for (size_type ie = 0, ne = omega.size(map_dim); ie < ne; ie++) {
 
275
    const geo_element& K = omega.get_geo_element (map_dim, ie);
 
276
    for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
 
277
      ball [map_dim][K[iloc]].insert(ie);
 
278
    }
 
279
  }
 
280
  // 2) scan sides, count adjacency and list bdry side indexes
 
281
  size_type sub_map_dim = map_dim-1;
 
282
  std::list<size_type, alloc_t> iside;
 
283
  for (size_type is = 0, ns = omega.size(sub_map_dim); is < ns; is++) {
 
284
    const geo_element& S = omega.get_geo_element (sub_map_dim, is);
 
285
    index_set adj = ball [map_dim][S[0]];
 
286
    for (size_type iloc = 1, nloc = S.size(); iloc < nloc; iloc++) {
 
287
      adj.inplace_intersection (ball [map_dim][S[iloc]]);
 
288
    }
 
289
    if (adj.size() != 1) continue;
 
290
    // here S is a boundary side:
 
291
    iside.push_back (is);
 
292
  }
 
293
 
 
294
  // 3) copy sides to a bdry array
 
295
  alloc_t alloc;
 
296
  element_t dummy_elt (reference_element::p, 0, alloc);
 
297
  e_array_t bdry_tmp (iside.size(), dummy_elt, alloc);
 
298
  size_type ibdr = 0;
 
299
  for (typename std::list<size_type, alloc_t>::const_iterator iter = iside.begin(), last = iside.end(); 
 
300
        iter != last; ++iter, ++ibdr) {
 
301
    size_type is = *iter;
 
302
    bdry_tmp [ibdr] = omega.get_geo_element (sub_map_dim, is);
 
303
  } 
 
304
  // 4) insert the bdry array into the mesh
 
305
  domain_indirect_basic<sequential> bdry (bdry_tmp, omega, ball);
 
306
  bdry.set_name ("boundary");
 
307
  omega.insert_domain_indirect (bdry);
 
308
}
 
309
int main(int argc, char**argv) {
 
310
    environment distributed(argc, argv);
 
311
    check_macro (communicator().size() == 1, "geo: command may be used as mono-process only");
 
312
    // ----------------------------
 
313
    // default options
 
314
    // ----------------------------
 
315
    bool do_upgrade  = false ;
 
316
    bool do_check    = false ;
 
317
    bool do_add_bdry = false ;
 
318
    string name = "output";
 
319
    dout.os() << verbose; bool bverbose = true;
 
320
    render_type render = no_render;
 
321
    geo_basic<Float,sequential> omega;
 
322
    dout.os() << lattice;
 
323
    dout.os() << nofill;
 
324
    dout.os() << grid;
 
325
#ifdef TO_CLEAN
 
326
    dout.os() << setorigin (point(0,0,0));
 
327
    dout.os() << setnormal (point(1,0,0));
 
328
#endif // TO_CLEAN
 
329
    string filename;
 
330
    string input_format = "geo";
 
331
    // ----------------------------
 
332
    // scan the command line
 
333
    // ----------------------------
 
334
    for (int i = 1; i < argc; i++) {
 
335
 
 
336
        // general options:
 
337
             if (strcmp (argv[i], "-clean") == 0)     dout.os() << clean;
 
338
        else if (strcmp (argv[i], "-noclean") == 0)   dout.os() << noclean;
 
339
        else if (strcmp (argv[i], "-execute") == 0)   dout.os() << execute;
 
340
        else if (strcmp (argv[i], "-noexecute") == 0) dout.os() << noexecute;
 
341
        else if (strcmp (argv[i], "-verbose") == 0)   { bverbose = true; dout.os() << verbose; }
 
342
        else if (strcmp (argv[i], "-noverbose") == 0) { bverbose = false; dout.os() << noverbose; }
 
343
        else if (strcmp (argv[i], "-add-boundary") == 0) { do_add_bdry = true; }
 
344
        else if (strcmp (argv[i], "-I") == 0)         {
 
345
            if (i == argc-1) { cerr << "geo -I: option argument missing" << endl; usage(); }
 
346
            append_dir_to_rheo_path (argv[++i]);
 
347
        }
 
348
        else if (argv [i][0] == '-' && argv [i][1] == 'I') { append_dir_to_rheo_path (argv[i]+2); }
 
349
        // output file option:
 
350
        else if (strcmp (argv[i], "-geo") == 0)       { dout.os() << rheo;     render = file_render; }
 
351
        else if (strcmp (argv[i], "-gmsh") == 0)      { dout.os() << gmsh;     render = gmsh_render; }
 
352
 
 
353
        // render spec:
 
354
        else if (strcmp (argv[i], "-gnuplot") == 0)   { dout.os() << gnuplot;  render = gnuplot_render; }
 
355
        else if (strcmp (argv[i], "-mayavi") == 0)   { dout.os() << mayavi;  render = mayavi_render; }
 
356
 
 
357
        // render option:
 
358
        else if (strcmp (argv[i], "-full") == 0)      { dout.os() << full; }
 
359
        else if (strcmp (argv[i], "-nofull") == 0)    { dout.os() << nofull; }
 
360
        else if (strcmp (argv[i], "-fill") == 0)      { dout.os() << fill; }
 
361
        else if (strcmp (argv[i], "-nofill") == 0)    { dout.os() << nofill; }
 
362
        else if (strcmp (argv[i], "-stereo") == 0)    { dout.os() << stereo << mayavi; render = mayavi_render; }
 
363
        else if (strcmp (argv[i], "-nostereo") == 0)  { dout.os() << nostereo; }
 
364
        else if (strcmp (argv[i], "-cut") == 0)       { dout.os() << cut; }
 
365
        else if (strcmp (argv[i], "-nocut") == 0)     { dout.os() << nocut; }
 
366
        else if (strcmp (argv[i], "-lattice") == 0)   { dout.os() << lattice; }
 
367
        else if (strcmp (argv[i], "-nolattice") == 0) { dout.os() << nolattice; }
 
368
        else if (strcmp (argv[i], "-subdivide") == 0) {
 
369
            if (i == argc-1) { cerr << "geo -subdivide: option argument missing" << endl; usage(); }
 
370
            size_t nsub = atoi(argv[++i]);
 
371
            dout.os() << setsubdivide (nsub);
 
372
        }
 
373
        else if (strcmp (argv[i], "-image-format") == 0) {
 
374
            if (i == argc-1) { cerr << "geo -image-format: option argument missing" << endl; usage(); }
 
375
            string format = argv[++i];
 
376
            if (format == "jpeg")       format = "jpg";
 
377
            if (format == "postscript") format = "ps";
 
378
            dout.os() << setimage_format(format);
 
379
        }
 
380
 
 
381
        // input options:
 
382
        else if (strcmp (argv[i], "-upgrade") == 0)   { do_upgrade = true; dout.os() << rheo; render = file_render; }
 
383
        else if (strcmp (argv[i], "-check") == 0)   { do_check = true; }
 
384
        else if (strcmp (argv[i], "-name") == 0) {
 
385
            if (i == argc-1) { cerr << "geo -name: option argument missing" << endl; usage(); }
 
386
            name = argv[++i];
 
387
        }
 
388
        else if (strcmp (argv[i], "-if") == 0 ||
 
389
                 strcmp (argv[i], "-input-format") == 0) {
 
390
            if (i == argc-1) { cerr << "geo "<<argv[i]<<": option argument missing" << endl; usage(); }
 
391
            input_format = argv[++i];
 
392
        }
 
393
        else if (strcmp (argv[i], "-") == 0) {
 
394
            filename = "-"; // geo on stdin
 
395
        }
 
396
        else if (argv[i][0] != '-') {
 
397
            filename = argv[i]; // input on file        
 
398
        }
 
399
        else { usage(); }
 
400
    }
 
401
    // ----------------------------
 
402
    // geo treatment
 
403
    // ----------------------------
 
404
    if (filename == "") {
 
405
      cerr << "geo: no input file specified" << endl;
 
406
      usage();
 
407
    } else if (filename == "-") {
 
408
      // geo on stdin
 
409
      if (do_upgrade) cin >> upgrade;
 
410
      string thename;
 
411
      if (name != "output") thename = name;
 
412
      cin >> setbasename(thename);
 
413
      set_input_format (din, input_format);
 
414
      din >> omega;
 
415
      dout.os() << setbasename(name) << reader_on_stdin;
 
416
    } else {
 
417
      // input geo on file
 
418
      string suffix = input_format;
 
419
      if (name == "output") {
 
420
        name = get_basename(delete_suffix (delete_suffix (filename, "gz"), suffix));
 
421
      }
 
422
      idiststream ips;
 
423
      ips.open (filename, suffix);
 
424
      check_macro(ips.good(), "\"" << filename << "[."<<suffix<<"[.gz]]\" not found.");
 
425
      if (do_upgrade) ips.is() >> upgrade;
 
426
      std::string root_name = delete_suffix (delete_suffix(filename, "gz"), suffix);
 
427
      std::string name = get_basename (root_name);
 
428
      ips.is() >> setbasename(name);
 
429
      set_input_format (ips, input_format);
 
430
      ips >> omega;
 
431
      dout.os() << setbasename(name);
 
432
    }
 
433
    if (do_add_bdry) {
 
434
      add_boundary (omega);
 
435
    }
 
436
    if (do_check) {
 
437
        check_macro (omega.check(bverbose), "geo check failed");
 
438
    }
 
439
    if (render == no_render) { dout.os() << gnuplot; }
 
440
    dout << omega;
 
441
}