From: wolf Date: Fri, 19 Apr 2002 14:23:08 +0000 (+0000) Subject: Preliminary import, to work on another computer. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6c9af7db445166624e08eef5e72ad70bb2d63d13;p=dealii-svn.git Preliminary import, to work on another computer. git-svn-id: https://svn.dealii.org/trunk@5694 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-14/step-14.cc b/deal.II/examples/step-14/step-14.cc index cc79c0902a..a4ab26478b 100644 --- a/deal.II/examples/step-14/step-14.cc +++ b/deal.II/examples/step-14/step-14.cc @@ -1,10 +1,10 @@ /* $Id$ */ -/* Author: Wolfgang Bangerth, University of Heidelberg, 2001, 2002 */ +/* Author: Wolfgang Bangerth, ETH Zurich, 2002 */ /* $Id$ */ /* Version: $Name$ */ /* */ -/* Copyright (C) 2001, 2002 by the deal.II authors */ +/* Copyright (C) 2002 by the deal.II authors */ /* */ /* This file is subject to QPL and may not be distributed */ /* without copyright and license information. Please refer */ @@ -508,6 +508,20 @@ namespace LaplaceSolver const Quadrature &quadrature, const Function &rhs_function, const Function &boundary_values); + + // XXX + virtual + void + solve_problem (); + + virtual + unsigned int + n_dofs () const; + + virtual + void + postprocess (const Evaluation::EvaluationBase &postprocessor) const; + protected: const SmartPointer > rhs_function; virtual void assemble_rhs (Vector &rhs) const; @@ -529,6 +543,32 @@ namespace LaplaceSolver {}; + template + void + PrimalSolver::solve_problem () + { + Solver::solve_problem (); + }; + + + + template + unsigned int + PrimalSolver::n_dofs() const + { + return Solver::n_dofs(); + }; + + + template + void + PrimalSolver:: + postprocess (const Evaluation::EvaluationBase &postprocessor) const + { + return Solver::postprocess(postprocessor); + }; + + template void @@ -729,6 +769,233 @@ RightHandSide::value (const Point &p, +namespace DualFunctional +{ + template + class DualFunctionalBase : public Subscriptor + { + public: + virtual + void + assemble_rhs (const DoFHandler &dof_handler, + Vector &rhs) const = 0; + }; + + + template + class PointValueEvaluation : public DualFunctionalBase + { + public: + PointValueEvaluation (const Point &evaluation_point); + + virtual + void + assemble_rhs (const DoFHandler &dof_handler, + Vector &rhs) const; + protected: + const Point evaluation_point; + }; + + + template + PointValueEvaluation:: + PointValueEvaluation (const Point &evaluation_point) + : + evaluation_point (evaluation_point) + {}; + + + template + void + PointValueEvaluation:: + assemble_rhs (const DoFHandler &dof_handler, + Vector &rhs) const + { + rhs.reinit (dof_handler.n_dofs()); + typename DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + bool evaluation_point_found = false; + for (; (cell!=endc) && !evaluation_point_found; ++cell) + for (unsigned int vertex=0; + vertex::vertices_per_cell; + ++vertex) + if (cell->vertex(vertex) == evaluation_point) + { + rhs(cell->vertex_dof_index(vertex,0)) = 1; + + evaluation_point_found = true; + break; + }; + + AssertThrow (evaluation_point_found, + ExcEvaluationPointNotFound(evaluation_point)); + }; + + +}; + + +namespace LaplaceSolver +{ + template + class DualSolver : public Solver + { + public: + DualSolver (Triangulation &triangulation, + const FiniteElement &fe, + const Quadrature &quadrature, + const DualFunctional::DualFunctionalBase &dual_functional); + + // XXX + virtual + void + solve_problem (); + + virtual + unsigned int + n_dofs () const; + + virtual + void + postprocess (const Evaluation::EvaluationBase &postprocessor) const; + + protected: + const SmartPointer > dual_functional; + virtual void assemble_rhs (Vector &rhs) const; + + static const ZeroFunction boundary_values; + }; + + + template + DualSolver:: + DualSolver (Triangulation &triangulation, + const FiniteElement &fe, + const Quadrature &quadrature, + const DualFunctional::DualFunctionalBase &dual_functional) + : + Base (triangulation), + Solver (triangulation, fe, + quadrature, boundary_values), + dual_functional (&dual_functional) + {}; + + + template + void + DualSolver::solve_problem () + { + Solver::solve_problem (); + }; + + + + template + unsigned int + DualSolver::n_dofs() const + { + return Solver::n_dofs(); + }; + + + template + void + DualSolver:: + postprocess (const Evaluation::EvaluationBase &postprocessor) const + { + return Solver::postprocess(postprocessor); + }; + + + + template + void + DualSolver:: + assemble_rhs (Vector &rhs) const + { + dual_functional->assemble_rhs (dof_handler, rhs); + }; + + + template + class WeightedResidual : public PrimalSolver, + public DualSolver + { + public: + WeightedResidual (Triangulation &coarse_grid, + const FiniteElement &fe, + const Quadrature &quadrature, + const Function &rhs_function, + const Function &boundary_values, + const DualFunctional::DualFunctionalBase &dual_functional); + + virtual + void + solve_problem (); + + virtual + void + postprocess (const Evaluation::EvaluationBase &postprocessor) const; + + virtual + unsigned int + n_dofs () const; + + virtual void refine_grid () {}; + }; + + + + template + WeightedResidual:: + WeightedResidual (Triangulation &coarse_grid, + const FiniteElement &fe, + const Quadrature &quadrature, + const Function &rhs_function, + const Function &boundary_values, + const DualFunctional::DualFunctionalBase &dual_functional) + : + Base (coarse_grid), + PrimalSolver (coarse_grid, fe, quadrature, + rhs_function, boundary_values), + DualSolver (coarse_grid, fe, quadrature, + dual_functional) + {}; + + + template + void + WeightedResidual::solve_problem () + { + PrimalSolver::solve_problem (); + DualSolver::solve_problem (); + }; + + + template + void + WeightedResidual::postprocess (const Evaluation::EvaluationBase &postprocessor) const + { + PrimalSolver::postprocess (postprocessor); + }; + + + template + unsigned int + WeightedResidual::n_dofs () const + { + return PrimalSolver::n_dofs(); + }; + + + template class WeightedResidual<2>; + +}; + + + + template void run_simulation (LaplaceSolver::Base &solver,