From: Wolfgang Bangerth Date: Tue, 28 Apr 1998 12:30:00 +0000 (+0000) Subject: Fix small problems. X-Git-Tag: v8.0.0~23043 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d1b212af4359b0670ae4a33d595c1df9b8e1e3f7;p=dealii.git Fix small problems. git-svn-id: https://svn.dealii.org/trunk@211 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/numerics/base.cc b/deal.II/deal.II/source/numerics/base.cc index cf4fad227e..da58fa61af 100644 --- a/deal.II/deal.II/source/numerics/base.cc +++ b/deal.II/deal.II/source/numerics/base.cc @@ -214,7 +214,7 @@ void ProblemBase::integrate_difference (const Function &exact_sol 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; + vector dof_values(fe.total_dofs, 0); cell->get_dof_values (solution, dof_values); vector psi; @@ -299,7 +299,7 @@ void ProblemBase::integrate_difference (const Function &exact_sol 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; + vector dof_values(fe.total_dofs, 0); cell->get_dof_values (solution, dof_values); vector > psi; @@ -359,7 +359,7 @@ void ProblemBase::fill_data (DataOut &out) const { template pair ProblemBase::get_solution_name () const { - return pair("solution", ""); + return pair("solution", ""); }; diff --git a/deal.II/deal.II/source/numerics/error_estimator.cc b/deal.II/deal.II/source/numerics/error_estimator.cc index 26250d2fe0..26f568a07c 100644 --- a/deal.II/deal.II/source/numerics/error_estimator.cc +++ b/deal.II/deal.II/source/numerics/error_estimator.cc @@ -59,20 +59,25 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, // to a face, for the present cell and its // neighbor. FEFaceValues fe_face_values_cell (fe, quadrature, - UpdateFlags(update_gradients | update_JxW_values | - update_jacobians | update_q_points | + UpdateFlags(update_gradients | + update_JxW_values | + update_jacobians | + update_q_points | update_normal_vectors)); FEFaceValues fe_face_values_neighbor (fe, quadrature, - UpdateFlags(update_gradients)); + UpdateFlags(update_gradients | + update_jacobians)); + + // loop variables DoFHandler::active_cell_iterator cell = dof.begin_active(), endc = dof.end(); - // loop over all cells for (unsigned int present_cell=0; cell!=endc; ++cell, ++present_cell) // loop over all faces of this cell for (unsigned int face_no=0; face_no<2*dim; ++face_no) { - const unsigned char boundary_indicator = cell->face(face_no)->boundary_indicator(); + const unsigned char boundary_indicator + = cell->face(face_no)->boundary_indicator(); if ((boundary_indicator != 255) && neumann_bc.find(boundary_indicator)==neumann_bc.end()) // this face is part of the boundary @@ -97,13 +102,14 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, // get a list of the values of // the solution for the ansatz // functions on this cell - vector dof_values; + 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()); + const vector > > &shape_grads(fe_face_values_cell. + get_shape_grads()); // compute the gradients of the solution // at the quadrature points by summing @@ -123,7 +129,8 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, Assert (cell->neighbor(face_no).state() == valid, ExcInternalError()); unsigned int neighbor_neighbor; - DoFHandler::active_cell_iterator neighbor = cell->neighbor(face_no); + DoFHandler::active_cell_iterator neighbor + = cell->neighbor(face_no); // find which number the current // face has relative to the neighboring @@ -137,7 +144,8 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, // get restriction of finite element // function of #neighbor# to the // common face. - fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor, fe, boundary); + fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor, + fe, boundary); // get a list of the values of // the solution for the ansatz @@ -178,7 +186,8 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, // // let phi be the name of the integrand vector phi(n_q_points,0); - const vector > &normal_vectors(fe_face_values_cell.get_normal_vectors()); + const vector > &normal_vectors(fe_face_values_cell. + get_normal_vectors()); for (unsigned int point=0; point::estimate_error (const DoFHandler &dof, error(present_cell) += sqrt(inner_product (phi.begin(), phi.end(), fe_face_values_cell.get_JxW_values().begin(), - 0.0)); + 0.0)) * cell->diameter(); }; };