From: Ralf Hartmann Date: Tue, 21 Sep 2004 07:19:16 +0000 (+0000) Subject: Make it compile also in 2d. X-Git-Tag: v8.0.0~14818 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fd97d5de3270f3a0b9df27f77fa1423bec27db46;p=dealii.git Make it compile also in 2d. git-svn-id: https://svn.dealii.org/trunk@9651 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/grid/grid_out.h b/deal.II/deal.II/include/grid/grid_out.h index 4163e4385d..a5e4519633 100644 --- a/deal.II/deal.II/include/grid/grid_out.h +++ b/deal.II/deal.II/include/grid/grid_out.h @@ -637,10 +637,17 @@ class GridOut * documentation of the * GridOutFlags::Gnuplot() class. */ - template - void write_gnuplot (const Triangulation &tria, - std::ostream &out, - const Mapping *mapping=0); + void write_gnuplot (const Triangulation<3> &tria, + std::ostream &out, + const Mapping<3> *mapping=0); + + /** + * Specialization of above + * function for 2d + */ + void write_gnuplot (const Triangulation<2> &tria, + std::ostream &out, + const Mapping<2> *mapping=0); /** * Declaration of the specialization diff --git a/deal.II/deal.II/source/grid/grid_out.cc b/deal.II/deal.II/source/grid/grid_out.cc index 04a1e7f964..043a1f27ec 100644 --- a/deal.II/deal.II/source/grid/grid_out.cc +++ b/deal.II/deal.II/source/grid/grid_out.cc @@ -602,23 +602,147 @@ void GridOut::write_gnuplot (const Triangulation<1> &tria, AssertThrow (out, ExcIO()); } +#endif -#else +#if deal_II_dimension==2 -template -void GridOut::write_gnuplot (const Triangulation &tria, - std::ostream &out, - const Mapping *mapping) +void GridOut::write_gnuplot (const Triangulation<2> &tria, + std::ostream &out, + const Mapping<2> *mapping) { AssertThrow (out, ExcIO()); + const unsigned int dim=2; const unsigned int n_additional_points= gnuplot_flags.n_boundary_face_points; const unsigned int n_points=2+n_additional_points; - typename Triangulation::active_cell_iterator cell=tria.begin_active(); - const typename Triangulation::active_cell_iterator endc=tria.end(); + Triangulation::active_cell_iterator cell=tria.begin_active(); + const Triangulation::active_cell_iterator endc=tria.end(); + + // if we are to treat curved + // boundaries, then generate a + // quadrature formula which will be + // used to probe boundary points at + // curved faces + Quadrature *q_projector=0; + std::vector > boundary_points; + if (mapping!=0) + { + boundary_points.resize(n_points); + boundary_points[0][0]=0; + boundary_points[n_points-1][0]=1; + for (unsigned int i=1; i dummy_weights(n_points, 1./n_points); + Quadrature<1> quadrature(boundary_points, dummy_weights); + + q_projector = new Quadrature (QProjector::project_to_all_faces(quadrature)); + } + + for (; cell!=endc; ++cell) + { + if (gnuplot_flags.write_cell_numbers) + out << "# cell " << cell << std::endl; + + if (mapping==0 || !cell->at_boundary()) + { + // write out the four sides + // of this cell by putting + // the four points (+ the + // initial point again) in + // a row and lifting the + // drawing pencil at the + // end + for (unsigned int i=0; i::vertices_per_cell; ++i) + out << cell->vertex(i) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + out << cell->vertex(0) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl + << std::endl // double new line for gnuplot 3d plots + << std::endl; + } + else + // cell is at boundary and we + // are to treat curved + // boundaries. so loop over + // all faces and draw them as + // small pieces of lines + { + for (unsigned int face_no=0; + face_no::faces_per_cell; ++face_no) + { + const Triangulation::face_iterator + face = cell->face(face_no); + if (face->at_boundary()) + { + // compute offset + // of quadrature + // points within + // set of projected + // points + const unsigned int offset=face_no*n_points; + for (unsigned int i=0; itransform_unit_to_real_cell + (cell, q_projector->point(offset+i))) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) + << std::endl; + + out << std::endl + << std::endl; + } + else + { + // if, however, the + // face is not at + // the boundary, + // then draw it as + // usual + out << face->vertex(0) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) + << std::endl + << face->vertex(1) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) + << std::endl + << std::endl + << std::endl; + } + } + } + } + + if (q_projector != 0) + delete q_projector; + + + AssertThrow (out, ExcIO()); +} + +#endif + + +#if deal_II_dimension==3 + +void GridOut::write_gnuplot (const Triangulation<3> &tria, + std::ostream &out, + const Mapping<3> *mapping) +{ + AssertThrow (out, ExcIO()); + + const unsigned int dim=3; + const unsigned int n_additional_points= + gnuplot_flags.n_boundary_face_points; + const unsigned int n_points=2+n_additional_points; + + Triangulation::active_cell_iterator cell=tria.begin_active(); + const Triangulation::active_cell_iterator endc=tria.end(); // if we are to treat curved // boundaries, then generate a @@ -649,224 +773,105 @@ void GridOut::write_gnuplot (const Triangulation &tria, if (gnuplot_flags.write_cell_numbers) out << "# cell " << cell << std::endl; - switch (dim) - { - case 1: - { - Assert(false, ExcInternalError()); - break; - }; - - case 2: + if (mapping==0 || n_points==2 || !cell->has_boundary_lines()) + for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) { - if (mapping==0 || !cell->at_boundary()) - { - // write out the four - // sides of this cell - // by putting the - // four points (+ the - // initial point - // again) in a row - // and lifting the - // drawing pencil at - // the end - out << cell->vertex(0) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << cell->vertex(1) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << cell->vertex(2) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << cell->vertex(3) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << cell->vertex(0) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << std::endl // double new line for gnuplot 3d plots - << std::endl; - } - else - // cell is at boundary - // and we are to treat - // curved - // boundaries. so loop - // over all faces and - // draw them as small - // pieces of lines - { - for (unsigned int face_no=0; - face_no::faces_per_cell; ++face_no) - { - const typename Triangulation::face_iterator - face = cell->face(face_no); - if (face->at_boundary()) - { - // compute - // offset of - // quadrature - // points - // within set - // of - // projected - // points. note - // that we - // need not - // care about - // reverted - // faces in - // 3d since - // boundary - // faces - // always - // have to be - // in - // standard - // orientation - // and the - // standard - // orientation - // quadrature - // points are - // first in - // the list - const unsigned int offset=face_no*n_points; - for (unsigned int i=0; itransform_unit_to_real_cell - (cell, q_projector->point(offset+i))) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) - << std::endl; - - out << std::endl - << std::endl; - } - else + const Triangulation::face_iterator + face = cell->face(face_no); + + out << "#face" << face_no << std::endl; + for (unsigned int v=0; v::vertices_per_face; ++v) + out << face->vertex(v) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + out << face->vertex(0) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + out << std::endl; + out << std::endl; + } + else + { + for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) + { + const Triangulation::face_iterator + face = cell->face(face_no); + + if (face->at_boundary()) + { + const unsigned int offset=face_no*n_points*n_points; + for (unsigned int i=0; ivertex(0) + const Point p0=mapping->transform_unit_to_real_cell( + cell, q_projector->point(offset+i*n_points+j)); + out << p0 << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) - << std::endl - << face->vertex(1) + << ' ' << static_cast(cell->material_id()) << std::endl; + out << (mapping->transform_unit_to_real_cell( + cell, q_projector->point(offset+(i+1)*n_points+j))) << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) - << std::endl - << std::endl - << std::endl; - }; - }; - }; - - break; - }; - - case 3: - { - if (mapping==0 || n_points==2 || !cell->has_boundary_lines()) - for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) + << ' ' << static_cast(cell->material_id()) << std::endl; + out << (mapping->transform_unit_to_real_cell( + cell, q_projector->point(offset+(i+1)*n_points+j+1))) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + out << (mapping->transform_unit_to_real_cell( + cell, q_projector->point(offset+i*n_points+j+1))) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + // and the + // first + // point + // again + out << p0 + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + out << std::endl << std::endl; + } + } + else { - const typename Triangulation::face_iterator - face = cell->face(face_no); + for (unsigned int l=0; l::lines_per_face; ++l) + { + const Triangulation::line_iterator + line=face->line(l); + + const Point &v0=line->vertex(0), + &v1=line->vertex(1); + if (line->at_boundary()) + { + // transform_real_to_unit_cell + // could be + // replaced + // by using + // QProjector<3>::project_to_line + // which is + // not yet + // implemented + const Point u0=mapping->transform_real_to_unit_cell(cell, v0), + u1=mapping->transform_real_to_unit_cell(cell, v1); + + for (unsigned int i=0; itransform_unit_to_real_cell + (cell, (1-boundary_points[i][0])*u0+boundary_points[i][0]*u1)) + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; + } + else + out << v0 + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl + << v1 + << ' ' << cell->level() + << ' ' << static_cast(cell->material_id()) << std::endl; - out << "#face" << face_no << std::endl; - for (unsigned int v=0; v::vertices_per_face; ++v) - out << face->vertex(v) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << face->vertex(0) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << std::endl; - out << std::endl; + out << std::endl << std::endl; + } } - else - { - for (unsigned int face_no=0; face_no::faces_per_cell; ++face_no) - { - const typename Triangulation::face_iterator - face = cell->face(face_no); - - if (face->at_boundary()) - { - const unsigned int offset=face_no*n_points*n_points; - for (unsigned int i=0; i p0=mapping->transform_unit_to_real_cell( - cell, q_projector->point(offset+i*n_points+j)); - out << p0 - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << (mapping->transform_unit_to_real_cell( - cell, q_projector->point(offset+(i+1)*n_points+j))) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << (mapping->transform_unit_to_real_cell( - cell, q_projector->point(offset+(i+1)*n_points+j+1))) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << (mapping->transform_unit_to_real_cell( - cell, q_projector->point(offset+i*n_points+j+1))) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << p0 - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << std::endl << std::endl; - } - } - else - { - for (unsigned int l=0; l::lines_per_face; ++l) - { - const typename Triangulation::line_iterator - line=face->line(l); - - const Point &v0=line->vertex(0), - &v1=line->vertex(1); - if (line->at_boundary()) - { - // transform_real_to_unit_cell - // could be replaced by using - // QProjector<3>::project_to_line - // which is not yet implemented - const Point u0=mapping->transform_real_to_unit_cell(cell, v0), - u1=mapping->transform_real_to_unit_cell(cell, v1); - - for (unsigned int i=0; itransform_unit_to_real_cell - (cell, (1-boundary_points[i][0])*u0+boundary_points[i][0]*u1)) - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - } - else - out << v0 - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl - << v1 - << ' ' << cell->level() - << ' ' << static_cast(cell->material_id()) << std::endl; - out << std::endl << std::endl; - } - } - } - } - - break; - }; - }; - }; + } + } + } if (q_projector != 0) delete q_projector; @@ -1373,10 +1378,6 @@ template void GridOut::write_ucd (const Triangulation &, std::ostream &); #if deal_II_dimension != 1 -template void GridOut::write_gnuplot -(const Triangulation &, - std::ostream &, - const Mapping *); template void GridOut::write_eps (const Triangulation &, std::ostream &,