From: Wolfgang Bangerth Date: Wed, 29 Apr 1998 12:09:38 +0000 (+0000) Subject: Use functionality to evaluate finite element functions and fix bugs in error estimator. X-Git-Tag: v8.0.0~23030 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8fa5356a4772f016c6fb602c9c6fb9494986c2e7;p=dealii.git Use functionality to evaluate finite element functions and fix bugs in error estimator. git-svn-id: https://svn.dealii.org/trunk@224 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/numerics/assembler.cc b/deal.II/deal.II/source/numerics/assembler.cc index a88b02e521..475545b464 100644 --- a/deal.II/deal.II/source/numerics/assembler.cc +++ b/deal.II/deal.II/source/numerics/assembler.cc @@ -101,7 +101,8 @@ void Assembler::assemble (const Equation &equation) { // re-init fe values for this cell fe_values.reinit (DoFHandler::cell_iterator (tria, present_level, - present_index), + present_index, + dof_handler), fe, boundary); const unsigned int n_dofs = dof_handler->get_selected_fe().total_dofs; diff --git a/deal.II/deal.II/source/numerics/base.cc b/deal.II/deal.II/source/numerics/base.cc index da58fa61af..b411e92bcd 100644 --- a/deal.II/deal.II/source/numerics/base.cc +++ b/deal.II/deal.II/source/numerics/base.cc @@ -211,12 +211,7 @@ void ProblemBase::integrate_difference (const Function &exact_sol // reference_function(x_j)-\psi_j // and assign that to the vector // \psi. - const unsigned int n_dofs = fe.total_dofs; const unsigned int n_q_points = q.n_quadrature_points; - const dFMatrix & shape_values = fe_values.get_shape_values(); - vector dof_values(fe.total_dofs, 0); - cell->get_dof_values (solution, dof_values); - vector psi; // in praxi: first compute @@ -225,9 +220,17 @@ void ProblemBase::integrate_difference (const Function &exact_sol psi); // then subtract finite element // solution - for (unsigned int j=0; j function_values (n_q_points, 0); + fe_values.get_function_values (solution, function_values); + + transform (psi.begin(), psi.end(), + function_values.begin(), + psi.begin(), + minus()); + }; + // for L1_norm and Linfty_norm: // take absolute @@ -296,12 +299,7 @@ void ProblemBase::integrate_difference (const Function &exact_sol // same procedure as above, but now // psi is a vector of gradients - const unsigned int n_dofs = fe.total_dofs; const unsigned int n_q_points = q.n_quadrature_points; - const vector > > & shape_grads = fe_values.get_shape_grads(); - vector dof_values(fe.total_dofs, 0); - cell->get_dof_values (solution, dof_values); - vector > psi; // in praxi: first compute @@ -311,10 +309,16 @@ void ProblemBase::integrate_difference (const Function &exact_sol // then subtract finite element // solution - for (unsigned int j=0; j > function_grads (n_q_points, Point()); + fe_values.get_function_grads (solution, function_grads); + transform (psi.begin(), psi.end(), + function_grads.begin(), + psi.begin(), + minus >()); + }; // take square of integrand vector psi_square (psi.size(), 0.0); for (unsigned int i=0; i::estimate_error (const DoFHandler &dof, // number of integration points per face const unsigned int n_q_points = quadrature.n_quadrature_points; - // number of dofs per cell - const unsigned int n_dofs = fe.total_dofs; // make up a fe face values object for the // restriction of the finite element function @@ -98,27 +96,8 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, // let psi be a short name for // [grad u_h] vector > psi(n_q_points); + fe_face_values_cell.get_function_grads (solution, psi); - // get a list of the values of - // the solution for the ansatz - // functions on this cell - vector dof_values(fe.total_dofs, 0); - cell->get_dof_values (solution, dof_values); - - // get a list of the gradients of - // the ansatz functions on this - // cell at the quadrature points - const vector > > &shape_grads(fe_face_values_cell. - get_shape_grads()); - - // compute the gradients of the solution - // at the quadrature points by summing - // over the ansatz functions. - for (unsigned int j=0; j::estimate_error (const DoFHandler &dof, fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor, fe, boundary); - // get a list of the values of - // the solution for the ansatz - // functions on this neighbor - neighbor->get_dof_values (solution, dof_values); - // get a list of the gradients of the - // - const vector > > &neighbor_grads (fe_face_values_cell. - get_shape_grads()); - // subtract the gradients of the - // solution on the neigbor cell - // at the quadrature points from - // those of the present cell - for (unsigned int j=0; j > neighbor_psi (n_q_points); + fe_face_values_neighbor.get_function_grads (solution, neighbor_psi); + + // compute the jump in the gradients + transform (psi.begin(), psi.end(), + neighbor_psi.begin(), + psi.begin(), + minus >()); }; @@ -188,6 +163,7 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, vector phi(n_q_points,0); const vector > &normal_vectors(fe_face_values_cell. get_normal_vectors()); + for (unsigned int point=0; point