From 73684b653e580ccc9185318a3d50bea2fae1d592 Mon Sep 17 00:00:00 2001 From: Alexander Grayver Date: Fri, 20 Mar 2020 19:57:42 +0100 Subject: [PATCH] Add mesh images --- examples/step-11/doc/results.dox | 21 +++++++++++++ examples/step-11/step-11.cc | 54 ++++++++++++++++++++++++++++++++ 2 files changed, 75 insertions(+) diff --git a/examples/step-11/doc/results.dox b/examples/step-11/doc/results.dox index ebb209867b..ec32c39377 100644 --- a/examples/step-11/doc/results.dox +++ b/examples/step-11/doc/results.dox @@ -42,3 +42,24 @@ of convergence but just to reduce the constant before the convergence order. On the other hand, using a cubic mapping only improves the result further insignificantly, except for the case of very coarse grids. + +We can also visualize the underlying meshes by using, for instance +ParaView. The image below shows initial meshes for different mapping +degrees: + + + +Clearly, the effect is most pronounced when we go from the linear to +quadratic mapping. This is also reflected in the error values given +in the table above. The effect of going from quadratic to cubic degree +is less dramatic, but still tangible owning to a more accurace +description of the circular boundary. + +Next, let's look at the meshes after three global refinements + + + +Here, the differences are much less visible, especially for higher order +mappings. Indeed, at this refinement level the error values reported +in the table are essentially identical between mappings of degree two +and three. diff --git a/examples/step-11/step-11.cc b/examples/step-11/step-11.cc index a97816dbeb..10625546f4 100644 --- a/examples/step-11/step-11.cc +++ b/examples/step-11/step-11.cc @@ -41,6 +41,7 @@ #include #include #include +#include // Just this one is new: it declares a class // DynamicSparsityPattern, which we will use and explain @@ -83,6 +84,7 @@ namespace Step11 void setup_system(); void assemble_and_solve(); void solve(); + void write_high_order_mesh(const unsigned cycle); Triangulation triangulation; FE_Q fe; @@ -397,6 +399,57 @@ namespace Step11 + // Next, we write the solution as well as the + // material ids to a VTU file. This is similar to what was done in many + // other tutorial programs. The new ingredient presented in this tutorial + // program is that we want to ensure that the data written to the file + // used for visualization is actually a faithful representation of what + // is used internally by deal.II. That is because most of the visualization + // data formats only represent cells by their vertex coordinates, but + // have no way of representing the curved boundaries that are used + // in deal.II when using higher order mappings -- in other words, what + // you see in the visualization tool is not actually what you are computing + // on. (The same, incidentally, is true when using higher order shape + // functions: Most visualization tools only render bilinear/trilinear + // representations. This is discussed in detail in DataOut::build_patches().) + // + // So we need to ensure that a high-order representation is written + // to the file. We need to consider two particular topics. Firstly, we tell + // the DataOut object via the DataOutBase::VtkFlags that we intend to + // interpret the subdivisions of the elements as a high-order Lagrange + // polynomial rather than a collection of bilinear patches. + // Recent visualization programs, like ParaView version 5.5 + // or newer, can then render a high-order solution (see a wiki + // page for more details). Secondly, we need to make sure that the mapping + // is passed to the DataOut::build_patches() method. Finally, the DataOut + // class only prints curved faces for boundary cells by default, so we + // need to ensure that also inner cells are printed in a curved representation + // via the mapping. + template + void LaplaceProblem::write_high_order_mesh(const unsigned cycle) + { + DataOut data_out; + + DataOutBase::VtkFlags flags; + flags.write_higher_order_cells = true; + data_out.set_flags(flags); + + data_out.attach_dof_handler(dof_handler); + data_out.add_data_vector(solution, "solution"); + + data_out.build_patches(mapping, + mapping.get_degree(), + DataOut::curved_inner_cells); + + std::ofstream file(("solution-c=" + std::to_string(cycle) + + ".p=" + std::to_string(mapping.get_degree()) + ".vtu") + .c_str()); + + data_out.write_vtu(file); + } + + // Finally the main function controlling the different steps to be // performed. Its content is rather straightforward, generating a // triangulation of a circle, associating a boundary to it, and then doing @@ -419,6 +472,7 @@ namespace Step11 { setup_system(); assemble_and_solve(); + write_high_order_mesh(cycle); triangulation.refine_global(); } -- 2.39.5