137 #include "rheolef/index_set_header.icc"
138 #include "rheolef/index_set_body.icc"
153 cerr <<
"bamg2geo: usage:" << endl
155 <<
"[-cartesian|-rz|-zr] "
156 <<
"[input.bamg [input.dmn]]"
161 struct node_type: std::array<double,2> {
163 struct element_type: std::array<double,4> {
165 size_t size()
const {
return nv; }
170 typedef vector<node_type> n_array_t;
171 typedef vector<element_type> e_array_t;
172 typedef vector<size_t> i_array_t;
175 i_array_t node_bamg_dom_id;
177 i_array_t element_bamg_dom_id;
178 e_array_t edge_boundary;
179 i_array_t edge_boundary_bamg_dom_id;
180 list<element_type> edge_list;
182 size_t ntri = 0, nqua = 0;
184 while (
bamg.good() && label !=
"End") {
186 if (label ==
"Vertices") {
195 cerr <<
"bamg2geo: error: empty bamg mesh file" << endl;
199 node_bamg_dom_id.resize (
nnod);
200 for (
size_t inod = 0; inod <
nnod; inod++) {
201 bamg >> node[inod][0]
203 >> node_bamg_dom_id[inod];
205 }
else if (label ==
"Triangles") {
212 cerr <<
"bamg2geo: error: triangles should precede quadrangles" << endl;
216 element.resize (ntri);
217 element_bamg_dom_id.resize (ntri);
218 for (
size_t it = 0; it < ntri; it++) {
219 bamg >> element[it][0]
222 >> element_bamg_dom_id[it];
228 }
else if (label ==
"Quadrilaterals") {
234 cerr <<
"nqua = " << nqua << endl;
235 element.resize (ntri+nqua);
236 element_bamg_dom_id.resize (ntri+nqua);
237 for (
size_t iq = 0; iq < nqua; iq++) {
238 bamg >> element[ntri+iq][0]
239 >> element[ntri+iq][1]
240 >> element[ntri+iq][2]
241 >> element[ntri+iq][3]
242 >> element_bamg_dom_id[ntri+iq];
243 element[ntri+iq][0]--;
244 element[ntri+iq][1]--;
245 element[ntri+iq][2]--;
246 element[ntri+iq][3]--;
247 element[ntri+iq].nv = 4;
249 }
else if (label ==
"Edges") {
257 edge_boundary.resize (nedg);
258 edge_boundary_bamg_dom_id.resize (nedg);
259 for (
size_t iedg = 0; iedg < nedg; iedg++) {
260 bamg >> edge_boundary[iedg][0]
261 >> edge_boundary[iedg][1]
262 >> edge_boundary_bamg_dom_id[iedg];
263 edge_boundary[iedg][0]--;
264 edge_boundary[iedg][1]--;
265 edge_boundary[iedg].nv = 2;
272 vector<string> node_domain_name;
273 vector<string> edge_domain_name;
274 vector<string> region_domain_name;
279 while (
c ==
'E' ||
c ==
'V' ||
c ==
'R') {
282 if (!
scatch(
dmn,
"RegionDomainNames",
true))
break;
285 region_domain_name.resize (n_dom_region);
286 for (
size_t k = 0; k < n_dom_region; k++) {
287 dmn >> region_domain_name[k];
289 }
else if (
c ==
'E') {
290 if (!
scatch(
dmn,
"EdgeDomainNames",
true))
break;
296 edge_domain_name.resize (n_dom_edge);
297 for (
size_t k = 0; k < n_dom_edge; k++) {
298 dmn >> edge_domain_name[k];
300 }
else if (
c ==
'V') {
301 if (!
scatch(
dmn,
"VertexDomainNames",
true))
break;
305 size_t n_dom_vertice;
306 dmn >> n_dom_vertice;
307 node_domain_name.resize (n_dom_vertice);
308 for (
size_t k = 0; k < n_dom_vertice; k++) {
309 dmn >> node_domain_name[k];
317 vector<index_set> ball_edge (node.size());
319 for (
size_t ie = 0; ie < element.size(); ++ie) {
320 for (
size_t iv0 = 0; iv0 < element[ie].nv; ++iv0) {
321 size_t iv1 = (iv0+1) % element[ie].nv;
322 size_t inod0 = element[ie][iv0];
323 size_t inod1 = element[ie][iv1];
326 if (iedge_set.size() > 1) {
327 cerr <<
"bamg2geo: error: connectivity problem (iedge.size="<<iedge_set.size()<<
")" << endl;
330 if (iedge_set.size() == 1)
continue;
331 ball_edge[inod0].insert (nedg);
332 ball_edge[inod1].insert (nedg);
333 element_type new_edge;
337 edge_list.push_back (new_edge);
342 std::copy (edge_list.begin(), edge_list.end(),
edge.begin());
347 typedef pair<size_t,size_t> pair_t;
348 typedef map<size_t, pair_t> map_t;
349 map_t edge_bamg_id2idom;
351 size_t edge_idom = 0;
352 for (
size_t ieb = 0, neb = edge_boundary_bamg_dom_id.size(); ieb < neb; ieb++) {
353 size_t bamg_id = edge_boundary_bamg_dom_id [ieb];
354 typename map_t::iterator iter = edge_bamg_id2idom.find (bamg_id);
355 if (iter != edge_bamg_id2idom.end()) {
357 ((*iter).second.second)++;
361 edge_bamg_id2idom.insert (pair<size_t,pair_t>(bamg_id, pair_t(edge_idom,1)));
364 size_t edge_ndom = edge_bamg_id2idom.size();
366 if (edge_ndom != edge_domain_name.size()) {
367 cerr <<
"bamg2geo: error: "
368 << edge_domain_name.size() <<
" domain name(s) founded while "
369 << edge_ndom <<
" bamg 1d domain(s) are defined" << endl
370 <<
"HINT: check domain name file (.dmn)" << endl;
374 struct edge_orient_type {
size_t index;
int orient; };
375 vector<list<edge_orient_type> > edge_domain (edge_ndom);
376 for (
size_t ieb = 0, neb = edge_boundary_bamg_dom_id.size(); ieb < neb; ieb++) {
377 size_t bamg_dom_id = edge_boundary_bamg_dom_id [ieb];
378 size_t edge_idom = edge_bamg_id2idom [bamg_dom_id].first;
379 size_t inod0 = edge_boundary[ieb][0];
380 size_t inod1 = edge_boundary[ieb][1];
383 if (iedge_set.size() != 1) {
384 cerr <<
"bamg2geo: error: connectivity problem (iedge.size="<<iedge_set.size()<<
")" << endl;
387 size_t ie = *(iedge_set.begin());
390 eo.orient = (inod0 ==
edge[ie][0]) ? 1 : -1;
391 edge_domain[edge_idom].push_back (eo);
396 vector<list<size_t> > node_domain;
397 if (node_domain_name.size() != 0) {
399 std::set<size_t> node_id;
400 for (
size_t iv = 0, nv = node_bamg_dom_id.size(); iv < nv; ++iv) {
401 size_t dom_id = node_bamg_dom_id [iv];
402 if (dom_id == 0)
continue;
403 node_id.insert (dom_id);
405 cerr <<
"node_id1.size = " << node_id.size() << endl;
408 for (
size_t iedg_bdr = 0, nedg_bdr = edge_boundary_bamg_dom_id.size(); iedg_bdr < nedg_bdr; ++iedg_bdr) {
409 size_t dom_id = edge_boundary_bamg_dom_id [iedg_bdr];
410 node_id.erase (dom_id);
412 cerr <<
"node_id2.size = " << node_id.size() << endl;
413 if (node_id.size() != node_domain_name.size()) {
414 cerr <<
"bamg2geo: error: unexpected VertexDomainNames with "<<node_domain_name.size()
415 <<
" domain names while the mesh provides " << node_id.size()
416 <<
" node labels" << endl;
420 node_domain.resize (node_id.size());
421 size_t node_idom = 0;
422 for (set<size_t>::const_iterator iter = node_id.begin(); iter != node_id.end(); ++iter, ++node_idom) {
423 size_t bamg_id = *iter;
424 string name = node_domain_name[node_idom];
425 for (
size_t i = 0,
n = node_bamg_dom_id.size(); i <
n; ++i) {
426 if (node_bamg_dom_id [i] != bamg_id)
continue;
427 node_domain[node_idom].push_back (i);
434 cout <<
"#!geo" << endl
439 <<
" dimension 2" << endl;
440 if (syscoord !=
"cartesian") {
441 cout <<
" coordinate_system " << syscoord << endl;
443 cout <<
" nodes " << node.size() << endl;
445 cout <<
" triangles " << ntri << endl;
448 cout <<
" quadrangles " << nqua << endl;
450 if (
edge.size() != 0) {
451 cout <<
" edges " <<
edge.size() << endl;
453 cout <<
"end header" << endl;
455 cout << setprecision(numeric_limits<double>::digits10);
456 for (
size_t i = 0; i < node.size(); ++i) {
457 cout << node[i][0] <<
" "
458 << node[i][1] << endl;
461 for (
size_t i = 0; i < element.size(); ++i) {
462 cout << ((element[i].nv == 3) ?
't' :
'q')
464 for (
size_t iv = 0; iv < element[i].nv; ++iv) {
465 cout << element[i][iv];
466 if (iv+1 < element[i].nv) { cout <<
" "; }
477 for (
size_t idom = 0; idom < edge_domain.size(); ++idom) {
478 cout <<
"domain" << endl
479 << edge_domain_name[idom] << endl
480 <<
"2 1 " << edge_domain[idom].size() << endl;
481 for (
auto e: edge_domain[idom]) {
482 cout << e.index*e.orient << endl;
486 for (
size_t idom = 0; idom < node_domain.size(); ++idom) {
487 cout <<
"domain" << endl
488 << node_domain_name[idom] << endl
489 <<
"2 0 " << node_domain[idom].size() << endl;
490 for (
auto i: node_domain[idom]) {
497 int main(
int argc,
char **argv) {
498 string syscoord =
"cartesian";
499 string bamg_name, dmn_name;
500 for (
int i = 1; i < argc; i++) {
501 if (strcmp (argv[i],
"-cartesian") == 0) syscoord =
"cartesian";
502 else if (strcmp (argv[i],
"-rz") == 0) syscoord =
"rz";
503 else if (strcmp (argv[i],
"-zr") == 0) syscoord =
"zr";
504 else if (argv[i][0] !=
'-') {
506 if (bamg_name ==
"") bamg_name = argv[i];
507 else if (dmn_name ==
"") dmn_name = argv[i];
509 cerr <<
"bamg2geo: too many file names `" << argv[i] <<
"'" << endl;
513 cerr <<
"bamg2geo: unknown option `" << argv[i] <<
"'" << endl;
517 if (bamg_name ==
"") {
521 ifstream
bamg (bamg_name.c_str());
523 cerr <<
"bamg2geo: unable to read file \""<<bamg_name<<
"\"" << endl; exit (1);
525 if (dmn_name ==
"") {
529 ifstream
dmn (dmn_name.c_str());
531 cerr <<
"bamg2geo: unable to read file \""<<dmn_name<<
"\"" << endl; exit (1);
void bamg2geo(istream &bamg, istream &dmn, string syscoord)
int main(int argc, char **argv)
see the edge page for the full documentation
void inplace_intersection(const index_set &b)
const size_t edge[n_edge][2]
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin color format format bamg
size_type nnod(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
This file is part of Rheolef.
bool scatch(std::istream &in, const std::string &ch, bool full_match=true)
scatch: see the rheostream page for the full documentation