From: Wolfgang Bangerth Date: Tue, 2 May 2006 23:07:24 +0000 (+0000) Subject: New test, but doesn't work right now. X-Git-Tag: v8.0.0~11799 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c8deffbadbd81f8a08783182c8f9633866621f9b;p=dealii.git New test, but doesn't work right now. git-svn-id: https://svn.dealii.org/trunk@12996 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/hp/step-3.cc b/tests/hp/step-3.cc new file mode 100644 index 0000000000..68a32f2a7e --- /dev/null +++ b/tests/hp/step-3.cc @@ -0,0 +1,226 @@ +//---------------------------- step-3.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-3.cc --------------------------- + + +// a hp-ified version of step-3 + + +#include +#include +#include +#include + +#include +#include +#include + +#include + +#include + +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include +#include +#include + +std::ofstream logfile("step-3/output"); + + + + +class LaplaceProblem +{ + public: + LaplaceProblem (); + + void run (); + + private: + void make_grid_and_dofs (); + void assemble_system (); + void solve (); + void output_results () const; + + Triangulation<2> triangulation; + hp::FECollection<2> fe; + hp::DoFHandler<2> dof_handler; + + SparsityPattern sparsity_pattern; + SparseMatrix system_matrix; + + Vector solution; + Vector system_rhs; +}; + + +LaplaceProblem::LaplaceProblem () : + fe (FE_Q<2>(1)), + dof_handler (triangulation) +{} + + + +void LaplaceProblem::make_grid_and_dofs () +{ + GridGenerator::hyper_cube (triangulation, -1, 1); + triangulation.refine_global (5); + deallog << "Number of active cells: " + << triangulation.n_active_cells() + << std::endl; + deallog << "Total number of cells: " + << triangulation.n_cells() + << std::endl; + + 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()); +} + + + +void LaplaceProblem::assemble_system () +{ + hp::QCollection<2> quadrature_formula(QGauss<2>(2)); + hp::FEValues<2> x_fe_values (fe, quadrature_formula, + update_values | update_gradients | 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); + + hp::DoFHandler<2>::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + { + x_fe_values.reinit (cell); + + const FEValues<2> &fe_values = x_fe_values.get_present_fe_values(); + + cell_matrix = 0; + cell_rhs = 0; + + for (unsigned int i=0; iget_dof_indices (local_dof_indices); + + for (unsigned int i=0; i boundary_values; + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ZeroFunction<2>(), + boundary_values); + MatrixTools::apply_boundary_values (boundary_values, + system_matrix, + solution, + system_rhs); +} + + + +void LaplaceProblem::solve () +{ + SolverControl solver_control (1000, 1e-12); + SolverCG<> cg (solver_control); + + cg.solve (system_matrix, solution, system_rhs, + PreconditionIdentity()); +} + + + +void LaplaceProblem::output_results () const +{ + DataOut<2,hp::DoFHandler<2> > data_out; + data_out.attach_dof_handler (dof_handler); + data_out.add_data_vector (solution, "solution"); + data_out.build_patches (); + + std::ofstream output ("solution.gpl"); + data_out.write_gnuplot (output); +} + + + +void LaplaceProblem::run () +{ + make_grid_and_dofs (); + assemble_system (); + solve (); + output_results (); +} + + + +int main () +{ + logfile.precision(2); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + LaplaceProblem laplace_problem; + laplace_problem.run (); + + return 0; +}