From: Wolfgang Bangerth Date: Sun, 24 May 1998 21:08:52 +0000 (+0000) Subject: Move code closer to the C++ standard and small performance changes. X-Git-Tag: v8.0.0~22912 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=50876aa6b64cb2c33974a108187d82fc45b9e6c3;p=dealii.git Move code closer to the C++ standard and small performance changes. git-svn-id: https://svn.dealii.org/trunk@342 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index 07cbbee068..d6ea8a9169 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -93,11 +93,13 @@ DoFLineAccessor::get_dof_indices (vector &dof_indices) const dof_handler->get_selected_fe().dofs_per_line), ExcVectorDoesNotMatch()); + const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex, + dofs_per_line = dof_handler->get_selected_fe().dofs_per_line; vector::iterator next = dof_indices.begin(); for (unsigned int vertex=0; vertex<2; ++vertex) - for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) + for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) + for (unsigned int d=0; d::get_dof_indices (vector &dof_indices) const dof_handler->get_selected_fe().dofs_per_quad), ExcVectorDoesNotMatch()); + const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex, + dofs_per_line = dof_handler->get_selected_fe().dofs_per_line, + dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad; vector::iterator next = dof_indices.begin(); for (unsigned int vertex=0; vertex<4; ++vertex) - for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) + for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) + for (unsigned int d=0; dline(line)->dof_index(d); - for (unsigned int d=0; dget_selected_fe().dofs_per_quad; ++d) + for (unsigned int d=0; d::child (const unsigned int i) const { +template <> DoFSubstructAccessor<1>::face_iterator DoFCellAccessor<1>::face (const unsigned int) const { Assert (false, ExcNotUsefulForThisDimension()); @@ -307,6 +313,7 @@ DoFCellAccessor<1>::face (const unsigned int) const { +template <> DoFSubstructAccessor<2>::face_iterator DoFCellAccessor<2>::face (const unsigned int i) const { return line(i); @@ -314,6 +321,7 @@ DoFCellAccessor<2>::face (const unsigned int i) const { +template <> void DoFCellAccessor<1>::get_dof_values (const dVector &values, vector &dof_values) const { @@ -324,17 +332,20 @@ DoFCellAccessor<1>::get_dof_values (const dVector &values, Assert (values.size() == dof_handler->n_dofs(), ExcVectorDoesNotMatch()); + const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex, + dofs_per_line = dof_handler->get_selected_fe().dofs_per_line; vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<2; ++vertex) - for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) + for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) + for (unsigned int d=0; d void DoFCellAccessor<2>::get_dof_values (const dVector &values, vector &dof_values) const { @@ -345,14 +356,17 @@ DoFCellAccessor<2>::get_dof_values (const dVector &values, Assert (values.size() == dof_handler->n_dofs(), ExcVectorDoesNotMatch()); + const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex, + dofs_per_line = dof_handler->get_selected_fe().dofs_per_line, + dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad; vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<4; ++vertex) - for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) + for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) + for (unsigned int d=0; dline(line)->dof_index(d)); - for (unsigned int d=0; dget_selected_fe().dofs_per_quad; ++d) + for (unsigned int d=0; d void KellyErrorEstimator<1>::estimate_error (const DoFHandler<1> &, const Quadrature<0> &, const FiniteElement<1> &, @@ -174,6 +175,7 @@ void KellyErrorEstimator::estimate_error (const DoFHandler &dof, +template <> void KellyErrorEstimator<1>::integrate_over_regular_face (const active_cell_iterator &, const unsigned int , const FiniteElement<1> &, @@ -334,6 +336,7 @@ integrate_over_regular_face (const active_cell_iterator &cell, +template <> void KellyErrorEstimator<1>:: integrate_over_irregular_face (const active_cell_iterator &, const unsigned int , diff --git a/deal.II/deal.II/source/numerics/matrices.cc b/deal.II/deal.II/source/numerics/matrices.cc index 2be140e98f..ad1723e1ef 100644 --- a/deal.II/deal.II/source/numerics/matrices.cc +++ b/deal.II/deal.II/source/numerics/matrices.cc @@ -78,6 +78,7 @@ void MatrixCreator::create_mass_matrix (const DoFHandler &dof, +template <> void MatrixCreator<1>::create_boundary_mass_matrix (const DoFHandler<1> &, const FiniteElement<1> &, const Quadrature<0> &, @@ -322,13 +323,14 @@ void MatrixTools::apply_boundary_values (const map &boundary_va for (; dof != endd; ++dof) { + const int dof_number = dof->first; // for each boundary dof: - + // set entries of this line // to zero - for (unsigned int j=sparsity_rowstart[dof->first]; - jfirst+1]; ++j) - if (sparsity_colnums[j] != dof->first) + const unsigned int last = sparsity_rowstart[dof_number+1]; + for (unsigned int j=sparsity_rowstart[dof->first]; j::apply_boundary_values (const map &boundary_va // store the new rhs entry to make // the gauss step more efficient double new_rhs; - if (matrix.diag_element(dof->first) != 0.0) - new_rhs = right_hand_side(dof->first) - = dof->second * matrix.diag_element(dof->first); + if (matrix.diag_element(dof_number) != 0.0) + new_rhs = right_hand_side(dof_number) + = dof->second * matrix.diag_element(dof_number); else { double first_diagonal_entry = 1; @@ -359,23 +361,23 @@ void MatrixTools::apply_boundary_values (const map &boundary_va break; }; - matrix.set(dof->first, dof->first, + matrix.set(dof_number, dof_number, first_diagonal_entry); - new_rhs = right_hand_side(dof->first) + new_rhs = right_hand_side(dof_number) = dof->second * first_diagonal_entry; }; // store the only nonzero entry // of this line for the Gauss // elimination step - const double diagonal_entry = matrix.diag_element(dof->first); + const double diagonal_entry = matrix.diag_element(dof_number); // do the Gauss step for (unsigned int row=0; rowfirst) && - ((signed int)row != dof->first)) + if ((sparsity_colnums[j] == dof_number) && + ((signed int)row != dof_number)) // this line has an entry // in the regarding column // but this is not the main @@ -391,7 +393,7 @@ void MatrixTools::apply_boundary_values (const map &boundary_va // preset solution vector - solution(dof->first) = dof->second; + solution(dof_number) = dof->second; }; }; @@ -420,23 +422,26 @@ void MassMatrix::assemble (dFMatrix &cell_matrix, const dFMatrix &values = fe_values.get_shape_values (); const vector &weights = fe_values.get_JxW_values (); + const unsigned int total_dofs = fe_values.total_dofs, + n_q_points = fe_values.n_quadrature_points; + if (coefficient != 0) { vector coefficient_values (fe_values.n_quadrature_points); coefficient->value_list (fe_values.get_quadrature_points(), coefficient_values); - for (unsigned int point=0; point::assemble (dFMatrix &cell_matrix, vector rhs_values (fe_values.n_quadrature_points); right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values); + const unsigned int total_dofs = fe_values.total_dofs, + n_q_points = fe_values.n_quadrature_points; + if (coefficient != 0) { - vector coefficient_values (fe_values.n_quadrature_points); + vector coefficient_values (n_q_points); coefficient->value_list (fe_values.get_quadrature_points(), coefficient_values); - for (unsigned int point=0; point::assemble (dFMatrix &cell_matrix, }; } else - for (unsigned int point=0; point::assemble (dVector &rhs, vector rhs_values(fe_values.n_quadrature_points); right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values); - for (unsigned int point=0; point::assemble (dFMatrix &cell_matrix, const vector &weights = fe_values.get_JxW_values (); right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values); + const unsigned int total_dofs = fe_values.total_dofs, + n_q_points = fe_values.n_quadrature_points; + if (coefficient != 0) { - vector coefficient_values(fe_values.n_quadrature_points); + vector coefficient_values(n_q_points); coefficient->value_list (fe_values.get_quadrature_points(), coefficient_values); - for (unsigned int point=0; point::assemble (dFMatrix &cell_matrix, }; } else - for (unsigned int point=0; point::assemble (dFMatrix &cell_matrix, const vector > >&gradients = fe_values.get_shape_grads (); const vector &weights = fe_values.get_JxW_values (); + const unsigned int total_dofs = fe_values.total_dofs, + n_q_points = fe_values.n_quadrature_points; + if (coefficient != 0) { - vector coefficient_values(fe_values.n_quadrature_points); + vector coefficient_values(n_q_points); coefficient->value_list (fe_values.get_quadrature_points(), coefficient_values); - for (unsigned int point=0; point::assemble (dVector &rhs, vector rhs_values(fe_values.n_quadrature_points); right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values); - for (unsigned int point=0; point::interpolate (const DoFHandler &dof, +template <> void VectorTools<1>::project (const DoFHandler<1> &, const ConstraintMatrix &, const FiniteElement<1> &, @@ -173,6 +174,7 @@ void VectorTools::project (const DoFHandler &dof, +template <> void VectorTools<1>::interpolate_boundary_values (const DoFHandler<1> &, const FunctionMap &, @@ -200,7 +202,7 @@ VectorTools::interpolate_boundary_values (const DoFHandler &dof, DoFHandler::active_face_iterator face = dof.begin_active_face(), endf = dof.end_face(); - FunctionMap::const_iterator function_ptr; + typename FunctionMap::const_iterator function_ptr; // field to store the indices of dofs vector face_dofs (fe.dofs_per_face);