From: wolf Date: Wed, 24 Apr 2002 11:56:38 +0000 (+0000) Subject: Check in for now, although only partly finished. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=be7cfc1ee4480c775227dcd1e5266c8e83506dc1;p=dealii-svn.git Check in for now, although only partly finished. git-svn-id: https://svn.dealii.org/trunk@5723 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 9f9105be06..06a5e79e7d 100644 --- a/deal.II/examples/step-14/step-14.cc +++ b/deal.II/examples/step-14/step-14.cc @@ -152,6 +152,9 @@ namespace Evaluation std::cout << " Point value=" << point_value << ", exact value=1.59492, error=" << 1.594915543-point_value << std::endl; +// std::cout << " Point value=" << point_value //TODO +// << ", exact value=1, error=" +// << 1.-point_value << std::endl; }; @@ -520,7 +523,7 @@ namespace LaplaceSolver const Function &rhs_function, const Function &boundary_values); - // XXX + //TODO!! virtual void solve_problem (); @@ -990,6 +993,7 @@ namespace Data q += sin(10*p(i)+5*p(0)*p(0)); const double exponential = exp(q); return exponential; +// return 0; // TODO! }; @@ -1017,6 +1021,8 @@ namespace Data t1 = t1*t1; return -u*(t1+t2+t3); +// const double pi = 3.1415926536; +// return 2.*pi*pi*sin(pi*p(0))*sin(pi*p(1)); //TODO!! }; @@ -1237,8 +1243,8 @@ namespace Data // And since we want that the // evaluation point (3/4,3/4) in // this example is a grid point, - // we refine once globally: - coarse_grid.refine_global (1); + // we refine twice globally: + coarse_grid.refine_global (2); }; }; @@ -1419,7 +1425,7 @@ namespace LaplaceSolver const Quadrature &face_quadrature, const DualFunctional::DualFunctionalBase &dual_functional); - // XXX + //TODO!! virtual void solve_problem (); @@ -1846,7 +1852,7 @@ namespace LaplaceSolver data_out.add_data_vector (DualSolver::solution, "dual_solution"); - data_out.build_patches (1); + data_out.build_patches (); #ifdef HAVE_STD_STRINGSTREAM std::ostringstream filename; @@ -1910,15 +1916,27 @@ namespace LaplaceSolver FETools::interpolate (PrimalSolver::dof_handler, PrimalSolver::solution, DualSolver::dof_handler, - primal_solution); + primal_solution); + //TODO!! + Vector tmp (PrimalSolver::dof_handler.n_dofs()); + Vector i_h_dual_solution (DualSolver::dof_handler.n_dofs()); + FETools::interpolate (DualSolver::dof_handler, + DualSolver::solution, + PrimalSolver::dof_handler, + tmp); + ConstraintMatrix primal_hanging_node_constraints; + DoFTools::make_hanging_node_constraints (PrimalSolver::dof_handler, + primal_hanging_node_constraints); + primal_hanging_node_constraints.close (); + primal_hanging_node_constraints.distribute (tmp); + FETools::interpolate (PrimalSolver::dof_handler, + tmp, + DualSolver::dof_handler, + i_h_dual_solution); + Vector dual_weights (DualSolver::dof_handler.n_dofs()); -// FETools::interpolation_difference (DualSolver::dof_handler, -// DualSolver::solution, -// *PrimalSolver::fe, -// dual_weights); dual_weights = DualSolver::solution; - abort (); // check Galerkin orthogonality, also for hanging nodes! - + dual_weights -= i_h_dual_solution; // Then we set up a map between // face iterators and their jump @@ -2606,18 +2624,18 @@ run_simulation (LaplaceSolver::Base &solver, template void solve_problem () { - Triangulation triangulation (Triangulation::maximum_smoothing); + Triangulation triangulation (Triangulation::smoothing_on_refinement); const FE_Q primal_fe(1); const FE_Q dual_fe(2); const QGauss4 quadrature; const QGauss4 face_quadrature; const Data::SetUpBase *data = - new Data::SetUp,dim> (); + new Data::SetUp,dim> (); data->create_coarse_grid (triangulation); - const Point evaluation_point(3./4.,3./4.); + const Point evaluation_point(0.5,0.5); const DualFunctional::PointValueEvaluation dual_functional (evaluation_point); @@ -2633,7 +2651,7 @@ void solve_problem () TableHandler results_table; Evaluation::PointValueEvaluation - postprocessor1 (Point(3./4.,3./4.), results_table); + postprocessor1 (Point(0.5,0.5), results_table); Evaluation::GridOutput postprocessor2 ("grid");