--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2015 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// a un-hp-ified version of hp/step-7
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <fstream>
+std::ofstream logfile("output");
+
+
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/function.h>
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/compressed_simple_sparsity_pattern.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/precondition.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_trace.h>
+#include <deal.II/numerics/matrix_tools.h>
+#include <deal.II/numerics/error_estimator.h>
+#include <deal.II/numerics/data_out.h>
+
+#include <deal.II/base/smartpointer.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/base/convergence_table.h>
+#include <deal.II/hp/fe_values.h>
+
+#include <typeinfo>
+#include <fstream>
+#include <iomanip>
+
+
+
+template <int dim>
+class SolutionBase
+{
+protected:
+ static const unsigned int n_source_centers = 3;
+ static const Point<dim> source_centers[n_source_centers];
+ static const double width;
+};
+
+
+template <>
+const Point<1>
+SolutionBase<1>::source_centers[SolutionBase<1>::n_source_centers]
+ = { Point<1>(-1.0 / 3.0),
+ Point<1>(0.0),
+ Point<1>(+1.0 / 3.0)
+ };
+
+template <>
+const Point<2>
+SolutionBase<2>::source_centers[SolutionBase<2>::n_source_centers]
+ = { Point<2>(-0.5, +0.5),
+ Point<2>(-0.5, -0.5),
+ Point<2>(+0.5, -0.5)
+ };
+
+template <>
+const Point<3>
+SolutionBase<3>::source_centers[SolutionBase<3>::n_source_centers]
+ = { Point<3>(-0.5, +0.5, 0.5),
+ Point<3>(-0.5, -0.5, -0.5),
+ Point<3>(+0.5, -0.5, 0)
+ };
+
+template <int dim>
+const double SolutionBase<dim>::width = 1./3.;
+
+
+
+template <int dim>
+class Solution : public Function<dim>,
+ protected SolutionBase<dim>
+{
+public:
+ Solution () : Function<dim>() {}
+
+ virtual double value (const Point<dim> &p,
+ const unsigned int component = 0) const;
+
+ virtual Tensor<1,dim> gradient (const Point<dim> &p,
+ const unsigned int component = 0) const;
+};
+
+
+template <int dim>
+double Solution<dim>::value (const Point<dim> &p,
+ const unsigned int) const
+{
+ double return_value = 0;
+ for (unsigned int i=0; i<this->n_source_centers; ++i)
+ {
+ const Tensor<1,dim> x_minus_xi = p - this->source_centers[i];
+ return_value += std::exp(-x_minus_xi.norm_square() /
+ (this->width * this->width));
+ }
+
+ return return_value;
+}
+
+
+template <int dim>
+Tensor<1,dim> Solution<dim>::gradient (const Point<dim> &p,
+ const unsigned int) const
+{
+ Tensor<1,dim> return_value;
+
+ for (unsigned int i=0; i<this->n_source_centers; ++i)
+ {
+ const Tensor<1,dim> x_minus_xi = p - this->source_centers[i];
+
+ return_value += (-2 / (this->width * this->width) *
+ std::exp(-x_minus_xi.norm_square() /
+ (this->width * this->width)) *
+ x_minus_xi);
+ }
+
+ return return_value;
+}
+
+
+
+template <int dim>
+class RightHandSide : public Function<dim>,
+ protected SolutionBase<dim>
+{
+public:
+ RightHandSide () : Function<dim>() {}
+
+ virtual double value (const Point<dim> &p,
+ const unsigned int component = 0) const;
+};
+
+
+template <int dim>
+double RightHandSide<dim>::value (const Point<dim> &p,
+ const unsigned int) const
+{
+ double return_value = 0;
+ for (unsigned int i=0; i<this->n_source_centers; ++i)
+ {
+ const Tensor<1,dim> x_minus_xi = p - this->source_centers[i];
+
+ return_value += ((2*dim - 4*x_minus_xi.norm_square()/
+ (this->width * this->width)) /
+ (this->width * this->width) *
+ std::exp(-x_minus_xi.norm_square() /
+ (this->width * this->width)));
+ return_value += std::exp(-x_minus_xi.norm_square() /
+ (this->width * this->width));
+ }
+
+ return return_value;
+}
+
+
+
+template <int dim>
+class HelmholtzProblem
+{
+public:
+ enum RefinementMode
+ {
+ global_refinement, adaptive_refinement
+ };
+
+ HelmholtzProblem (const unsigned int fe_degree,
+ const RefinementMode refinement_mode);
+
+ void run ();
+
+private:
+ void setup_system ();
+ void assemble_system (const bool do_reconstruct);
+ void solve ();
+ void refine_grid ();
+ void process_solution (const unsigned int cycle);
+
+ Triangulation<dim> triangulation;
+ FE_Q<dim> fe;
+ DoFHandler<dim> dof_handler;
+ ConstraintMatrix constraints;
+ SparsityPattern sparsity_pattern;
+ SparseMatrix<double> system_matrix;
+ Vector<double> solution;
+ Vector<double> system_rhs;
+
+ FE_TraceQ<dim> fe_trace;
+ DoFHandler<dim> dof_handler_trace;
+ ConstraintMatrix constraints_trace;
+ SparsityPattern sparsity_pattern_trace;
+ SparseMatrix<double> system_matrix_trace;
+ Vector<double> solution_trace_full;
+ Vector<double> solution_trace;
+ Vector<double> system_rhs_trace;
+
+ const RefinementMode refinement_mode;
+
+ ConvergenceTable convergence_table;
+};
+
+
+
+
+template <int dim>
+HelmholtzProblem<dim>::HelmholtzProblem (const unsigned int fe_degree,
+ const RefinementMode refinement_mode) :
+ fe (fe_degree),
+ dof_handler (triangulation),
+ fe_trace (fe_degree),
+ dof_handler_trace (triangulation),
+ refinement_mode (refinement_mode)
+{
+ deallog << "Solving with Q" << fe_degree << " elements, "
+ << (refinement_mode == global_refinement ? "global" : "adaptive")
+ << " refinement" << std::endl
+ << "==========================================="
+ << (refinement_mode == global_refinement ? "" : "==") << std::endl
+ << std::endl;
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::setup_system ()
+{
+ dof_handler.distribute_dofs (fe);
+
+ constraints.clear ();
+ DoFTools::make_hanging_node_constraints (dof_handler,
+ constraints);
+ VectorTools::interpolate_boundary_values (dof_handler,
+ 0,
+ Solution<dim>(),
+ constraints);
+ constraints.close ();
+
+ {
+ CompressedSimpleSparsityPattern csp (dof_handler.n_dofs());
+ DoFTools::make_sparsity_pattern (dof_handler, csp, constraints, false);
+ sparsity_pattern.copy_from(csp);
+ }
+
+ system_matrix.reinit (sparsity_pattern);
+
+ solution.reinit (dof_handler.n_dofs());
+ system_rhs.reinit (dof_handler.n_dofs());
+
+
+ dof_handler_trace.distribute_dofs (fe_trace);
+
+ constraints_trace.clear ();
+ DoFTools::make_hanging_node_constraints (dof_handler_trace,
+ constraints_trace);
+ VectorTools::interpolate_boundary_values (dof_handler_trace,
+ 0,
+ Solution<dim>(),
+ constraints_trace);
+ constraints_trace.close ();
+
+ {
+ CompressedSimpleSparsityPattern csp (dof_handler_trace.n_dofs());
+ DoFTools::make_sparsity_pattern (dof_handler_trace, csp, constraints_trace, false);
+ sparsity_pattern_trace.copy_from(csp);
+ }
+
+ system_matrix_trace.reinit (sparsity_pattern_trace);
+
+ solution_trace.reinit (dof_handler_trace.n_dofs());
+ system_rhs_trace.reinit (dof_handler_trace.n_dofs());
+ solution_trace_full.reinit (dof_handler.n_dofs());
+
+ deallog << "Number of DoFs: " << dof_handler.n_dofs() << " / "
+ << dof_handler_trace.n_dofs() << std::endl;
+ deallog << "Number of matrix nnz: " << sparsity_pattern.n_nonzero_elements() << " / "
+ << sparsity_pattern_trace.n_nonzero_elements() << std::endl;
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::assemble_system (const bool do_reconstruct)
+{
+ QGauss<dim> quadrature_formula (fe.degree+1);
+ QGauss<dim-1> face_quadrature_formula (fe.degree+1);
+
+ const unsigned int n_q_points = quadrature_formula.size();
+ const unsigned int n_face_q_points = face_quadrature_formula.size();
+
+ const unsigned int dofs_per_cell = fe.dofs_per_cell;
+
+ FullMatrix<double> cell_matrix (dofs_per_cell, dofs_per_cell);
+ Vector<double> cell_rhs (dofs_per_cell);
+
+ FullMatrix<double> trace_matrix (fe_trace.dofs_per_cell, fe_trace.dofs_per_cell);
+ Vector<double> trace_rhs (fe_trace.dofs_per_cell);
+
+ FullMatrix<double> eliminate_matrix(cell_matrix.m()-trace_matrix.m(),
+ cell_matrix.m()-trace_matrix.m());
+ FullMatrix<double> temp_matrix(trace_matrix.m(), eliminate_matrix.m());
+
+ std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
+ std::vector<types::global_dof_index> dof_indices_trace (fe_trace.dofs_per_cell);
+ Vector<double> local_trace(fe_trace.dofs_per_cell);
+ Vector<double> local_interior(dofs_per_cell-fe_trace.dofs_per_cell);
+
+ FEValues<dim> x_fe_values (fe, quadrature_formula,
+ update_values | update_gradients |
+ update_q_points | update_JxW_values);
+
+ FEFaceValues<dim> x_fe_face_values (fe, face_quadrature_formula,
+ update_values | update_q_points |
+ update_normal_vectors | update_JxW_values);
+
+ const RightHandSide<dim> right_hand_side;
+ std::vector<double> rhs_values (n_q_points);
+
+ const Solution<dim> exact_solution;
+
+ typename DoFHandler<dim>::active_cell_iterator
+ cell = dof_handler.begin_active(),
+ endc = dof_handler.end(),
+ tracec = dof_handler_trace.begin_active();
+ for (; cell!=endc; ++cell, ++tracec)
+ {
+ cell_rhs = 0;
+
+ x_fe_values.reinit (cell);
+ const FEValues<dim> &fe_values = x_fe_values.get_present_fe_values();
+
+ right_hand_side.value_list (fe_values.get_quadrature_points(),
+ rhs_values);
+
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ {
+ for (unsigned int j=0; j<dofs_per_cell; ++j)
+ {
+ double sum = 0;
+ const Tensor<1,dim> *shape_grad_i = &fe_values.shape_grad(i,0);
+ const Tensor<1,dim> *shape_grad_j = &fe_values.shape_grad(j,0);
+ const double *shape_value_i = &fe_values.shape_value(i,0);
+ const double *shape_value_j = &fe_values.shape_value(j,0);
+ for (unsigned int q_index=0; q_index<n_q_points; ++q_index)
+ sum += (shape_grad_i[q_index] * shape_grad_j[q_index] +
+ shape_value_i[q_index] * shape_value_j[q_index]
+ ) * fe_values.JxW(q_index);
+ cell_matrix(i,j) = sum;
+ }
+
+ for (unsigned int q_point = 0; q_point < n_q_points; ++q_point)
+ cell_rhs(i) += (fe_values.shape_value(i,q_point) *
+ rhs_values [q_point] *
+ fe_values.JxW(q_point));
+ }
+
+ for (unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
+ if (cell->face(face)->at_boundary()
+ &&
+ (cell->face(face)->boundary_indicator() == 1))
+ {
+ x_fe_face_values.reinit (cell, face);
+ const FEFaceValues<dim> &fe_face_values
+ = x_fe_face_values.get_present_fe_values ();
+
+ for (unsigned int q_point=0; q_point<n_face_q_points; ++q_point)
+ {
+ const double neumann_value
+ = (exact_solution.gradient (fe_face_values.quadrature_point(q_point)) *
+ fe_face_values.normal_vector(q_point));
+
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ cell_rhs(i) += (neumann_value *
+ fe_face_values.shape_value(i,q_point) *
+ fe_face_values.JxW(q_point));
+ }
+ }
+
+ const unsigned int shift = fe_trace.dofs_per_cell;
+ const unsigned int sizes = fe.dofs_per_cell - fe_trace.dofs_per_cell;
+ for (unsigned int i=0; i<sizes; ++i)
+ for (unsigned int j=0; j<sizes; ++j)
+ eliminate_matrix(i,j) = cell_matrix(i+shift, j+shift);
+ if (sizes > 0)
+ eliminate_matrix.gauss_jordan();
+
+ cell->get_dof_indices (local_dof_indices);
+ if (do_reconstruct == false)
+ {
+ cell->get_dof_indices (local_dof_indices);
+ constraints.distribute_local_to_global(cell_matrix, cell_rhs,
+ local_dof_indices,
+ system_matrix, system_rhs);
+
+ for (unsigned int i=0; i<shift; ++i)
+ for (unsigned int j=0; j<sizes; ++j)
+ {
+ double sum = 0;
+ for (unsigned int k=0; k<sizes; ++k)
+ sum += cell_matrix(i,k+shift) * eliminate_matrix(k,j);
+ temp_matrix(i,j) = sum;
+ }
+ for (unsigned int i=0; i<shift; ++i)
+ for (unsigned int j=0; j<shift; ++j)
+ {
+ double sum = 0;
+ for (unsigned int k=0; k<sizes; ++k)
+ sum += temp_matrix(i,k) * cell_matrix(shift+k,j);
+ trace_matrix(i,j) = cell_matrix(i,j) - sum;
+ }
+ for (unsigned int i=0; i<shift; ++i)
+ {
+ double sum = 0;
+ for (unsigned int k=0; k<sizes; ++k)
+ sum += temp_matrix(i,k) * cell_rhs(shift+k);
+ trace_rhs(i) = cell_rhs(i) - sum;
+ }
+
+ tracec->get_dof_indices (dof_indices_trace);
+ constraints_trace.distribute_local_to_global(trace_matrix, trace_rhs,
+ dof_indices_trace,
+ system_matrix_trace,
+ system_rhs_trace);
+ }
+ else
+ {
+ tracec->get_interpolated_dof_values(solution_trace, local_trace);
+ for (unsigned int i=0; i<shift; ++i)
+ solution_trace_full(local_dof_indices[i]) = local_trace(i);
+ for (unsigned int i=0; i<sizes; ++i)
+ {
+ double sum = 0;
+ for (unsigned int k=0; k<shift; ++k)
+ sum += cell_matrix(shift+i,k) * local_trace(k);
+ local_interior(i) = cell_rhs(shift+i) - sum;
+ }
+ for (unsigned int i=0; i<sizes; ++i)
+ {
+ double sum = 0;
+ for (unsigned int j=0; j<sizes; ++j)
+ sum += eliminate_matrix(i,j) * local_interior(j);
+ solution_trace_full(local_dof_indices[shift+i]) = sum;
+ }
+ }
+ }
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::solve ()
+{
+ {
+ SolverControl solver_control (1000, 1e-12);
+ SolverCG<> cg (solver_control);
+
+ PreconditionSSOR<> preconditioner;
+ preconditioner.initialize(system_matrix, 1.2);
+
+ cg.solve (system_matrix, solution, system_rhs,
+ preconditioner);
+
+ constraints.distribute (solution);
+ }
+ {
+ SolverControl solver_control (1000, 1e-12);
+ SolverCG<> cg (solver_control);
+
+ PreconditionSSOR<> preconditioner;
+ preconditioner.initialize(system_matrix_trace, 1.2);
+
+ cg.solve (system_matrix_trace, solution_trace, system_rhs_trace,
+ preconditioner);
+
+ constraints_trace.distribute (solution_trace);
+ }
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::refine_grid ()
+{
+ switch (refinement_mode)
+ {
+ case global_refinement:
+ {
+ triangulation.refine_global (1);
+ break;
+ }
+
+ case adaptive_refinement:
+ {
+ Vector<float> estimated_error_per_cell (triangulation.n_active_cells());
+
+ typename FunctionMap<dim>::type neumann_boundary;
+ KellyErrorEstimator<dim>::estimate (dof_handler,
+ QGauss<dim-1>(3),
+ neumann_boundary,
+ solution,
+ estimated_error_per_cell);
+
+ GridRefinement::refine_and_coarsen_fixed_number (triangulation,
+ estimated_error_per_cell,
+ 0.3, 0.03);
+
+ triangulation.execute_coarsening_and_refinement ();
+
+ break;
+ }
+
+ default:
+ {
+ Assert (false, ExcNotImplemented());
+ }
+ }
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::process_solution (const unsigned int cycle)
+{
+ {
+ Vector<float> difference_per_cell (triangulation.n_active_cells());
+ VectorTools::integrate_difference (dof_handler,
+ solution,
+ Solution<dim>(),
+ difference_per_cell,
+ QGauss<dim>(fe.degree+2),
+ VectorTools::L2_norm);
+ const double L2_error = difference_per_cell.l2_norm();
+
+ VectorTools::integrate_difference (dof_handler,
+ solution,
+ Solution<dim>(),
+ difference_per_cell,
+ QGauss<dim>(fe.degree+2),
+ VectorTools::H1_seminorm);
+ const double H1_error = difference_per_cell.l2_norm();
+
+ const unsigned int n_active_cells=triangulation.n_active_cells();
+ const unsigned int n_dofs=dof_handler.n_dofs();
+
+ convergence_table.add_value("cycle", cycle);
+ convergence_table.add_value("cells", n_active_cells);
+ convergence_table.add_value("dofs", n_dofs);
+ convergence_table.add_value("L2", L2_error);
+ convergence_table.add_value("H1", H1_error);
+ }
+ {
+ Vector<float> difference_per_cell (triangulation.n_active_cells());
+ VectorTools::integrate_difference (dof_handler,
+ solution_trace_full,
+ Solution<dim>(),
+ difference_per_cell,
+ QGauss<dim>(fe.degree+2),
+ VectorTools::L2_norm);
+ const double L2_error = difference_per_cell.l2_norm();
+
+ VectorTools::integrate_difference (dof_handler,
+ solution_trace_full,
+ Solution<dim>(),
+ difference_per_cell,
+ QGauss<dim>(fe.degree+2),
+ VectorTools::H1_seminorm);
+ const double H1_error = difference_per_cell.l2_norm();
+
+ const unsigned int n_active_cells=triangulation.n_active_cells();
+ const unsigned int n_dofs=dof_handler.n_dofs();
+
+ convergence_table.add_value("L2 sc", L2_error);
+ convergence_table.add_value("H1 sc", H1_error);
+ }
+}
+
+
+
+template <int dim>
+void HelmholtzProblem<dim>::run ()
+{
+ for (unsigned int cycle=0; cycle<6-dim; ++cycle)
+ {
+ if (cycle == 0)
+ {
+ GridGenerator::hyper_cube (triangulation, -1, 1);
+ triangulation.refine_global (1);
+
+ typename Triangulation<dim>::cell_iterator
+ cell = triangulation.begin (),
+ endc = triangulation.end();
+ for (; cell!=endc; ++cell)
+ for (unsigned int face=0;
+ face<GeometryInfo<dim>::faces_per_cell;
+ ++face)
+ if ((cell->face(face)->center()(0) == -1)
+ ||
+ (cell->face(face)->center()(1) == -1))
+ cell->face(face)->set_boundary_indicator (1);
+ }
+ else
+ refine_grid ();
+
+
+ setup_system ();
+
+ assemble_system (false);
+ solve ();
+ assemble_system (true);
+
+ process_solution (cycle);
+ }
+
+ convergence_table.set_precision("L2", 3);
+ convergence_table.set_precision("H1", 3);
+
+ convergence_table.set_scientific("L2", true);
+ convergence_table.set_scientific("H1", true);
+
+ convergence_table.set_precision("L2 sc", 3);
+ convergence_table.set_precision("H1 sc", 3);
+
+ convergence_table.set_scientific("L2 sc", true);
+ convergence_table.set_scientific("H1 sc", true);
+
+ deallog << std::endl;
+ convergence_table.write_text(deallog.get_file_stream());
+}
+
+
+int main ()
+{
+ deallog << std::setprecision(2);
+ logfile << std::setprecision(2);
+
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ try
+ {
+ deallog.depth_console (0);
+
+
+ for (unsigned int deg = 1; deg < 5; ++deg)
+ {
+ HelmholtzProblem<2>
+ helmholtz_problem_2d (deg, HelmholtzProblem<2>::adaptive_refinement);
+
+ helmholtz_problem_2d.run ();
+
+ deallog << std::endl;
+ }
+ for (unsigned int deg = 1; deg < 4; ++deg)
+ {
+ HelmholtzProblem<2>
+ helmholtz_problem_2d (deg, HelmholtzProblem<2>::global_refinement);
+
+ helmholtz_problem_2d.run ();
+
+ deallog << std::endl;
+ }
+
+ {
+ HelmholtzProblem<3>
+ helmholtz_problem_3d (2, HelmholtzProblem<3>::adaptive_refinement);
+
+ helmholtz_problem_3d.run ();
+
+ deallog << std::endl;
+ }
+ }
+ catch (std::exception &exc)
+ {
+ deallog << std::endl << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ deallog << "Exception on processing: " << std::endl
+ << exc.what() << std::endl
+ << "Aborting!" << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ return 1;
+ }
+ catch (...)
+ {
+ deallog << std::endl << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ deallog << "Unknown exception!" << std::endl
+ << "Aborting!" << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ return 1;
+ }
+
+ return 0;
+}
+
+
+template const double SolutionBase<2>::width;
--- /dev/null
+
+DEAL::Solving with Q1 elements, adaptive refinement
+DEAL::=============================================
+DEAL::
+DEAL::Number of DoFs: 9 / 9
+DEAL::Number of matrix nnz: 21 / 21
+DEAL:cg::Starting value 5.5
+DEAL:cg::Convergence step 4 value 0
+DEAL:cg::Starting value 5.5
+DEAL:cg::Convergence step 4 value 0
+DEAL::Number of DoFs: 22 / 22
+DEAL::Number of matrix nnz: 86 / 86
+DEAL:cg::Starting value 6.0
+DEAL:cg::Convergence step 11 value 0
+DEAL:cg::Starting value 6.0
+DEAL:cg::Convergence step 11 value 0
+DEAL::Number of DoFs: 46 / 46
+DEAL::Number of matrix nnz: 222 / 222
+DEAL:cg::Starting value 5.3
+DEAL:cg::Convergence step 14 value 0
+DEAL:cg::Starting value 5.3
+DEAL:cg::Convergence step 14 value 0
+DEAL::Number of DoFs: 84 / 84
+DEAL::Number of matrix nnz: 444 / 444
+DEAL:cg::Starting value 3.7
+DEAL:cg::Convergence step 17 value 0
+DEAL:cg::Starting value 3.7
+DEAL:cg::Convergence step 17 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 9 4.037e+00 4.590e+00 4.037e+00 4.590e+00
+1 13 22 5.358e-01 1.430e+00 5.358e-01 1.430e+00
+2 31 46 3.087e-01 1.258e+00 3.087e-01 1.258e+00
+3 61 84 7.744e-02 1.096e+00 7.744e-02 1.096e+00
+DEAL::
+DEAL::Solving with Q2 elements, adaptive refinement
+DEAL::=============================================
+DEAL::
+DEAL::Number of DoFs: 25 / 21
+DEAL::Number of matrix nnz: 153 / 107
+DEAL:cg::Starting value 11.
+DEAL:cg::Convergence step 11 value 0
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 11 value 0
+DEAL::Number of DoFs: 71 / 58
+DEAL::Number of matrix nnz: 643 / 460
+DEAL:cg::Starting value 3.9
+DEAL:cg::Convergence step 19 value 0
+DEAL:cg::Starting value 3.9
+DEAL:cg::Convergence step 16 value 0
+DEAL::Number of DoFs: 148 / 120
+DEAL::Number of matrix nnz: 1492 / 1062
+DEAL:cg::Starting value 3.4
+DEAL:cg::Convergence step 27 value 0
+DEAL:cg::Starting value 3.6
+DEAL:cg::Convergence step 21 value 0
+DEAL::Number of DoFs: 265 / 210
+DEAL::Number of matrix nnz: 3101 / 2240
+DEAL:cg::Starting value 2.6
+DEAL:cg::Convergence step 39 value 0
+DEAL:cg::Starting value 2.6
+DEAL:cg::Convergence step 28 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 25 1.471e+00 2.988e+00 1.471e+00 2.988e+00
+1 13 71 9.942e-02 1.127e+00 9.942e-02 1.127e+00
+2 28 148 7.311e-02 9.526e-01 7.311e-02 9.526e-01
+3 55 265 1.091e-02 2.577e-01 1.091e-02 2.577e-01
+DEAL::
+DEAL::Solving with Q3 elements, adaptive refinement
+DEAL::=============================================
+DEAL::
+DEAL::Number of DoFs: 49 / 33
+DEAL::Number of matrix nnz: 589 / 261
+DEAL:cg::Starting value 3.9
+DEAL:cg::Convergence step 20 value 0
+DEAL:cg::Starting value 0.30
+DEAL:cg::Convergence step 14 value 0
+DEAL::Number of DoFs: 146 / 94
+DEAL::Number of matrix nnz: 2378 / 1122
+DEAL:cg::Starting value 3.7
+DEAL:cg::Convergence step 32 value 0
+DEAL:cg::Starting value 3.3
+DEAL:cg::Convergence step 22 value 0
+DEAL::Number of DoFs: 342 / 218
+DEAL::Number of matrix nnz: 6126 / 2902
+DEAL:cg::Starting value 3.0
+DEAL:cg::Convergence step 47 value 0
+DEAL:cg::Starting value 3.0
+DEAL:cg::Convergence step 30 value 0
+DEAL::Number of DoFs: 630 / 386
+DEAL::Number of matrix nnz: 12314 / 5938
+DEAL:cg::Starting value 1.8
+DEAL:cg::Convergence step 65 value 0
+DEAL:cg::Starting value 2.0
+DEAL:cg::Convergence step 39 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 49 3.566e-01 1.675e+00 3.566e-01 1.675e+00
+1 13 146 1.336e-02 2.093e-01 1.336e-02 2.093e-01
+2 31 342 7.824e-03 1.526e-01 7.824e-03 1.526e-01
+3 61 630 1.181e-03 3.811e-02 1.181e-03 3.811e-02
+DEAL::
+DEAL::Solving with Q4 elements, adaptive refinement
+DEAL::=============================================
+DEAL::
+DEAL::Number of DoFs: 81 / 45
+DEAL::Number of matrix nnz: 1617 / 483
+DEAL:cg::Starting value 5.9
+DEAL:cg::Convergence step 30 value 0
+DEAL:cg::Starting value 0.75
+DEAL:cg::Convergence step 17 value 0
+DEAL::Number of DoFs: 247 / 130
+DEAL::Number of matrix nnz: 6311 / 2072
+DEAL:cg::Starting value 2.9
+DEAL:cg::Convergence step 49 value 0
+DEAL:cg::Starting value 2.8
+DEAL:cg::Convergence step 26 value 0
+DEAL::Number of DoFs: 583 / 304
+DEAL::Number of matrix nnz: 16071 / 5334
+DEAL:cg::Starting value 2.5
+DEAL:cg::Convergence step 73 value 0
+DEAL:cg::Starting value 2.5
+DEAL:cg::Convergence step 35 value 0
+DEAL::Number of DoFs: 1113 / 564
+DEAL::Number of matrix nnz: 32325 / 10914
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 98 value 0
+DEAL:cg::Starting value 1.6
+DEAL:cg::Convergence step 46 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 81 4.573e-02 5.591e-01 4.573e-02 5.591e-01
+1 13 247 3.694e-03 8.184e-02 3.694e-03 8.184e-02
+2 31 583 2.238e-03 5.906e-02 2.238e-03 5.906e-02
+3 61 1113 1.549e-04 6.438e-03 1.549e-04 6.438e-03
+DEAL::
+DEAL::Solving with Q1 elements, global refinement
+DEAL::===========================================
+DEAL::
+DEAL::Number of DoFs: 9 / 9
+DEAL::Number of matrix nnz: 21 / 21
+DEAL:cg::Starting value 5.5
+DEAL:cg::Convergence step 4 value 0
+DEAL:cg::Starting value 5.5
+DEAL:cg::Convergence step 4 value 0
+DEAL::Number of DoFs: 25 / 25
+DEAL::Number of matrix nnz: 109 / 109
+DEAL:cg::Starting value 6.0
+DEAL:cg::Convergence step 11 value 0
+DEAL:cg::Starting value 6.0
+DEAL:cg::Convergence step 11 value 0
+DEAL::Number of DoFs: 81 / 81
+DEAL::Number of matrix nnz: 501 / 501
+DEAL:cg::Starting value 3.8
+DEAL:cg::Convergence step 16 value 0
+DEAL:cg::Starting value 3.8
+DEAL:cg::Convergence step 16 value 0
+DEAL::Number of DoFs: 289 / 289
+DEAL::Number of matrix nnz: 2149 / 2149
+DEAL:cg::Starting value 2.3
+DEAL:cg::Convergence step 26 value 0
+DEAL:cg::Starting value 2.3
+DEAL:cg::Convergence step 26 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 9 4.037e+00 4.590e+00 4.037e+00 4.590e+00
+1 16 25 4.926e-01 1.409e+00 4.926e-01 1.409e+00
+2 64 81 7.984e-02 1.130e+00 7.984e-02 1.130e+00
+3 256 289 2.101e-02 5.828e-01 2.101e-02 5.828e-01
+DEAL::
+DEAL::Solving with Q2 elements, global refinement
+DEAL::===========================================
+DEAL::
+DEAL::Number of DoFs: 25 / 21
+DEAL::Number of matrix nnz: 153 / 107
+DEAL:cg::Starting value 11.
+DEAL:cg::Convergence step 11 value 0
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 11 value 0
+DEAL::Number of DoFs: 81 / 65
+DEAL::Number of matrix nnz: 801 / 575
+DEAL:cg::Starting value 3.8
+DEAL:cg::Convergence step 20 value 0
+DEAL:cg::Starting value 3.8
+DEAL:cg::Convergence step 17 value 0
+DEAL::Number of DoFs: 289 / 225
+DEAL::Number of matrix nnz: 3633 / 2639
+DEAL:cg::Starting value 2.6
+DEAL:cg::Convergence step 40 value 0
+DEAL:cg::Starting value 2.6
+DEAL:cg::Convergence step 29 value 0
+DEAL::Number of DoFs: 1089 / 833
+DEAL::Number of matrix nnz: 15441 / 11279
+DEAL:cg::Starting value 1.4
+DEAL:cg::Convergence step 79 value 0
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 53 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 25 1.471e+00 2.988e+00 1.471e+00 2.988e+00
+1 16 81 9.523e-02 1.122e+00 9.523e-02 1.122e+00
+2 64 289 9.409e-03 2.484e-01 9.409e-03 2.484e-01
+3 256 1089 1.194e-03 6.228e-02 1.194e-03 6.228e-02
+DEAL::
+DEAL::Solving with Q3 elements, global refinement
+DEAL::===========================================
+DEAL::
+DEAL::Number of DoFs: 49 / 33
+DEAL::Number of matrix nnz: 589 / 261
+DEAL:cg::Starting value 3.9
+DEAL:cg::Convergence step 20 value 0
+DEAL:cg::Starting value 0.30
+DEAL:cg::Convergence step 14 value 0
+DEAL::Number of DoFs: 169 / 105
+DEAL::Number of matrix nnz: 2941 / 1397
+DEAL:cg::Starting value 3.7
+DEAL:cg::Convergence step 34 value 0
+DEAL:cg::Starting value 3.3
+DEAL:cg::Convergence step 23 value 0
+DEAL::Number of DoFs: 625 / 369
+DEAL::Number of matrix nnz: 13045 / 6381
+DEAL:cg::Starting value 1.8
+DEAL:cg::Convergence step 63 value 0
+DEAL:cg::Starting value 2.0
+DEAL:cg::Convergence step 38 value 0
+DEAL::Number of DoFs: 2401 / 1377
+DEAL::Number of matrix nnz: 54853 / 27197
+DEAL:cg::Starting value 1.1
+DEAL:cg::Convergence step 120 value 0
+DEAL:cg::Starting value 1.2
+DEAL:cg::Convergence step 69 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 4 49 3.566e-01 1.675e+00 3.566e-01 1.675e+00
+1 16 169 9.770e-03 1.813e-01 9.770e-03 1.813e-01
+2 64 625 9.013e-04 3.472e-02 9.013e-04 3.472e-02
+3 256 2401 6.648e-05 5.054e-03 6.648e-05 5.054e-03
+DEAL::
+DEAL::Solving with Q2 elements, adaptive refinement
+DEAL::=============================================
+DEAL::
+DEAL::Number of DoFs: 125 / 117
+DEAL::Number of matrix nnz: 1085 / 893
+DEAL:cg::Starting value 7.8
+DEAL:cg::Convergence step 11 value 0
+DEAL:cg::Starting value 3.6
+DEAL:cg::Convergence step 10 value 0
+DEAL::Number of DoFs: 407 / 378
+DEAL::Number of matrix nnz: 6443 / 5454
+DEAL:cg::Starting value 5.2
+DEAL:cg::Convergence step 18 value 0
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 15 value 0
+DEAL::Number of DoFs: 1198 / 1099
+DEAL::Number of matrix nnz: 30428 / 26217
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 26 value 0
+DEAL:cg::Starting value 1.5
+DEAL:cg::Convergence step 21 value 0
+DEAL::
+cycle cells dofs L2 H1 L2 sc H1 sc
+0 8 125 5.433e-01 2.365e+00 5.433e-01 2.365e+00
+1 29 407 2.488e-01 1.565e+00 2.488e-01 1.565e+00
+2 99 1198 5.615e-02 7.982e-01 5.615e-02 7.982e-01
+DEAL::