-The first method uses a triangulated approximation of the circle with
-unit radius and integrates the function that is constant one over it. Of course, if
-the domain were the exact unit circle, then the area would be $\pi$, but
-since we only use an approximation by piecewise polynomial segments,
-the value of the area we integrate over is not exactly $\pi$. However, it is known that as
-we refine the triangulation, a $Q_p$ mapping approximates the boundary
-with an order $h^{p+1}$, where $h$ is the mesh
-size. We will check the values of the computed area of the circle and
-their convergence towards $\pi$ under mesh refinement for different
-mappings. We will also find a convergence behavior that is surprising
-at first, but has a good explanation.
+The first method uses a triangulated approximation of the circle with unit
+radius and integrates a unit magnitude constant function ($f = 1$) over it. Of
+course, if the domain were the exact unit circle, then the area would be $\pi$,
+but since we only use an approximation by piecewise polynomial segments, the
+value of the area we integrate over is not exactly $\pi$. However, it is known
+that as we refine the triangulation, a $Q_p$ mapping approximates the boundary
+with an order $h^{p+1}$, where $h$ is the mesh size. We will check the values
+of the computed area of the circle and their convergence towards $\pi$ under
+mesh refinement for different mappings. We will also find a convergence
+behavior that is surprising at first, but has a good explanation.
// of the library). Rather, we just pack the functionality into separate
// functions. We make these functions templates on the number of space
// dimensions to conform to usual practice when using deal.II, although we
- // will only use them for two space dimensions.
+ // will only use them for two space dimensions and throw an exception when
+ // attempted to use for any other spatial dimension.
//
// The first of these functions just generates a triangulation of a circle
// (hyperball) and outputs the $Q_p$ mapping of its cells for different values
// the triangulation, with $w(x_i)$ being the weight of quadrature point
// $x_i$. The integrals on each cell are approximated by numerical
// quadrature, hence the only additional ingredient we need is to set up a
- // FEValues object that provides the corresponding `JxW' values of each
- // cell. (Note that `JxW' is meant to abbreviate <i>Jacobian determinant
+ // FEValues object that provides the corresponding `JxW` values of each
+ // cell. (Note that `JxW` is meant to abbreviate <i>Jacobian determinant
// times weight</i>; since in numerical quadrature the two factors always
// occur at the same places, we only offer the combined quantity, rather
// than two separate ones.) We note that here we won't use the FEValues
// object in its original purpose, i.e. for the computation of values of
// basis functions of a specific finite element at certain quadrature
- // points. Rather, we use it only to gain the `JxW' at the quadrature
+ // points. Rather, we use it only to gain the `JxW` at the quadrature
// points, irrespective of the (dummy) finite element we will give to the
// constructor of the FEValues object. The actual finite element given to
// the FEValues object is not used at all, so we could give any.
// We now create a finite element. Unlike the rest of the example
// programs, we do not actually need to do any computations with shape
- // functions; we only need the `JxW' values from an FEValues
+ // functions; we only need the `JxW` values from an FEValues
// object. Hence we use the special finite element class FE_Nothing
// which has exactly zero degrees of freedom per cell (as the name
// implies, the local basis on each cell is the empty set). A more
const FE_Nothing<dim> fe;
// Likewise, we need to create a DoFHandler object. We do not actually
- // use it, but it will provide us with `active_cell_iterators' that
+ // use it, but it will provide us with `active_cell_iterators` that
// are needed to reinitialize the FEValues object on each cell of the
// triangulation.
DoFHandler<dim> dof_handler(triangulation);
// Now we set up the FEValues object, giving the Mapping, the dummy
// finite element and the quadrature object to the constructor,
- // together with the update flags asking for the `JxW' values at the
+ // together with the update flags asking for the `JxW` values at the
// quadrature points only. This tells the FEValues object that it
// needs not compute other quantities upon calling the
// <code>reinit</code> function, thus saving computation time.
{
// In this loop we first add the number of active cells of the
// current triangulation to the table. This function automatically
- // creates a table column with superscription `cells', in case
+ // creates a table column with superscription `cells`, in case
// this column was not created before.
table.add_value("cells", triangulation.n_active_cells());
// function below.
dof_handler.distribute_dofs(fe);
- // We define the variable area as `long double' like we did for
- // the pi variable before.
+ // We define the variable area as `long double` like we did for
+ // the `pi` variable before.
long double area = 0;
// Now we loop over all cells, reinitialize the FEValues object
- // for each cell, and add up all the `JxW' values for this cell to
- // `area'...
+ // for each cell, and add up all the `JxW` values for this cell to
+ // `area`...
for (const auto &cell : dof_handler.active_cell_iterators())
{
fe_values.reinit(cell);
table.add_value("error", static_cast<double>(std::fabs(area - pi)));
}
- // We want to compute the convergence rates of the `error'
+ // We want to compute the convergence rates of the `error`
// column. Therefore we need to omit the other columns from the
// convergence rate evaluation before calling
- // `evaluate_all_convergence_rates'
+ // `evaluate_all_convergence_rates`
table.omit_column_from_convergence_rate_evaluation("cells");
table.omit_column_from_convergence_rate_evaluation("eval.pi");
table.evaluate_all_convergence_rates(
std::cout << "Computation of Pi by the perimeter:" << std::endl
<< "===================================" << std::endl;
- // We take the same order of quadrature but this time a `dim-1'
+ // We take the same order of quadrature but this time a `dim-1`
// dimensional quadrature as we will integrate over (boundary) lines
// rather than over cells.
const QGauss<dim - 1> quadrature(4);
dof_handler.distribute_dofs(fe);
// Now we run over all cells and over all faces of each cell. Only
- // the contributions of the `JxW' values on boundary faces are
- // added to the long double variable `perimeter'.
+ // the contributions of the `JxW` values on boundary faces are
+ // added to the long double variable `perimeter`.
long double perimeter = 0;
for (const auto &cell : dof_handler.active_cell_iterators())
for (const auto &face : cell->face_iterators())
{
std::cout.precision(16);
- Step10::gnuplot_output<2>();
+ const unsigned int dim = 2;
- Step10::compute_pi_by_area<2>();
- Step10::compute_pi_by_perimeter<2>();
+ Step10::gnuplot_output<dim>();
+
+ Step10::compute_pi_by_area<dim>();
+ Step10::compute_pi_by_perimeter<dim>();
}
catch (std::exception &exc)
{
projects / dealii.git / commitdiff