From: Wolfgang Bangerth Date: Wed, 3 May 2006 02:26:26 +0000 (+0000) Subject: Add another test X-Git-Tag: v8.0.0~11792 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=698081da3b321d4197dce32ed716847132330430;p=dealii.git Add another test git-svn-id: https://svn.dealii.org/trunk@13009 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/hp/step-5.cc b/tests/hp/step-5.cc new file mode 100644 index 0000000000..de9ac2b34f --- /dev/null +++ b/tests/hp/step-5.cc @@ -0,0 +1,339 @@ +//---------------------------- step-5.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2005, 2006 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- step-5.cc --------------------------- + + +// a hp-ified version of step-5 + + +#include +#include +std::ofstream logfile("step-5/output"); + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + +#include +#include + + + +template +class LaplaceProblem +{ + public: + LaplaceProblem (); + void run (); + + private: + void setup_system (); + void assemble_system (); + void solve (); + void output_results (const unsigned int cycle) const; + + Triangulation triangulation; + hp::FECollection fe; + hp::DoFHandler dof_handler; + + SparsityPattern sparsity_pattern; + SparseMatrix system_matrix; + + Vector solution; + Vector system_rhs; +}; + + + +template +class Coefficient : public Function +{ + public: + Coefficient () : Function() {}; + + virtual double value (const Point &p, + const unsigned int component = 0) const; + + virtual void value_list (const std::vector > &points, + std::vector &values, + const unsigned int component = 0) const; +}; + + + +template +double Coefficient::value (const Point &p, + const unsigned int /*component*/) const +{ + if (p.square() < 0.5*0.5) + return 20; + else + return 1; +} + + + + +template +void Coefficient::value_list (const std::vector > &points, + std::vector &values, + const unsigned int component) const +{ + Assert (values.size() == points.size(), + ExcDimensionMismatch (values.size(), points.size())); + + Assert (component == 0, + ExcIndexRange (component, 0, 1)); + + const unsigned int n_points = points.size(); + + for (unsigned int i=0; i +LaplaceProblem::LaplaceProblem () : + fe (FE_Q(1)), + dof_handler (triangulation) +{} + + + + +template +void LaplaceProblem::setup_system () +{ + dof_handler.distribute_dofs (fe); + + deallog << " Number of degrees of freedom: " + << dof_handler.n_dofs() + << std::endl; + + sparsity_pattern.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern); + sparsity_pattern.compress(); + + system_matrix.reinit (sparsity_pattern); + + solution.reinit (dof_handler.n_dofs()); + system_rhs.reinit (dof_handler.n_dofs()); +} + + + + +template +void LaplaceProblem::assemble_system () +{ + hp::QCollection quadrature_formula(QGauss(2)); + + hp::FEValues x_fe_values (fe, quadrature_formula, + update_values | update_gradients | + update_q_points | update_JxW_values); + + const unsigned int dofs_per_cell = fe[0].dofs_per_cell; + const unsigned int n_q_points = quadrature_formula[0].n_quadrature_points; + + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + Vector cell_rhs (dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + const Coefficient coefficient; + std::vector coefficient_values (n_q_points); + + typename hp::DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + { + cell_matrix = 0; + cell_rhs = 0; + + x_fe_values.reinit (cell); + const FEValues<2> &fe_values = x_fe_values.get_present_fe_values(); + + coefficient.value_list (fe_values.get_quadrature_points(), + coefficient_values); + + for (unsigned int q_point=0; q_pointget_dof_indices (local_dof_indices); + for (unsigned int i=0; i boundary_values; + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ZeroFunction(), + boundary_values); + MatrixTools::apply_boundary_values (boundary_values, + system_matrix, + solution, + system_rhs); +} + + + +template +void LaplaceProblem::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); + + deallog << " " << solver_control.last_step() + << " CG iterations needed to obtain convergence." + << std::endl; +} + + + +template +void LaplaceProblem::output_results (const unsigned int cycle) const +{ + DataOut > data_out; + + data_out.attach_dof_handler (dof_handler); + data_out.add_data_vector (solution, "solution"); + + data_out.build_patches (); + + DataOutBase::EpsFlags eps_flags; + eps_flags.z_scaling = 4; + eps_flags.azimut_angle = 40; + eps_flags.turn_angle = 10; + data_out.set_flags (eps_flags); + + data_out.write_eps (deallog.get_file_stream()); +} + + + + +template +void LaplaceProblem::run () +{ + for (unsigned int cycle=0; cycle<6; ++cycle) + { + deallog << "Cycle " << cycle << ':' << std::endl; + + if (cycle != 0) + triangulation.refine_global (1); + else + { + GridIn grid_in; + grid_in.attach_triangulation (triangulation); + std::ifstream input_file("step-5/circle-grid.inp"); + Assert (dim==2, ExcInternalError()); + + grid_in.read_ucd (input_file); + + static const HyperBallBoundary boundary; + triangulation.set_boundary (0, boundary); + } + + deallog << " Number of active cells: " + << triangulation.n_active_cells() + << std::endl + << " Total number of cells: " + << triangulation.n_cells() + << std::endl; + + setup_system (); + assemble_system (); + solve (); + output_results (cycle); + } +} + + + +int main () +{ + logfile.precision(2); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + deallog.depth_console (0); + + LaplaceProblem<2> laplace_problem_2d; + laplace_problem_2d.run (); + +/* + Coefficient<2> coefficient; + std::vector > points (2); + std::vector coefficient_values (1); + coefficient.value_list (points, coefficient_values); +*/ + + return 0; +} diff --git a/tests/hp/step-5/circle-grid.inp b/tests/hp/step-5/circle-grid.inp new file mode 100644 index 0000000000..f28a7a238a --- /dev/null +++ b/tests/hp/step-5/circle-grid.inp @@ -0,0 +1,46 @@ +25 20 0 0 0 +1 -0.7071 -0.7071 0 +2 0.7071 -0.7071 0 +3 -0.2668 -0.2668 0 +4 0.2668 -0.2668 0 +5 -0.2668 0.2668 0 +6 0.2668 0.2668 0 +7 -0.7071 0.7071 0 +8 0.7071 0.7071 0 +9 0 -1 0 +10 0.5 -0.5 0 +11 0 -0.3139 0 +12 -0.5 -0.5 0 +13 0 -0.6621 0 +14 -0.3139 0 0 +15 -0.5 0.5 0 +16 -1 0 0 +17 -0.6621 0 0 +18 0.3139 0 0 +19 0 0.3139 0 +20 0 0 0 +21 1 0 0 +22 0.5 0.5 0 +23 0.6621 0 0 +24 0 1 0 +25 0 0.6621 0 +1 0 quad 1 9 13 12 +2 0 quad 9 2 10 13 +3 0 quad 13 10 4 11 +4 0 quad 12 13 11 3 +5 0 quad 1 12 17 16 +6 0 quad 12 3 14 17 +7 0 quad 17 14 5 15 +8 0 quad 16 17 15 7 +9 0 quad 3 11 20 14 +10 0 quad 11 4 18 20 +11 0 quad 20 18 6 19 +12 0 quad 14 20 19 5 +13 0 quad 2 21 23 10 +14 0 quad 21 8 22 23 +15 0 quad 23 22 6 18 +16 0 quad 10 23 18 4 +17 0 quad 7 15 25 24 +18 0 quad 15 5 19 25 +19 0 quad 25 19 6 22 +20 0 quad 24 25 22 8 diff --git a/tests/hp/step-5/cmp/generic b/tests/hp/step-5/cmp/generic new file mode 100644 index 0000000000..309c3cb828 --- /dev/null +++ b/tests/hp/step-5/cmp/generic @@ -0,0 +1,43 @@ + +DEAL::Cycle 0: +DEAL:: Number of active cells: 20 +DEAL:: Total number of cells: 20 +DEAL:: Number of degrees of freedom: 25 +DEAL:cg::Starting value 0.51 +DEAL:cg::Convergence step 13 value 0 +DEAL:: 13 CG iterations needed to obtain convergence. +DEAL::Cycle 1: +DEAL:: Number of active cells: 80 +DEAL:: Total number of cells: 100 +DEAL:: Number of degrees of freedom: 89 +DEAL:cg::Starting value 0.31 +DEAL:cg::Convergence step 18 value 0 +DEAL:: 18 CG iterations needed to obtain convergence. +DEAL::Cycle 2: +DEAL:: Number of active cells: 320 +DEAL:: Total number of cells: 420 +DEAL:: Number of degrees of freedom: 337 +DEAL:cg::Starting value 0.17 +DEAL:cg::Convergence step 29 value 0 +DEAL:: 29 CG iterations needed to obtain convergence. +DEAL::Cycle 3: +DEAL:: Number of active cells: 1280 +DEAL:: Total number of cells: 1700 +DEAL:: Number of degrees of freedom: 1313 +DEAL:cg::Starting value 0.092 +DEAL:cg::Convergence step 52 value 0 +DEAL:: 52 CG iterations needed to obtain convergence. +DEAL::Cycle 4: +DEAL:: Number of active cells: 5120 +DEAL:: Total number of cells: 6820 +DEAL:: Number of degrees of freedom: 5185 +DEAL:cg::Starting value 0.047 +DEAL:cg::Convergence step 95 value 0 +DEAL:: 95 CG iterations needed to obtain convergence. +DEAL::Cycle 5: +DEAL:: Number of active cells: 20480 +DEAL:: Total number of cells: 27300 +DEAL:: Number of degrees of freedom: 20609 +DEAL:cg::Starting value 0.024 +DEAL:cg::Convergence step 182 value 0 +DEAL:: 182 CG iterations needed to obtain convergence.